OOPSE/libmdtools/SimInfo.cpp

#include <cstdlib>
#include <cstring>


#include "SimInfo.hpp"
#define __C
#include "fSimulation.h"
#include "simError.h"

#include "fortranWrappers.hpp"

#ifdef IS_MPI
#include "mpiSimulation.hpp"
#endif

SimInfo* currentInfo;

SimInfo::SimInfo(){
  excludes = NULL;
  n_constraints = 0;
  n_oriented = 0;
  n_dipoles = 0;
  ndf = 0;
  ndfRaw = 0;
  the_integrator = NULL;
  setTemp = 0;
  thermalTime = 0.0;
  rCut = 0.0;

  usePBC = 0;
  useLJ = 0; 
  useSticky = 0;
  useDipole = 0;
  useReactionField = 0;
  useGB = 0;
  useEAM = 0;

  wrapMeSimInfo( this );
}

void SimInfo::setBox(double newBox[3]) {
  double smallestBox, maxCutoff;
  int status;
  box_x = newBox[0];
  box_y = newBox[1];
  box_z = newBox[2];
  setFortranBoxSize(newBox);

  smallestBox = box_x;
  if (box_y < smallestBox) smallestBox = box_y;
  if (box_z < smallestBox) smallestBox = box_z;

  maxCutoff = smallestBox / 2.0;

  if (rList > maxCutoff) {
    sprintf( painCave.errMsg,
             "New Box size is forcing neighborlist radius down to %lf\n",
             maxCutoff );
    painCave.isFatal = 0;
    simError();

    rList = maxCutoff;

    sprintf( painCave.errMsg,
             "New Box size is forcing cutoff radius down to %lf\n",
             maxCutoff - 1.0 );
    painCave.isFatal = 0;
    simError();

    rCut = rList - 1.0;

    // list radius changed so we have to refresh the simulation structure.
    refreshSim();
  }

  if (rCut > maxCutoff) {
    sprintf( painCave.errMsg,
             "New Box size is forcing cutoff radius down to %lf\n",
             maxCutoff );
    painCave.isFatal = 0;
    simError();

    status = 0;
    LJ_new_rcut(&rCut, &status);
    if (status != 0) {
      sprintf( painCave.errMsg,
               "Error in recomputing LJ shifts based on new rcut\n");
      painCave.isFatal = 1;
      simError();
    }
  }
}

void SimInfo::getBox(double theBox[3]) {
  theBox[0] = box_x;
  theBox[1] = box_y;
  theBox[2] = box_z;
}
 
int SimInfo::getNDF(){
  int ndf_local, ndf;
  
  ndf_local = 3 * n_atoms + 3 * n_oriented - n_constraints;

#ifdef IS_MPI
  MPI_Allreduce(&ndf_local,&ndf,1,MPI_INT,MPI_SUM, MPI_COMM_WORLD);
#else
  ndf = ndf_local;
#endif

  ndf = ndf - 3;

  return ndf;
}

int SimInfo::getNDFraw() {
  int ndfRaw_local, ndfRaw;

  // Raw degrees of freedom that we have to set
  ndfRaw_local = 3 * n_atoms + 3 * n_oriented;
  
#ifdef IS_MPI
  MPI_Allreduce(&ndfRaw_local,&ndfRaw,1,MPI_INT,MPI_SUM, MPI_COMM_WORLD);
#else
  ndfRaw = ndfRaw_local;
#endif

  return ndfRaw;
}
 
