rp/field/FieldIo.h

#ifndef RP_FIELD_IO_H
#define RP_FIELD_IO_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 <prdc/crystal/UnitCell.h>     // nested LatticeSystem enum
#include <pscf/math/IntVec.h>          // template with default
#include <util/containers/DArray.h>    // member
 
// Forward declarations for classes used in interfaces
namespace Util {
   class FileMaster;
   template <typename T> class DMatrix;
}
namespace Pscf {
   template <int D> class Mesh;
   namespace Prdc {
      template <int D> class Basis;
      template <int D> class SpaceGroup;
   }
}
 
namespace Pscf {
namespace Rp {
 
   using namespace Util;
   using namespace Prdc;

   template <int D, class RFT, class KFT, class FFT>
   class FieldIo
   {
 
   public:


      void associate(Mesh<D> const & mesh,
                     FFT const & fft,
                     typename UnitCell<D>::LatticeSystem const & lattice,
                     bool const & hasGroup,
                     std::string const & groupName,
                     SpaceGroup<D> const & group,
                     Basis<D> & basis);

      void setFileMaster(FileMaster const & fileMaster);

      void setNMonomer(int nMonomer);


      void
      readFieldsBasis(std::istream& in,
            DArray< DArray<double> >& fields,
                      UnitCell<D> & unitCell) const;

      void readFieldsBasis(std::string filename,
                           DArray< DArray<double> >& fields,
                           UnitCell<D> & unitCell) const;

      void readFieldBasis(std::istream& in,
                DArray<double>& field,
                          UnitCell<D> & unitCell) const;

      void readFieldBasis(std::string filename,
                DArray<double>& field,
                          UnitCell<D> & unitCell) const;

      void writeFieldsBasis(std::ostream& out,
                            DArray< DArray<double> > const & fields,
                            UnitCell<D> const & unitCell) const;

      void writeFieldsBasis(std::string filename,
                            DArray< DArray<double> > const & fields,
                            UnitCell<D> const & unitCell) const;

      void writeFieldBasis(std::ostream& out,
                           DArray<double> const & field,
                           UnitCell<D> const & unitCell) const;

      void writeFieldBasis(std::string filename,
                           DArray<double> const & field,
                           UnitCell<D> const & unitCell) const;


      virtual
      bool readFieldsRGrid(std::istream& in,
                           DArray<RFT>& fields,
                           UnitCell<D> & unitCell) const = 0;

      bool readFieldsRGrid(std::string filename,
                           DArray<RFT>& fields,
                           UnitCell<D> & unitCell) const;

      virtual
      void readFieldsRGridData(std::istream& in,
                               DArray<RFT>& fields,
                               int nMonomer) const = 0;

      virtual
      bool readFieldRGrid(std::istream &in,
                          RFT & field,
                          UnitCell<D>& unitCell) const = 0;

      bool readFieldRGrid(std::string filename,
                          RFT & field,
                          UnitCell<D>& unitCell) const;

      virtual
      void writeFieldsRGrid(std::ostream& out,
                            DArray<RFT> const & fields,
                            UnitCell<D> const & unitCell,
                            bool writeHeader = true,
                            bool isSymmetric = true,
                            bool writeMeshSize = true) const = 0;

      void writeFieldsRGrid(std::string filename,
                            DArray<RFT> const & fields,
                            UnitCell<D> const & unitCell,
                            bool isSymmetric = true) const;

      virtual
      void writeFieldRGrid(std::ostream &out,
                           RFT const & field,
                           UnitCell<D> const & unitCell,
                           bool writeHeader = true,
                           bool isSymmetric = true) const = 0;

      void writeFieldRGrid(std::string filename,
                           RFT const & field,
                           UnitCell<D> const & unitCell,
                           bool isSymmetric = true) const;


      virtual
      void readFieldsKGrid(std::istream& in,
                           DArray<KFT>& fields,
                           UnitCell<D> & unitCell) const = 0;

