rpg/fts/simulator/Simulator.h

#ifndef RPG_SIMULATOR_H
#define RPG_SIMULATOR_H
 
/*
* PSCF - Polymer Self-Consistent Field
*
* Copyright 2015 - 2025, The Regents of the University of Minnesota
* Distributed under the terms of the GNU General Public License.
*/
 
#include <util/param/ParamComposite.h>     // base class
 
#include <rpg/fts/simulator/SimState.h>    // member
#include <prdc/cuda/RField.h>              // memmber (template arg)
#include <util/random/Random.h>            // member
#include <pscf/cuda/CudaRandom.h>          // member
#include <util/containers/DArray.h>        // member (template)
#include <util/containers/DMatrix.h>       // member (template)
 
namespace Pscf {
namespace Rpg {
 
   template <int D> class System;
   template <int D> class Compressor;
   template <int D> class CompressorFactory;
   template <int D> class Perturbation;
   template <int D> class PerturbationFactory;
   template <int D> class Ramp;
   template <int D> class RampFactory;
 
   using namespace Util;
   using namespace Pscf::Prdc::Cuda;

   template <int D>
   class Simulator : public ParamComposite
   {
 
   public:

      Simulator(System<D>& system);

      ~Simulator();

      void allocate();

      virtual void readParameters(std::istream &in);


      virtual void simulate(int nStep);

      virtual void analyze(int min, int max,
                           std::string classname,
                           std::string filename);

      void clearData();

      virtual void outputTimers(std::ostream& out) const;

      virtual void outputMdeCounter(std::ostream& out) const;

      virtual void clearTimers();

      long iStep();

      long iTotalStep();


      void analyzeChi();

      DArray<double> const & chiEvals() const;

      double chiEval(int a) const;

      DMatrix<double> const & chiEvecs() const;

      double chiEvecs(int a, int i) const;

      DArray<double> const & sc() const;

      double sc(int a) const;


      void computeHamiltonian();

      double hamiltonian() const;

      double idealHamiltonian() const;

      double fieldHamiltonian() const;

      double perturbationHamiltonian() const;

      bool hasHamiltonian() const;


      void computeWc();

      DArray< RField<D> > const & wc() const;

      RField<D> const & wc(int a) const;

      bool hasWc() const;


      void computeCc();

      DArray< RField<D> > const & cc() const;

      RField<D> const & cc(int a) const;

      bool hasCc() const;


      void computeDc();

      DArray< RField<D> > const & dc() const;

      RField<D> const & dc(int i) const;

      bool hasDc() const;


      void saveState();

      void restoreState();

      void clearState();


      System<D>& system();

      Compressor<D>& compressor();

      Compressor<D> const & compressor() const;

      bool hasCompressor() const;

      Random& random();

      CudaRandom& cudaRandom();

      bool hasPerturbation() const;

      Perturbation<D> const & perturbation() const;

      Perturbation<D>& perturbation();

      bool hasRamp() const;

      Ramp<D> const & ramp() const;

      Ramp<D>& ramp();

 
   protected:
 
      // Protected member functions
 
      using Util::ParamComposite::setClassName;

      void readRandomSeed(std::istream& in);

      CompressorFactory<D>& compressorFactory();

      void readCompressor(std::istream& in, bool& isEnd);

      PerturbationFactory<D>& perturbationFactory();

      void readPerturbation(std::istream& in, bool& isEnd);

      void setPerturbation(Perturbation<D>* ptr);

      RampFactory<D>& rampFactory();

      void readRamp(std::istream& in, bool& isEnd);

      void setRamp(Ramp<D>* ptr);
 
      // Protected data members

      Random random_;

      CudaRandom cudaRandom_;

      DArray< RField<D> > wc_;

      DArray< RField<D> > cc_;

      DArray< RField<D> > dc_;

      mutable SimState<D> state_;

      double hamiltonian_;

      double idealHamiltonian_;

      double fieldHamiltonian_;

      double perturbationHamiltonian_;

      long iStep_;

      long iTotalStep_;

      long seed_;

      bool hasHamiltonian_;

      bool hasWc_;

