Gasparini_1_Potential Class Reference

Computes potential energy for Gasparini potential. More...

#include <potential.h>

Inheritance diagram for Gasparini_1_Potential:

Collaboration diagram for Gasparini_1_Potential:

Public Member Functions
	Gasparini_1_Potential (double, double, const Container *)
	Constructor.
	~Gasparini_1_Potential ()
	Destructor.
double	V (const dVec &r)
	The potential.
dVec	gradV (const dVec &)
	The gradient of the potential.
double	grad2V (const dVec &r)
	Laplacian of the potential.
Array< dVec, 1 >	initialConfig (const Container *, MTRand &, const int)
	Initial configuration corresponding to FixedAziz potential. More...
Array< double, 1 >	getExcLen ()
	get the exclusion lengths ay and az More...
Public Member Functions inherited from PotentialBase
	PotentialBase ()
	Constructor.
virtual	~PotentialBase ()
	Destructor.
virtual double	V (const dVec &, const dVec &)
	The effective potential for the pair product approximation.
virtual double	dVdlambda (const dVec &, const dVec &)
	The derivative of the effective potential with respect to lambda and tau.
virtual double	dVdtau (const dVec &, const dVec &)
void	output (const double)
	A debug method that output's the potential to a supplied separation. More...

Data Fields
const double	excZ
const double	excY
const double	V0
Data Fields inherited from PotentialBase
double	tailV
	Tail correction factor.

Additional Inherited Members
Protected Member Functions inherited from PotentialBase
double	deltaSeparation (double sep1, double sep2) const
	Return the minimum image difference for 1D separations.

Detailed Description

Computes potential energy for Gasparini potential.

This is the case of gradV = 0 for all positions in cell.

Definition at line 990 of file potential.h.

Member Function Documentation

getExcLen()

Array< double, 1 > Gasparini_1_Potential::getExcLen ( )

virtual

get the exclusion lengths ay and az

Returns the exclusion lengths ay and az.

Reimplemented from PotentialBase.

Definition at line 2033 of file potential.cpp.

                                                {
 
    Array<double, 1> excLens(2);
    excLens(0) = excY; // was ay
    excLens(1) = excZ; // was az
    
    return (excLens);
}

initialConfig()

Array< dVec, 1 > Gasparini_1_Potential::initialConfig ( const Container *  boxPtr, MTRand &  random, const int  numParticles  )

virtual

Initial configuration corresponding to FixedAziz potential.

Return initial particle configuration.

Return initial positions to exclude volume in the simulation cell.

Reimplemented from PotentialBase.