void SimInfo::refreshSim(){

  simtype fInfo;
  int isError;
  int n_global;
  int* excl;
  
  fInfo.rrf = 0.0;
  fInfo.rt = 0.0;
  fInfo.dielect = 0.0;

  fInfo.box[0] = box_x;
  fInfo.box[1] = box_y;
  fInfo.box[2] = box_z;

  fInfo.rlist = rList;
  fInfo.rcut = rCut;

  if( useDipole ){
    fInfo.rrf = ecr;
    fInfo.rt = ecr - est;
    if( useReactionField )fInfo.dielect = dielectric;
  }

  fInfo.SIM_uses_PBC = usePBC;
  //fInfo.SIM_uses_LJ = 0;
  fInfo.SIM_uses_LJ = useLJ;
  fInfo.SIM_uses_sticky = useSticky;
  //fInfo.SIM_uses_sticky = 0;
  fInfo.SIM_uses_dipoles = useDipole;
  //fInfo.SIM_uses_dipoles = 0;
  //fInfo.SIM_uses_RF = useReactionField;
  fInfo.SIM_uses_RF = 0;
  fInfo.SIM_uses_GB = useGB;
  fInfo.SIM_uses_EAM = useEAM;

  excl = Exclude::getArray();

#ifdef IS_MPI
  n_global = mpiSim->getTotAtoms();
#else
  n_global = n_atoms;
#endif

  isError = 0;

  setFsimulation( &fInfo, &n_global, &n_atoms, identArray, &n_exclude, excl, 
                  &nGlobalExcludes, globalExcludes, molMembershipArray, 
                  &isError );

  if( isError ){

    sprintf( painCave.errMsg,
             "There was an error setting the simulation information in fortran.\n" );
    painCave.isFatal = 1;
    simError();
  }

#ifdef IS_MPI
  sprintf( checkPointMsg,
           "succesfully sent the simulation information to fortran.\n");
  MPIcheckPoint();
#endif // is_mpi

  this->ndf = this->getNDF();
  this->ndfRaw = this->getNDFraw();

}