      void readFieldsKGrid(std::string filename,
                           DArray<KFT>& fields,
                           UnitCell<D> & unitCell) const;

      virtual
      void writeFieldsKGrid(std::ostream& out,
                            DArray<KFT> const & fields,
                            UnitCell<D> const & unitCell,
                            bool isSymmetric = true) const = 0;

      void writeFieldsKGrid(std::string filename,
                           DArray<KFT> const & fields,
                           UnitCell<D> const & unitCell,
                           bool isSymmetric = true) const;


      virtual
      void convertBasisToKGrid(DArray<double> const & components,
                               KFT& dft) const = 0;

      void convertBasisToKGrid(DArray< DArray<double> > const & in,
                               DArray<KFT>& out) const;

      virtual
      void convertKGridToBasis(KFT const & in,
                               DArray<double> & out,
                               bool checkSymmetry = true,
                               double epsilon = 1.0e-8) const = 0;

      void convertKGridToBasis(DArray<KFT> const & in,
                               DArray< DArray<double> > & out,
                               bool checkSymmetry = true,
                               double epsilon = 1.0e-8) const;

      void convertBasisToRGrid(DArray<double> const & in,
                               RFT & out) const;

      void convertBasisToRGrid(DArray< DArray<double> > const & in,
                               DArray<RFT> & out) const ;

      void convertRGridToBasis(RFT const & in,
                               DArray<double> & out,
                               bool checkSymmetry = true,
                               double epsilon = 1.0e-8) const;

      void convertRGridToBasis(DArray<RFT> const & in,
                               DArray< DArray<double> > & out,
                               bool checkSymmetry = true,
                               double epsilon = 1.0e-8) const;

      void convertKGridToRGrid(DArray<KFT> const & in,
                               DArray<RFT> & out) const;

      void convertKGridToRGrid(KFT const & in,
                               RFT & out) const;

      void convertRGridToKGrid(DArray<RFT> const & in,
                               DArray<KFT> & out) const;

      void convertRGridToKGrid(RFT const & in,
                               KFT & out) const;


      void convertBasisToRGrid(std::string const & inFileName,
                               std::string const & outFileName) const;

      void convertRGridToBasis(std::string const & inFileName,
                               std::string const & outFileName) const;

      void convertKGridToRGrid(std::string const & inFileName,
                               std::string const & outFileName) const;

      void convertRGridToKGrid(std::string const & inFileName,
                               std::string const & outFileName) const;

      void convertKGridToBasis(std::string const & inFileName,
                               std::string const & outFileName) const;

      void convertBasisToKGrid(std::string const & inFileName,
                               std::string const & outFileName) const;


      virtual
      bool hasSymmetry(KFT const & in,
                       double epsilon = 1.0e-8,
                       bool verbose = true) const = 0;

      bool hasSymmetry(RFT const & in,
                       double epsilon = 1.0e-8,
                       bool verbose = true) const;

      bool hasSymmetry(std::string const & inFileName,
                       double epsilon = 1.0E-8) const;


      void compareFieldsBasis(DArray< DArray<double> > const & field1,
                              DArray< DArray<double> > const & field2)
      const;

      void compareFieldsBasis(std::string const & filename1,
                              std::string const & filename2) const;

      virtual
      void compareFieldsRGrid(DArray< RFT > const & field1,
                              DArray< RFT > const & field2) const = 0;

      void compareFieldsRGrid(std::string const & filename1,
                              std::string const & filename2) const;


      virtual
      void scaleFieldBasis(DArray<double>& field,
                           double factor) const;

      void scaleFieldsBasis(DArray< DArray<double> >& fields,
                            double factor) const;

      void scaleFieldsBasis(std::string const & inFileName,
                            std::string const & outFileName,
                            double factor) const;

      virtual
      void scaleFieldRGrid(RFT& field, double factor) const = 0;

      void scaleFieldsRGrid(DArray<RFT> & fields, double factor) const;

