PSCF v1.3
Pscf::Rpg::Block< D > Class Template Reference

Block within a branched polymer. More...

#include <Block.h>

Inheritance diagram for Pscf::Rpg::Block< D >:
Pscf::BlockTmpl< Propagator< D > > Pscf::Edge

Public Types

using Base = BlockTmpl< Propagator<D> >
 Base class.
using PropagatorT = Propagator<D>
 Propagator type.
using FFTT = FFT<D>
 Fast Fourier Transform (FFT) type.
using WaveListT = WaveList<D>
 Wavelist type.
using FieldIoT = FieldIo<D>
 FieldIo type.
Public Types inherited from Pscf::BlockTmpl< Propagator< D > >
using PropagatorT
 Modified diffusion equation solver (propagator) type.

Public Member Functions

 Block ()
 Constructor.
 ~Block ()
 Destructor.
void associate (Mesh< D > const &mesh, FFT< D > const &fft, UnitCell< D > const &cell, WaveList< D > &wavelist)
 Create permanent associations with related objects.
void allocate (double ds, bool useBatchedFFT=true)
 Allocate memory and set contour step size for thread model.
void clearUnitCellData ()
 Clear all internal data that depends on lattice parameters.
void setLength (double newLength)
 Set or reset block length (only used in thread model).
void setKuhn (double kuhn)
 Set or reset monomer statistical segment length.
void setupSolver (RField< D > const &w)
 Set solver for this block.
void stepThread (RField< D > const &qin, RField< D > &qout)
 Compute step of integration loop, from i to i+1.
void stepBead (RField< D > const &qin, RField< D > &qout)
 Compute one step of solution of MDE for the bead model.
void stepFieldBead (RField< D > &q)
 Apply the exponential field operator for the bead model.
void stepBondBead (RField< D > const &qin, RField< D > &qout)
 Compute a bond operator for the bead model.
void stepHalfBondBead (RField< D > const &qin, RField< D > &qout)
 Compute a half-bond operator for the bead model.
void computeConcentrationThread (double prefactor)
 Compute unnormalized concentration for block by integration.
void computeConcentrationBead (double prefactor)
 Compute the concentration for this block, using the bead model.
double averageProduct (RField< D > const &q0, RField< D > const &q1)
 Compute the spatial average of the product used to compute Q.
double averageProductBead (RField< D > const &q0, RField< D > const &q1)
 Compute the spatial average of the product used to compute Q.
void computeStressThread (double prefactor)
 Compute stress contribution for this block, using thread model.
void computeStressBead (double prefactor)
 Compute stress contribution for this block, using bead model.
double stress (int n)
 Get derivative of free energy w/ respect to a unit cell parameter.
Mesh< D > const & mesh () const
 Return associated spatial Mesh by const reference.
FFT< D > const & fft () const
 Return associated FFT<D> object by const reference.
double ds () const
 Contour length step size.
int ns () const
 Number of contour length steps.
void setId (int id)
 Set the id for this block.
void setVertexIds (int vertexId0, int vertexId1)
 Set indices of associated vertices.
void setMonomerId (int monomerId)
 Set the monomer type id.
int id () const
 Get the id of this block (unique within the polymer).
int monomerId () const
 Get the monomer type id for this block.
const Pair< int > & vertexIds () const
 Get the pair of associated vertex ids.
int vertexId (int i) const
 Get the id of one associated vertex.
double length () const
 Get the length of this block, in the thread model.
int nBead () const
 Get the number of beads in this block, in the bead model.
Public Member Functions inherited from Pscf::BlockTmpl< Propagator< D > >
 BlockTmpl ()
 Constructor.
virtual ~BlockTmpl ()
 Destructor.
Propagator< D > & propagator (int directionId)
 Get a Propagator for a specified direction.
Propagator< D > const & propagator (int directionId) const
 Get a const Propagator for a specified direction.
Propagator< D >::FieldT & cField ()
 Get the associated monomer concentration field.
Propagator< D >::FieldT const & cField () const
 Get the associated const monomer concentration field.
double kuhn () const
 Get monomer statistical segment length.
Public Member Functions inherited from Pscf::Edge
 Edge ()
 Constructor.
virtual ~Edge ()
 Destructor.
template<class Archive>
void serialize (Archive &ar, unsigned int versionId)
 Serialize to/from archive.
void setId (int id)
 Set the id for this block.
void setMonomerId (int monomerId)
 Set the monomer type id.
void setVertexIds (int vertexId0, int vertexId1)
 Set indices of associated vertices.