Revision:	490
Committed:	Fri Apr 11 15:16:59 2003 UTC (21 years, 3 months ago) by gezelter
File size:	3955 byte(s)
Log Message:	Bug fix in progress for NPT
#	User	Rev	Content
1	mmeineke	377	#include <cstdlib>
2			#include <cstring>
3
4
5			#include "SimInfo.hpp"
6			#define __C
7			#include "fSimulation.h"
8			#include "simError.h"
9
10			#include "fortranWrappers.hpp"
11
12	gezelter	490	#ifdef IS_MPI
13			#include "mpiSimulation.hpp"
14			#endif
15
16	mmeineke	377	SimInfo* currentInfo;
17
18			SimInfo::SimInfo(){
19			excludes = NULL;
20			n_constraints = 0;
21			n_oriented = 0;
22			n_dipoles = 0;
23	gezelter	458	ndf = 0;
24			ndfRaw = 0;
25	mmeineke	377	the_integrator = NULL;
26			setTemp = 0;
27			thermalTime = 0.0;
28	mmeineke	420	rCut = 0.0;
29	mmeineke	377
30			usePBC = 0;
31			useLJ = 0;
32			useSticky = 0;
33			useDipole = 0;
34			useReactionField = 0;
35			useGB = 0;
36			useEAM = 0;
37
38	gezelter	457	wrapMeSimInfo( this );
39			}
40	mmeineke	377
41	gezelter	457	void SimInfo::setBox(double newBox[3]) {
42	gezelter	463	double smallestBox, maxCutoff;
43			int status;
44	gezelter	457	box_x = newBox[0];
45			box_y = newBox[1];
46			box_z = newBox[2];
47			setFortranBoxSize(newBox);
48	gezelter	463
49			smallestBox = box_x;
50			if (box_y < smallestBox) smallestBox = box_y;
51			if (box_z < smallestBox) smallestBox = box_z;
52
53			maxCutoff = smallestBox / 2.0;
54
55			if (rList > maxCutoff) {
56			sprintf( painCave.errMsg,
57			"New Box size is forcing neighborlist radius down to %lf\n",
58			maxCutoff );
59			painCave.isFatal = 0;
60			simError();
61
62			rList = maxCutoff;
63
64			sprintf( painCave.errMsg,
65			"New Box size is forcing cutoff radius down to %lf\n",
66			maxCutoff - 1.0 );
67			painCave.isFatal = 0;
68			simError();
69
70			rCut = rList - 1.0;
71
72			// list radius changed so we have to refresh the simulation structure.
73			refreshSim();
74			}
75
76			if (rCut > maxCutoff) {
77			sprintf( painCave.errMsg,
78			"New Box size is forcing cutoff radius down to %lf\n",
79			maxCutoff );
80			painCave.isFatal = 0;
81			simError();
82
83			status = 0;
84			LJ_new_rcut(&rCut, &status);
85			if (status != 0) {
86			sprintf( painCave.errMsg,
87			"Error in recomputing LJ shifts based on new rcut\n");
88			painCave.isFatal = 1;
89			simError();
90			}
91			}
92	gezelter	457	}
93	mmeineke	377
94	gezelter	457	void SimInfo::getBox(double theBox[3]) {
95			theBox[0] = box_x;
96			theBox[1] = box_y;
97			theBox[2] = box_z;
98	mmeineke	377	}
99	gezelter	458
100			int SimInfo::getNDF(){
101			int ndf_local, ndf;
102	gezelter	457
103	gezelter	458	ndf_local = 3 * n_atoms + 3 * n_oriented - n_constraints;
104
105			#ifdef IS_MPI
106			MPI_Allreduce(&ndf_local,&ndf,1,MPI_INT,MPI_SUM, MPI_COMM_WORLD);
107			#else
108			ndf = ndf_local;
109			#endif
110
111			ndf = ndf - 3;
112
113			return ndf;
114			}
115
116			int SimInfo::getNDFraw() {
117			int ndfRaw_local, ndfRaw;
118
119			// Raw degrees of freedom that we have to set
120			ndfRaw_local = 3 * n_atoms + 3 * n_oriented;
121
122			#ifdef IS_MPI
123			MPI_Allreduce(&ndfRaw_local,&ndfRaw,1,MPI_INT,MPI_SUM, MPI_COMM_WORLD);
124			#else
125			ndfRaw = ndfRaw_local;
126			#endif
127
128			return ndfRaw;
129			}
130
131	mmeineke	377	void SimInfo::refreshSim(){
132
133			simtype fInfo;
134			int isError;
135	gezelter	490	int n_global;
136	mmeineke	424	int* excl;
137	mmeineke	469
138			fInfo.rrf = 0.0;
139			fInfo.rt = 0.0;
140			fInfo.dielect = 0.0;
141	mmeineke	377
142			fInfo.box[0] = box_x;
143			fInfo.box[1] = box_y;
144			fInfo.box[2] = box_z;
145
146			fInfo.rlist = rList;
147			fInfo.rcut = rCut;
148
149	mmeineke	469	if( useDipole ){
150			fInfo.rrf = ecr;
151			fInfo.rt = ecr - est;
152			if( useReactionField )fInfo.dielect = dielectric;
153			}
154
155	mmeineke	377	fInfo.SIM_uses_PBC = usePBC;
156	mmeineke	443	//fInfo.SIM_uses_LJ = 0;
157	chuckv	439	fInfo.SIM_uses_LJ = useLJ;
158	mmeineke	443	fInfo.SIM_uses_sticky = useSticky;
159			//fInfo.SIM_uses_sticky = 0;
160	chuckv	482	fInfo.SIM_uses_dipoles = useDipole;
161			//fInfo.SIM_uses_dipoles = 0;
162	mmeineke	443	//fInfo.SIM_uses_RF = useReactionField;
163			fInfo.SIM_uses_RF = 0;
164	mmeineke	377	fInfo.SIM_uses_GB = useGB;
165			fInfo.SIM_uses_EAM = useEAM;
166
167	mmeineke	424	excl = Exclude::getArray();
168	mmeineke	377
169	gezelter	490	#ifdef IS_MPI
170			n_global = mpiSim->getTotAtoms();
171			#else
172			n_global = n_atoms;
173			#endif
174
175	mmeineke	377	isError = 0;
176
177	gezelter	490	setFsimulation( &fInfo, &n_global, &n_atoms, identArray, &n_exclude, excl,
178	gezelter	483	&nGlobalExcludes, globalExcludes, molMembershipArray,
179			&isError );
180	mmeineke	377
181			if( isError ){
182
183			sprintf( painCave.errMsg,
184			"There was an error setting the simulation information in fortran.\n" );
185			painCave.isFatal = 1;
186			simError();
187			}
188
189			#ifdef IS_MPI
190			sprintf( checkPointMsg,
191			"succesfully sent the simulation information to fortran.\n");
192			MPIcheckPoint();
193			#endif // is_mpi
194	gezelter	458
195	gezelter	474	this->ndf = this->getNDF();
196			this->ndfRaw = this->getNDFraw();
197	gezelter	458
198	mmeineke	377	}
199