1#ifndef RPC_INTRACORRELATION_TPP
2#define RPC_INTRACORRELATION_TPP
11#include "IntraCorrelation.h"
12#include <rpc/System.h>
14#include <pscf/mesh/MeshIterator.h>
15#include <pscf/chem/Debye.h>
16#include <pscf/chem/EdgeIterator.h>
17#include <prdc/crystal/shiftToMinimum.h>
27 : systemPtr_(&system),
41 Mixture<D> const & mixture = system().mixture();
42 UnitCell<D> const & unitCell = system().domain().unitCell();
43 IntVec<D> const & dimensions = system().domain().mesh().dimensions();
44 const int nPolymer = mixture.
nPolymer();
45 const double vMonomer = mixture.
vMonomer();
48 for (
int i = 0; i < D; ++i) {
50 kMeshDimensions_[i] = dimensions[i];
51 kSize_ *= dimensions[i];
53 kMeshDimensions_[i] = dimensions[i]/2 + 1;
54 kSize_ *= (dimensions[i]/2 + 1);
69 Gmin = shiftToMinimum(G, dimensions, unitCell);
70 Gsq[iter.
rank()] = unitCell.
ksq(Gmin);
73 double phi, cPolymer, polymerLength;
74 double length, lengthA, lengthB, ksq;
75 double kuhn, kuhnA, kuhnB, eA, eB;
76 int monomerId, monomerIdA, monomerIdB, rank;
79 for (
int i = 0; i < correlations.
capacity(); ++i){
80 correlations[i] = 0.0;
84 for (
int i = 0; i < nPolymer; i++){
88 const int nBlock = polymer.
nBlock();
93 for (
int j = 0; j < nBlock; j++) {
95 polymerLength += length;
99 cPolymer = phi/(polymerLength*vMonomer);
102 for (
int j = 0; j < nBlock; j++) {
112 correlations[rank] += cPolymer *
Debye::d(ksq, length, kuhn);
122 for (
int ia = 1; ia < nBlock; ++ia) {
126 lengthA = blockA.
length();
131 for (
int ib = 0; ib < ia; ++ib) {
135 lengthB = blockB.
length();
142 EdgeItr.
begin(ia, ib);
143 while (EdgeItr.
notEnd()) {
145 if (edgeId != ia && edgeId != ib){
147 double lengthK = blockK.
length();
150 d += lengthK * kuhnK * kuhnK;
157 double prefactor = 2.0*cPolymer;
161 x = std::exp( -d * ksq / 6.0);
164 correlations[rank] += prefactor * x * eA * eB;
Edge iterator for graph associated with a polymer.
int currentEdgeId() const
Get index of the current edge.
bool notEnd() const
Return true iff currentId != targetId.
void begin(int sourceId, int targetId)
Initialize iterator.
An IntVec<D, T> is a D-component vector of elements of integer type T.
Iterator over points in a Mesh<D>.
int rank() const
Get the rank of current element.
void begin()
Set iterator to the first point in the mesh.
bool atEnd() const
Is this the end (i.e., one past the last point)?
void setDimensions(const IntVec< D > &dimensions)
Set the grid dimensions in all directions.
IntVec< D > position() const
Get current position in the grid, as integer vector.
Polymer & polymer(int id)
Get a polymer object.
double vMonomer() const
Get monomer reference volume (set to 1.0 by default).
Monomer const & monomer(int id) const
Get a Monomer type descriptor (const reference).
int nPolymer() const
Get number of polymer species.
double kuhn() const
Statistical segment length (random walk step size).
Block & block(int id)
Get a specified Block (block solver and descriptor)
double phi() const
Get the overall volume fraction for this species.
int nBlock() const
Number of blocks.
Field of real double precision values on an FFT mesh.
virtual double ksq(IntVec< D > const &k) const
Compute square magnitude of reciprocal lattice vector.
Base template for UnitCell<D> classes, D=1, 2 or 3.
Block within a branched polymer.
int monomerId() const
Get the monomer type id.
double length() const
Get the length (number of monomers) in this block.
void computeIntraCorrelations(RField< D > &correlations)
Compute and modify intramolecular correlations.
IntraCorrelation(System< D > &system)
Constructor.
~IntraCorrelation()
Destructor.
Solver for a mixture of polymers and solvents.
Descriptor and solver for one polymer species.
Main class for SCFT or PS-FTS simulation of one system.
int capacity() const
Return allocated size.
Dynamically allocatable contiguous array template.
void allocate(int capacity)
Allocate the underlying C array.
void setClassName(const char *className)
Set class name string.
File containing preprocessor macros for error handling.
#define UTIL_CHECK(condition)
Assertion macro suitable for serial or parallel production code.
double d(double ksq, double length, double kuhn)
Compute and return a homopolymer correlation function.
double e(double ksq, double length, double kuhn)
Compute and return one-sided correlation factor for one block.
PSCF package top-level namespace.
Utility classes for scientific computation.