      void scaleFieldsRGrid(std::string const & inFileName,
                            const std::string & outFileName,
                            double factor) const;


      void estimateWBasis(DMatrix<double> const & chi,
                          DArray< DArray<double> > & fields) const;

      void estimateWBasis(std::string const & inFileName,
                          std::string const & outFileName,
                          DMatrix<double> const & chi) const;


      virtual
      void replicateUnitCell(std::ostream& out,
                             DArray<RFT> const & fields,
                             UnitCell<D> const & unitCell,
                             IntVec<D> const & replicas) const = 0;

      void replicateUnitCell(std::string filename,
                             DArray<RFT> const & fields,
                             UnitCell<D> const & unitCell,
                             IntVec<D> const & replicas) const;

      void replicateUnitCell(std::string const & inFileName,
                             std::string const & outFileName,
                             IntVec<D> const & replicas) const;


      virtual
      void expandRGridDimension(std::ostream &out,
                                DArray<RFT > const & fields,
                                UnitCell<D> const & unitCell,
                                int d,
                                DArray<int> const& newGridDimensions)
      const = 0;

      void expandRGridDimension(std::string filename,
                                DArray<RFT > const & fields,
                                UnitCell<D> const & unitCell,
                                int d,
                                DArray<int> newGridDimensions) const;

      void expandRGridDimension(std::string const & inFileName,
                                std::string const & outFileName,
                                int d,
                                DArray<int> newGridDimensions) const;


      void readFieldHeader(std::istream& in,
                           int& nMonomer,
                           UnitCell<D> & unitCell,
                           bool & isSymmetric) const;

      void writeFieldHeader(std::ostream& out,
                            int nMonomer,
                            UnitCell<D> const & unitCell,
                            bool isSymmetric = true) const;


      Mesh<D> const & mesh() const
      {
         UTIL_ASSERT(meshPtr_);
         return *meshPtr_;
      }

      Basis<D> const & basis() const
      {
         UTIL_ASSERT(basisPtr_);
         return *basisPtr_;
      }

      FileMaster const & fileMaster() const
      {
         UTIL_ASSERT(fileMasterPtr_);
         return *fileMasterPtr_;
      }

 
   protected:

      FieldIo();

      ~FieldIo() = default;

      typename UnitCell<D>::LatticeSystem const & lattice() const
      {
         UTIL_ASSERT(latticePtr_);
         return *latticePtr_;
      }

      bool hasGroup() const
      {
         UTIL_ASSERT(hasGroupPtr_);
         return *hasGroupPtr_;
      }

      std::string const & groupName() const
      {
         UTIL_ASSERT(groupNamePtr_);
         return *groupNamePtr_;
      }

      SpaceGroup<D> const & group() const
      {
         UTIL_ASSERT(groupPtr_);
         return *groupPtr_;
      }

      FFT const & fft() const
      {
         UTIL_ASSERT(fftPtr_);
         return *fftPtr_;
      }

      void checkAllocateRGrid() const;

      void checkAllocateKGrid() const;

      void checkAllocateBasis(std::string const & inFileName) const;
 
   private:

      Mesh<D> const * meshPtr_;

      FFT const * fftPtr_;

      typename UnitCell<D>::LatticeSystem const * latticePtr_;

      bool const * hasGroupPtr_;

      std::string const * groupNamePtr_;

      SpaceGroup<D> const * groupPtr_;

      Basis<D> * basisPtr_;

      FileMaster const * fileMasterPtr_;

      int nMonomer_;
 
      // Mutable work space used by functions that read field files

      mutable DArray< DArray<double> > tmpFieldsBasis_;

      mutable DArray<RFT> tmpFieldsRGrid_;

      mutable DArray<KFT> tmpFieldsKGrid_;

      mutable KFT workDft_;

      mutable bool isAllocatedBasis_;

      mutable bool isAllocatedRGrid_;

      mutable bool isAllocatedKGrid_;
 
   };
 
} // namespace Rp
} // namespace Pscf
#endif