MixtureTmpl.h

#ifndef PSCF_MIXTURE_TMPL_H
#define PSCF_MIXTURE_TMPL_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 <pscf/chem/MixtureBase.h>
#include <util/param/ParamComposite.h>
#include <util/containers/DArray.h>
 
namespace Pscf
{
 
   using namespace Util;

   template <class PT, class ST>
   class MixtureTmpl : public MixtureBase, public ParamComposite
   {
   public:
 
      // Public typedefs

      using SolventT = ST;

      using PolymerT = PT;

      using BlockT = typename PT::BlockT;

      using PropagatorT = typename BlockT::PropagatorT;
 
      // Public member functions

      MixtureTmpl();

      ~MixtureTmpl();

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

 
      /*
      * Public member functions inherited from MixtureBase:
      *
      *   int nMonomer() const;
      *   int nPolymer() const;
      *   int nSolvent() const;
      *   int nBlock() const;
      *   double vMonomer() const;
      *   Monomer const & monomer(int id) const;
      *
      *   void setVmonomer(double);
      */

      PolymerT& polymer(int id);

      PolymerT const & polymer(int id) const;

      PolymerSpecies const & polymerSpecies(int id) const final;

      SolventT& solvent(int id);

      SolventT const & solvent(int id) const;

      SolventSpecies const & solventSpecies(int id) const final;

 
   private:

      DArray<PolymerT> polymers_;

      DArray<SolventT> solvents_;
 
      // Restrict access to inherited protected data
      using MixtureBase::monomers_;
      using MixtureBase::nMonomer_;
      using MixtureBase::nPolymer_;
      using MixtureBase::nSolvent_;
      using MixtureBase::nBlock_;
      using MixtureBase::vMonomer_;
 
   };
 
   // Inline member functions
 
   template <class PT, class ST>
   inline PT& MixtureTmpl<PT,ST>::polymer(int id)
   {  
      UTIL_CHECK(id < nPolymer_);
      return polymers_[id];
   }
 
   template <class PT, class ST>
   inline PT const & MixtureTmpl<PT,ST>::polymer(int id) const
   {  
      UTIL_CHECK(id < nPolymer_);
      return polymers_[id];
   }
 
   template <class PT, class ST>
   inline ST& MixtureTmpl<PT,ST>::solvent(int id)
   {  
      UTIL_CHECK(id < nSolvent_);
      return solvents_[id]; 
   }
 
   template <class PT, class ST>
   inline ST const & MixtureTmpl<PT,ST>::solvent(int id) const
   {  
      UTIL_CHECK(id < nSolvent_);
      return solvents_[id]; 
   }
 
   // Non-inline member functions
 
   /*
   * Constructor.
   */
   template <class PT, class ST>
   MixtureTmpl<PT,ST>::MixtureTmpl()
    : MixtureBase(),
      ParamComposite(),
      polymers_(),
      solvents_()
   {}
 
   /*
   * Destructor.
   */
   template <class PT, class ST>
   MixtureTmpl<PT,ST>::~MixtureTmpl()
   {}
 
   /*
   * Get a PolymerSpecies descriptor by non-const reference.
   */
   template <class PT, class ST>
   PolymerSpecies const & MixtureTmpl<PT,ST>::polymerSpecies(int id) const
   {  
      UTIL_CHECK(id < nPolymer_);
      return polymers_[id];
   }
 
   /*
   * Get a SolventSpecies descriptor by const reference.
   */
   template <class PT, class ST>
   SolventSpecies const & MixtureTmpl<PT,ST>::solventSpecies(int id) const
   {  
      UTIL_CHECK(id < nSolvent_);
      return solvents_[id]; 
   }
 
   /*
   * Read all parameters and initialize.
   */
   template <class PT, class ST>
   void MixtureTmpl<PT,ST>::readParameters(std::istream& in)
   {
      // Read nMonomer and monomers array
      read<int>(in, "nMonomer", nMonomer_);
      monomers_.allocate(nMonomer_);
      for (int i = 0; i < nMonomer_; ++i) {
         monomers_[i].setId(i);
      }
      readDArray< Monomer >(in, "monomers", monomers_, nMonomer_);
 
      /*
      * The input format for a single monomer is defined in the istream
      * extraction operation (operator >>) for a Pscf::Monomer, in file
      * pscf/chem/Monomer.cpp. The text representation contains only the
      * value for the monomer statistical segment Monomer::kuhn.
      */
 
      // Read nPolymer (required parameter, must be > 0)
      read<int>(in, "nPolymer", nPolymer_);
      UTIL_CHECK(nPolymer_ > 0);
 
      // Optionally read nSolvent, with nSolvent=0 by default
      nSolvent_ = 0;
      readOptional<int>(in, "nSolvent", nSolvent_);
 
      // Read polymers and compute nBlock
      nBlock_ = 0;
      if (nPolymer_ > 0) {
 
         polymers_.allocate(nPolymer_);
         for (int i = 0; i < nPolymer_; ++i) {
            readParamComposite(in, polymer(i));
            nBlock_ = nBlock_ + polymer(i).nBlock();
         }
   
         // Set statistical segment lengths for all blocks
         double kuhn;
         int monomerId;
         for (int i = 0; i < nPolymer_; ++i) {
            for (int j = 0; j < polymer(i).nBlock(); ++j) {
               monomerId = polymer(i).block(j).monomerId();
               kuhn = monomer(monomerId).kuhn();
               polymer(i).block(j).setKuhn(kuhn);
            }
         }
 
      }
 
      // Read solvents
      if (nSolvent_ > 0) {
 
         solvents_.allocate(nSolvent_);
         for (int i = 0; i < nSolvent_; ++i) {
            readParamComposite(in, solvent(i));
         }
 
      }
 
      // Optionally read monomer reference value
      vMonomer_ = 1.0; // Default value
      readOptional(in, "vMonomer", vMonomer_);
 
   }
 
}
#endif