virtual void setNBead (int nBead)
 Set the number of beads in this block (only valid for bead model).
void setPolymerType (PolymerType::Enum type)
 Set the type of the parent polymer (branched or linear).
int id () const
 Get the id of this block (unique within the polymer).
int monomerId () const
 Get the monomer type id for this block.
const Pair< int > & vertexIds () const
 Get the pair of associated vertex ids.
int vertexId (int i) const
 Get the id of one associated vertex.
double length () const
 Get the length of this block, in the thread model.
int nBead () const
 Get the number of beads in this block, in the bead model.
PolymerType::Enum polymerType () const
 Get the type of the parent polymer (branched or linear).

Detailed Description

template<int D>
class Pscf::Rpg::Block< D >

Block within a branched polymer.

Derived from BlockTmpl<Propagator<D>>. A BlockTmpl<Propagator<D>> has two Propagator<D> members and is derived from Edge.

Definition at line 54 of file rpg/solvers/Block.h.

Member Typedef Documentation

◆ Base

template<int D>
using Pscf::Rpg::Block< D >::Base = BlockTmpl< Propagator<D> >

Base class.

Definition at line 62 of file rpg/solvers/Block.h.

◆ PropagatorT

template<int D>
using Pscf::Rpg::Block< D >::PropagatorT = Propagator<D>

Propagator type.

Definition at line 65 of file rpg/solvers/Block.h.

◆ FFTT

template<int D>
using Pscf::Rpg::Block< D >::FFTT = FFT<D>

Fast Fourier Transform (FFT) type.

Definition at line 68 of file rpg/solvers/Block.h.

◆ WaveListT

template<int D>
using Pscf::Rpg::Block< D >::WaveListT = WaveList<D>

Wavelist type.

Definition at line 71 of file rpg/solvers/Block.h.

◆ FieldIoT

template<int D>
using Pscf::Rpg::Block< D >::FieldIoT = FieldIo<D>

FieldIo type.

Definition at line 74 of file rpg/solvers/Block.h.

Constructor & Destructor Documentation

◆ Block()

template<int D>
Pscf::Rpg::Block< D >::Block ( )

Constructor.

Definition at line 215 of file rpg/solvers/Block.tpp.

References Pscf::BlockTmpl< Propagator< D > >::propagator().

◆ ~Block()

template<int D>
Pscf::Rpg::Block< D >::~Block ( )

Destructor.

Definition at line 238 of file rpg/solvers/Block.tpp.

Member Function Documentation

◆ associate()

template<int D>
void Pscf::Rpg::Block< D >::associate ( Mesh< D > const & mesh,
FFT< D > const & fft,
UnitCell< D > const & cell,
WaveList< D > & wavelist )

Create permanent associations with related objects.

This function creates associations of this block with the Mesh, FFT, UnitCell, and WaveList objects by storing their addresses. It must be called before allocate().

Parameters
meshMesh<D> object - spatial discretization mesh
fftFFT<D> object - Fourier transforms
cellUnitCell<D> object - crystallographic unit cell
wavelistWaveList<D> object - properties of wavevectors

Definition at line 242 of file rpg/solvers/Block.tpp.

References fft(), mesh(), Pscf::Prdc::UnitCellBase< D >::nParameter(), and UTIL_CHECK.

◆ allocate()

template<int D>
void Pscf::Rpg::Block< D >::allocate ( double ds,
bool useBatchedFFT = true )

Allocate memory and set contour step size for thread model.

This function choses a value for the number ns of contour variable grid points for this block, sets the step size, and allocates memory for several private arrays. Spatial grid dimensions are obtained from a pointers to the associated mesh.

For the thread model, if PolymerModel::isThread() is true, the value for the number ns of contour variable grid points for this block is chosen to yield a value for the the actual step size length/(ns-1) as close as possible to the input parameter ds (the target step size), consistent with the requirements that ns be an odd integer and ns > 1. These requirements allow use of Simpson's rule for integration with respect to the contour variable s to compute monomer concentration fields and stress contributions.

For the bead model, if PolymerModel::isThread() is true, the value of ns is given by nBead + 2.

Precondition: associate() must be called before this function

Parameters
dstarget (optimal) value for contour length step size
useBatchedFFTFlag indicating whether to use batched FFTs

Definition at line 263 of file rpg/solvers/Block.tpp.