      bool hasCc_;

      bool hasDc_;
 
   private:

      DMatrix<double> chiP_;

      DMatrix<double> chiEvecs_;

      DArray<double>  chiEvals_;

      DArray<double>  sc_;

      RField<D> wcs_;

      System<D>* systemPtr_;

      CompressorFactory<D>* compressorFactoryPtr_;

      Compressor<D>* compressorPtr_;

      PerturbationFactory<D>* perturbationFactoryPtr_;

      Perturbation<D>* perturbationPtr_;

      RampFactory<D>* rampFactoryPtr_;

      Ramp<D>* rampPtr_;

      bool isAllocated_;
 
   };
 
   // Inline functions
 
   // Get the parent System.
   template <int D>
   inline System<D>& Simulator<D>::system()
   {
      UTIL_CHECK(systemPtr_);
      return *systemPtr_;
   }
 
   // Get the CPU random number generator.
   template <int D>
   inline Random& Simulator<D>::random()
   {  return random_; }
 
   // Get the GPU random number generator.
   template <int D>
   inline CudaRandom& Simulator<D>::cudaRandom()
   {  return cudaRandom_; }
 
   // Does this Simulator have a Compressor?
   template <int D>
   inline bool Simulator<D>::hasCompressor() const
   {  return (bool)compressorPtr_; }
 
   // Get the Compressor by non-const reference.
   template <int D>
   inline Compressor<D>& Simulator<D>::compressor()
   {
      UTIL_CHECK(compressorPtr_);
      return *compressorPtr_;
   }
 
   // Get the Compressor by const reference.
   template <int D>
   inline Compressor<D> const & Simulator<D>::compressor() const
   {
      UTIL_CHECK(compressorPtr_);
      return *compressorPtr_;
   }
 
   // Get the Compressor factory.
   template <int D>
   inline CompressorFactory<D>& Simulator<D>::compressorFactory()
   {
      UTIL_CHECK(compressorFactoryPtr_);
      return *compressorFactoryPtr_;
   }
 
   // Does this Simulator have an associated Perturbation?
   template <int D>
   inline bool Simulator<D>::hasPerturbation() const
   {  return (bool)perturbationPtr_; }
 
   // Get the perturbation (if any) by const reference.
   template <int D>
   inline Perturbation<D> const & Simulator<D>::perturbation() const
   {
      UTIL_CHECK(perturbationPtr_);
      return *perturbationPtr_;
   }
 
   // Get the perturbation (if any) by non-const reference.
   template <int D>
   inline Perturbation<D>& Simulator<D>::perturbation()
   {
      UTIL_CHECK(perturbationPtr_);
      return *perturbationPtr_;
   }
 
   // Get the perturbation factory.
   template <int D>
   inline PerturbationFactory<D>& Simulator<D>::perturbationFactory()
   {
      UTIL_CHECK(perturbationFactoryPtr_);
      return *perturbationFactoryPtr_;
   }
 
   // Does this Simulator have an associated Ramp?
   template <int D>
   inline bool Simulator<D>::hasRamp() const
   {  return (bool)rampPtr_; }
 
   // Get the ramp by const reference.
   template <int D>
   inline Ramp<D> const & Simulator<D>::ramp() const
   {
      UTIL_CHECK(rampPtr_);
      return *rampPtr_;
   }
 
   // Get the ramp by non-const reference.
   template <int D>
   inline Ramp<D>& Simulator<D>::ramp()
   {
      UTIL_CHECK(rampPtr_);
      return *rampPtr_;
   }
 
   // Get the ramp factory.
   template <int D>
   inline RampFactory<D>& Simulator<D>::rampFactory()
   {
      UTIL_CHECK(rampFactoryPtr_);
      return *rampFactoryPtr_;
   }
 
   // Return an array of eigenvalues of projected chi matrix.
   template <int D>
   inline double Simulator<D>::chiEval(int a) const
   {  return chiEvals_[a]; }
 
   // Return an array of eigenvalues of projected chi matrix.
   template <int D>
   inline DArray<double> const & Simulator<D>::chiEvals() const
   {  return chiEvals_; }
 