Definition at line 1874 of file potential.cpp.

                                                {
 
    /* label the lengths of the sides of the simulation cell */
    dVec lside;
    lside[0] = boxPtr->side[0];
    lside[1] = boxPtr->side[1];
    lside[2] = boxPtr->side[NDIM-1];
 
    /* calculate actual volume */
    double volTot = product(lside);
 
    /* calculate density */
    double density = (1.0*numParticles/volTot);
 
    /* calculate excluded volume */
    double volExc = lside[0]*(2.0*excY)*(2.0*excZ);
 
    /* calculate actual number of particles */
    int correctNum = int(numParticles-density*volExc);
 
    /* The particle configuration */
    Array<dVec,1> initialPos(correctNum);
 
    /* get linear size per particle  */
    double initSide = pow((1.0*correctNum/(volTot-volExc)),-1.0/(1.0*NDIM));
 
    /* For accessible volume, determine the number of 
     * initial grid boxes there are in each dimension and compute
     * their size. */
    int totNumGridBoxes1 = 1;
    int totNumGridBoxes2 = 1;
    iVec numNNGrid1;
    iVec numNNGrid2;
    dVec sizeNNGrid1;
    dVec sizeNNGrid2;
 
    /* divide space into two regions, insert particles appropriately */
    double V1 = (lside[1]-2.0*excY)*(2.0*excZ)*lside[0];
    double V2 = (lside[2]-2.0*excZ)*lside[1]*lside[0];
 
    double fracV1 = V1/(V1+V2);
 
    int numIn1 = int(correctNum*fracV1);
    int numIn2 = (correctNum-numIn1);
    
    /* grid space in volume 1 */
    /* x */
    numNNGrid1[0] = static_cast<int>(ceil((1.0*lside[0]/initSide)-EPS));
    if (numNNGrid1[0] < 1)
        numNNGrid1[0] = 1;
    sizeNNGrid1[0] = 1.0*lside[0]/(1.0*numNNGrid1[0]);
    totNumGridBoxes1 *= numNNGrid1[0];
    /* y */
    numNNGrid1[1] = static_cast<int>(ceil(((lside[1]-2.0*excY)/initSide)-EPS));
    if (numNNGrid1[1] < 1)
        numNNGrid1[1] = 1;
    sizeNNGrid1[1] = (lside[1]-2.0*excY)/(1.0*numNNGrid1[1]);
    totNumGridBoxes1 *= numNNGrid1[1];
    /* z */
    numNNGrid1[2] = static_cast<int>(ceil(((2.0*excZ)/initSide)-EPS));
    if (numNNGrid1[2] < 1)
        numNNGrid1[2] = 1;
    sizeNNGrid1[2] = (2.0*excZ)/(1.0*numNNGrid1[2]);
    totNumGridBoxes1 *= numNNGrid1[2];
    
    /* grid space in volume 2 */
    /* x */
    numNNGrid2[0] = static_cast<int>(ceil((1.0*lside[0]/initSide)-EPS));
    if (numNNGrid2[0] < 1)
        numNNGrid2[0] = 1;
    sizeNNGrid2[0] = 1.0*lside[0]/(1.0*numNNGrid2[0]);
    totNumGridBoxes2 *= numNNGrid2[0];
    /* y */
    numNNGrid2[1] = static_cast<int>(ceil((1.0*lside[1]/initSide)-EPS));
    if (numNNGrid2[1] < 1)
        numNNGrid2[1] = 1;
    sizeNNGrid2[1] = 1.0*lside[1]/(1.0*numNNGrid2[1]);
    totNumGridBoxes2 *= numNNGrid2[1];
    /* z */
    numNNGrid2[2] = static_cast<int>(ceil(((lside[2]-2.0*excZ)/initSide)-EPS));
    if (numNNGrid2[2] < 1)
        numNNGrid2[2] = 1;
    sizeNNGrid2[2] = (lside[2]-2.0*excZ)/(1.0*numNNGrid2[2]);
    totNumGridBoxes2 *= numNNGrid2[2];
    
    /* Place particles in the middle of the boxes -- volume 1 */
    PIMC_ASSERT(totNumGridBoxes1>=numIn1);
    dVec pos1;
 
    for (int n = 0; n < totNumGridBoxes1; n++) {
        iVec gridIndex1;
        /* update grid index */
        for (int i = 0; i < NDIM; i++) {
            int scale = 1;
            for (int j = i+1; j < NDIM; j++) 
                scale *= numNNGrid1[j];
            gridIndex1[i] = (n/scale) % numNNGrid1[i];
        }
        /* place particle in position vector, skipping over excluded volume */
        pos1[0] = (gridIndex1[0]+0.5)*sizeNNGrid1[0] + +0.5*lside[0] + 2.0*EPS;
        pos1[1] = (gridIndex1[1]+0.5)*sizeNNGrid1[1] + excY + 2.0*EPS;
        pos1[2] = (gridIndex1[2]+0.5)*sizeNNGrid1[2] - excZ + 2.0*EPS;
 
        if ((pos1[1]<-excY || pos1[1]>excY) || (pos1[2]<-excZ || pos1[2]>excZ))
            boxPtr->putInside(pos1);
 
        if (n < numIn1){
            initialPos(n) = pos1;
        }
        else 
            break;
    }
 
    /* Place particles in the middle of the boxes -- volume 2 */
    PIMC_ASSERT(totNumGridBoxes2>=numIn2);
    dVec pos2;
 
    for (int n = 0; n < totNumGridBoxes2; n++) {
        iVec gridIndex2;
        /* update grid index */
        for (int i = 0; i < NDIM; i++) {
            int scale = 1;
            for (int j = i+1; j < NDIM; j++) 
                scale *= numNNGrid2[j];
            gridIndex2[i] = (n/scale) % numNNGrid2[i];
        }
        /* place particles in position vectors */
        pos2[0] = (gridIndex2[0]+0.5)*sizeNNGrid2[0] + 0.5*lside[0] + 2.0*EPS;
        pos2[1] = (gridIndex2[1]+0.5)*sizeNNGrid2[1] + 0.5*lside[1] + 2.0*EPS;
        pos2[2] = (gridIndex2[2]+0.5)*sizeNNGrid2[2] + excZ + 2.0*EPS;
 
        if ((pos2[1]<-excY || pos2[1]>excY) || (pos2[2]<-excZ || pos2[2]>excZ))
            boxPtr->putInside(pos2);
        
        if (n < numIn2){
            initialPos(n+numIn1) = pos2;
        }
        else 
            break;
    }
    /* do we want to output the initial config to disk? */
    bool outToDisk = 1;
    ofstream OF;
    if (outToDisk){
        OF.open("./OUTPUT/initialConfig.dat");
        OF<<"# Cartesian Coordinates of initial Positions (X-Y-Z)"<<endl;
        OF<<"# "<<lside[0]<<"\t"<<lside[1]<<"\t"<<lside[2]<<endl;
        OF<<"# "<< excY <<"\t"<< excZ <<endl;
        for (int i=0; i< int(initialPos.size()); i++)
            OF<<initialPos(i)(0)<< "\t"<<initialPos(i)(1)<<"\t"<<initialPos(i)(2)<<endl;
        OF.close();
    }
    
    return initialPos;
}

Here is the call graph for this function:

---

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

• include/potential.h
• src/potential.cpp