References Pscf::BlockTmpl< Propagator< D > >::cField(), Pscf::Prdc::Cpu::FFT< D >::computeKMesh(), ds(), Pscf::PolymerModel::isBead(), Pscf::PolymerModel::isThread(), length(), mesh(), nBead(), Pscf::BlockTmpl< Propagator< D > >::propagator(), and UTIL_CHECK.

◆ clearUnitCellData()

template<int D>
void Pscf::Rpg::Block< D >::clearUnitCellData ( )

Clear all internal data that depends on lattice parameters.

This method changes the internal hasExpKsq_ flag to false, so that the expKsq arrays will need to be recalculated before a step function can be called.

Definition at line 393 of file rpg/solvers/Block.tpp.

References UTIL_CHECK.

◆ setLength()

template<int D>
void Pscf::Rpg::Block< D >::setLength ( double newLength)
virtual

Set or reset block length (only used in thread model).

Precondition: PolymerModel::isThread()

Parameters
newLengthnew block length

Reimplemented from Pscf::Edge.

Definition at line 339 of file rpg/solvers/Block.tpp.

References Pscf::PolymerModel::isThread(), length(), Pscf::BlockTmpl< Propagator< D > >::propagator(), Pscf::Edge::setLength(), and UTIL_CHECK.

◆ setKuhn()

template<int D>
void Pscf::Rpg::Block< D >::setKuhn ( double kuhn)
virtual

Set or reset monomer statistical segment length.

Parameters
kuhnnew monomer statistical segment length.

Reimplemented from Pscf::BlockTmpl< Propagator< D > >.

Definition at line 383 of file rpg/solvers/Block.tpp.

References Pscf::BlockTmpl< Propagator< D > >::BlockTmpl(), Pscf::BlockTmpl< Propagator< D > >::kuhn(), and setKuhn().

Referenced by averageProductBead(), and setKuhn().

◆ setupSolver()

template<int D>
void Pscf::Rpg::Block< D >::setupSolver ( RField< D > const & w)

Set solver for this block.

This should be called once after every change in w fields, the unit cell parameters, the block length or kuhn length, before entering the loop used to solve the MDE for either propagator. This function is called by Polymer<D>::compute.

Parameters
wchemical potential field for this monomer type

Definition at line 433 of file rpg/solvers/Block.tpp.

References Pscf::Prdc::Cuda::VecOp::divSV(), Pscf::Prdc::Cuda::VecOp::expVc(), Pscf::PolymerModel::isThread(), mesh(), and UTIL_CHECK.

◆ stepThread()

template<int D>
void Pscf::Rpg::Block< D >::stepThread ( RField< D > const & qin,
RField< D > & qout )

Compute step of integration loop, from i to i+1.

This function is called internally by the Propagator::solve function within a loop over steps. It is implemented in the Block class because the same private data structures are needed for the two propagators associated with a Block.

Parameters
qinslice i of propagator q (input)
qoutslice i+1 of propagator q (output)

Definition at line 460 of file rpg/solvers/Block.tpp.

References fft(), mesh(), Pscf::Prdc::Cuda::VecOp::mulEqV(), Pscf::Prdc::Cuda::VecOp::mulEqVPair(), Pscf::Prdc::Cuda::VecOp::mulVVPair(), and UTIL_CHECK.

◆ stepBead()

template<int D>
void Pscf::Rpg::Block< D >::stepBead ( RField< D > const & qin,
RField< D > & qout )

Compute one step of solution of MDE for the bead model.

This function is called internally by the Propagator::solve function within a loop over steps. It is implemented in the Block class because the same private data structures are needed for the two propagators associated with a Block.

Parameters
qinslice i of propagator q (input)
qoutslice i+1 of propagator q (output)

Definition at line 501 of file rpg/solvers/Block.tpp.

References stepBondBead(), and stepFieldBead().

◆ stepFieldBead()

template<int D>
void Pscf::Rpg::Block< D >::stepFieldBead ( RField< D > & q)

Apply the exponential field operator for the bead model.

This function applies exp( -w(r)) in-place, where w(r) is the w-field for the monomer type of this block.

Parameters
qslice of propagator q, modified in placde

Definition at line 511 of file rpg/solvers/Block.tpp.

References Util::Array< Data >::capacity(), mesh(), Pscf::Prdc::Cuda::VecOp::mulEqV(), and UTIL_CHECK.

Referenced by stepBead().

◆ stepBondBead()

template<int D>
void Pscf::Rpg::Block< D >::stepBondBead ( RField< D > const & qin,
RField< D > & qout )

Compute a bond operator for the bead model.