   // Return a matrix of eigenvectors of the projected chi matrix.
   template <int D>
   inline DMatrix<double> const & Simulator<D>::chiEvecs() const
   {  return chiEvecs_; }
 
   // Return a matrix of eigenvectors of the projected chi matrix.
   template <int D>
   inline double Simulator<D>::chiEvecs(int a, int i) const
   {  return chiEvecs_(a, i); }
 
   // Return array of values of vector S.
   template <int D>
   inline DArray<double> const & Simulator<D>::sc() const
   {  return sc_; }
 
   // Return one component of vector S.
   template <int D>
   inline double Simulator<D>::sc(int a) const
   {  return sc_[a]; }
 
   // Has the Hamiltonian been computed for the current w fields ?
   template <int D>
   inline bool Simulator<D>::hasHamiltonian() const
   {  return hasHamiltonian_; }
 
   // Get the precomputed Hamiltonian
   template <int D>
   inline double Simulator<D>::hamiltonian() const
   {
      UTIL_CHECK(hasHamiltonian_);
      return hamiltonian_;
   }
 
   // Get the ideal gas component of the precomputed Hamiltonian
   template <int D>
   inline double Simulator<D>::idealHamiltonian() const
   {
      UTIL_CHECK(hasHamiltonian_);
      return idealHamiltonian_;
   }
 
   // Get the W field component of the precomputed Hamiltonian.
   template <int D>
   inline double Simulator<D>::fieldHamiltonian() const
   {
      UTIL_CHECK(hasHamiltonian_);
      return fieldHamiltonian_;
   }
 
   // Get the perturbation component of the precomputed Hamiltonian.
   template <int D>
   inline double Simulator<D>::perturbationHamiltonian() const
   {
      UTIL_CHECK(hasHamiltonian_);
      return perturbationHamiltonian_;
   }
 
   // Return all eigencomponents of the w fields.
   template <int D>
   inline DArray< RField<D> > const & Simulator<D>::wc() const
   {  return wc_; }
 
   // Return a single eigenvector component of the w fields.
   template <int D>
   inline RField<D> const & Simulator<D>::wc(int a) const
   {  return wc_[a]; }
 
   // Have eigenvector components of current w fields been computed?
   template <int D>
   inline bool Simulator<D>::hasWc() const
   {  return hasWc_; }
 
   // Return all eigenvector components of the current c fields.
   template <int D>
   inline DArray< RField<D> > const & Simulator<D>::cc() const
   {  return cc_; }
 
   // Return a single eigenvector component of the current c fields.
   template <int D>
   inline RField<D> const & Simulator<D>::cc(int a) const
   {  return cc_[a]; }
 
   // Have eigenvector components of current c fields been computed?
   template <int D>
   inline bool Simulator<D>::hasCc() const
   {  return hasCc_; }
 
   // Return all eigenvector components of the current d fields.
   template <int D>
   inline DArray< RField<D> > const & Simulator<D>::dc() const
   {  return dc_; }
 
   // Return a single eigenvector component of the current d fields.
   template <int D>
   inline RField<D> const & Simulator<D>::dc(int a) const
   {  return dc_[a]; }
 
   // Have eigenvector components of current d fields been computed?
   template <int D>
   inline bool Simulator<D>::hasDc() const
   {  return hasDc_; }
 
   // Clear all data (eigen-components of w field and Hamiltonian)
   template <int D>
   inline void Simulator<D>::clearData()
   {
      hasHamiltonian_ = false;
      hasWc_ = false;
      hasCc_ = false;
      hasDc_ = false;
   }
 
   // Return the current converged simulation step index.
   template <int D>
   inline long Simulator<D>::iStep()
   {  return iStep_; }
 
   // Return the current simulation step index.
   template <int D>
   inline long Simulator<D>::iTotalStep()
   {  return iTotalStep_; }
 
   #ifndef RPG_SIMULATOR_TPP
   // Suppress implicit instantiation
   extern template class Simulator<1>;
   extern template class Simulator<2>;
   extern template class Simulator<3>;
   #endif
 
}
}
#endif