This function applies exp( -|G|^2 b^2 / 6) in Fourier space. To do so, it peforms a forward FFT, multplication in k-space, and an inverse FFT.

Parameters
qinslice i of propagator q (input)
qoutslice i+1 of propagator q (output)

Definition at line 525 of file rpg/solvers/Block.tpp.

References Util::Array< Data >::capacity(), fft(), mesh(), Pscf::Prdc::Cuda::VecOp::mulEqV(), and UTIL_CHECK.

Referenced by stepBead().

◆ stepHalfBondBead()

template<int D>
void Pscf::Rpg::Block< D >::stepHalfBondBead ( RField< D > const & qin,
RField< D > & qout )

Compute a half-bond operator for the bead model.

This function applies exp( -|G|^2 b^2 / 12) in Fourier space. To do so, it peforms a forward FFT, multplication in k-space, and an inverse FFT.

Parameters
qinslice i of propagator q (input)
qoutslice i+1 of propagator q (output)

Definition at line 551 of file rpg/solvers/Block.tpp.

References Util::Array< Data >::capacity(), fft(), mesh(), Pscf::Prdc::Cuda::VecOp::mulEqV(), and UTIL_CHECK.

◆ computeConcentrationThread()

template<int D>
void Pscf::Rpg::Block< D >::computeConcentrationThread ( double prefactor)

Compute unnormalized concentration for block by integration.

The "prefactor" parameter must equal phi/(L q), where L is the total length of all blocks in the polymer species and q is the species partition function.

Upon return, grid point r of array cField() contains the integral int ds q(r,s)q^{*}(r,L-s) times the prefactor, where q(r,s) is the solution obtained from propagator(0), and q^{*} is the solution of propagator(1), and s is a contour variable that is integrated over the domain 0 < s < length(), where length() is the block length.

Parameters
prefactorconstant multiplying integral over s

Definition at line 574 of file rpg/solvers/Block.tpp.

References Pscf::BlockTmpl< Propagator< D > >::cField(), Pscf::Prdc::Cuda::VecOp::eqS(), mesh(), Pscf::Prdc::Cuda::VecOp::mulEqS(), Pscf::BlockTmpl< Propagator< D > >::propagator(), Pscf::Rpg::Propagator< D >::q(), and UTIL_CHECK.

◆ computeConcentrationBead()

template<int D>
void Pscf::Rpg::Block< D >::computeConcentrationBead ( double prefactor)

Compute the concentration for this block, using the bead model.

This function is called by Polymer::compute if a bead model is is used.

The "prefactor" parameter must equal phi/(N q), where N is the total number of beads owned by all blocks of the polymer, and q is the species partition function.

Upon return, grid point r of array cField() contains the sum sum_s ds q(r,s) q^{*}(r,N-s) exp(W(r)*ds) over beads owned by this block times the "prefactor" parameter, where q(r,s) and q^{*}(r,s) are propagators associated with different directions.

Parameters
prefactorconstant multiplying sum over beads

Definition at line 610 of file rpg/solvers/Block.tpp.

References Pscf::BlockTmpl< Propagator< D > >::cField(), Pscf::Prdc::Cuda::VecOp::eqS(), mesh(), Pscf::Prdc::Cuda::VecOp::mulEqS(), Pscf::BlockTmpl< Propagator< D > >::propagator(), Pscf::Rpg::Propagator< D >::q(), and UTIL_CHECK.

◆ averageProduct()

template<int D>
double Pscf::Rpg::Block< D >::averageProduct ( RField< D > const & q0,
RField< D > const & q1 )

Compute the spatial average of the product used to compute Q.

This function computes the spatial average of the product q0[i]*q1[i], where q0 and q1 and are complementary propagator slices, and i is a spatial mesh rank.

◆ averageProductBead()

template<int D>
double Pscf::Rpg::Block< D >::averageProductBead ( RField< D > const & q0,
RField< D > const & q1 )

Compute the spatial average of the product used to compute Q.

This computes the spatial average of the product q0[i]*q1[i]/exp(-W[i]), where q0 and q1 and are complementary propagator slices for a bead model, and i is mesh rank. This is used in the bead model for computation of Q from propagator slices associated with a bead that is owned by the propagator.

References Pscf::BlockTmpl< Propagator< D > >::BlockTmpl(), Pscf::BlockTmpl< Propagator< D > >::cField(), Pscf::Edge::id(), Pscf::BlockTmpl< Propagator< D > >::kuhn(), Pscf::Edge::length(), Pscf::Edge::monomerId(), Pscf::Edge::nBead(), Pscf::BlockTmpl< Propagator< D > >::propagator(), Pscf::Edge::setId(), setKuhn(), Pscf::Edge::setLength(), Pscf::Edge::setMonomerId(), Pscf::Edge::setVertexIds(), Pscf::Edge::vertexId(), and Pscf::Edge::vertexIds().

◆ computeStressThread()

template<int D>
void Pscf::Rpg::Block< D >::computeStressThread ( double prefactor)

Compute stress contribution for this block, using thread model.

This function is called by Polymer<D>::computeStress. The prefactor is equal to that used in computeConcentrationThread.

Parameters
prefactorconstant multiplying integral over s

Definition at line 640 of file rpg/solvers/Block.tpp.

References Util::FSArray< Data, Capacity >::append(), fft(), Pscf::Rpg::Propagator< D >::isSolved(), Pscf::PolymerModel::isThread(), Pscf::BlockTmpl< Propagator< D > >::kuhn(), mesh(), Pscf::BlockTmpl< Propagator< D > >::propagator(), Pscf::Rpg::Propagator< D >::q(), Pscf::Rpg::Propagator< D >::qAll(), and UTIL_CHECK.

◆ computeStressBead()

template<int D>
void Pscf::Rpg::Block< D >::computeStressBead ( double prefactor)

◆ stress()

template<int D>
double Pscf::Rpg::Block< D >::stress ( int n)
inline

Get derivative of free energy w/ respect to a unit cell parameter.

Parameters
nunit cell parameter index

Definition at line 479 of file rpg/solvers/Block.h.

◆ mesh()

template<int D>
Mesh< D > const & Pscf::Rpg::Block< D >::mesh ( ) const
inline

◆ fft()

template<int D>
FFT< D > const & Pscf::Rpg::Block< D >::fft ( ) const
inline

Return associated FFT<D> object by const reference.

Definition at line 492 of file rpg/solvers/Block.h.

References UTIL_ASSERT.

Referenced by associate(), computeStressBead(), computeStressThread(), stepBondBead(), stepHalfBondBead(), and stepThread().

◆ ds()

template<int D>
double Pscf::Rpg::Block< D >::ds ( ) const
inline

Contour length step size.

Definition at line 474 of file rpg/solvers/Block.h.

Referenced by allocate().

◆ ns()

template<int D>
int Pscf::Rpg::Block< D >::ns ( ) const
inline

Number of contour length steps.

Definition at line 469 of file rpg/solvers/Block.h.

◆ setId()

template<int D>
void Pscf::Edge::setId ( int id)

Set the id for this block.

Parameters
idinteger index for this block

Definition at line 89 of file Edge.cpp.

◆ setVertexIds()

template<int D>
void Pscf::Edge::setVertexIds ( int vertexId0,
int vertexId1 )

Set indices of associated vertices.

Parameters
vertexId0integer id of vertex 0
vertexId1integer id of vertex 1

Definition at line 104 of file Edge.cpp.

◆ setMonomerId()

template<int D>
void Pscf::Edge::setMonomerId ( int monomerId)

Set the monomer type id.

Parameters
monomerIdinteger id of monomer type

Definition at line 96 of file Edge.cpp.

◆ id()

template<int D>
int Pscf::Edge::id ( ) const
inline

Get the id of this block (unique within the polymer).

Definition at line 147 of file Edge.h.

◆ monomerId()

template<int D>
int Pscf::Edge::monomerId ( ) const
inline

Get the monomer type id for this block.

Definition at line 152 of file Edge.h.

◆ vertexIds()

template<int D>
const Pair< int > & Pscf::Edge::vertexIds ( ) const
inline

Get the pair of associated vertex ids.

Definition at line 157 of file Edge.h.

◆ vertexId()

template<int D>
int Pscf::Edge::vertexId ( int i) const
inline

Get the id of one associated vertex.

Parameters
iindex of vertex (0 or 1)

Definition at line 164 of file Edge.h.

◆ length()

template<int D>
double Pscf::Edge::length ( ) const
inline

Get the length of this block, in the thread model.

Precondition: PolymerModel::isThread()

Definition at line 171 of file Edge.h.

Referenced by allocate(), and setLength().

◆ nBead()

template<int D>
int Pscf::Edge::nBead ( ) const
inline

Get the number of beads in this block, in the bead model.

Precondition: PolymerModel::isBead()

Definition at line 178 of file Edge.h.

Referenced by allocate().


The documentation for this class was generated from the following files: