OOPSE/libmdtools/SimInfo.cpp

#include <cstdlib>
#include <cstring>


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

#include "fortranWrappers.hpp"

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* 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();

  isError = 0;

//   fInfo;
//   n_atoms; 
//   identArray;
//   n_exclude;
//   excludes;
//   nGlobalExcludes;
//   globalExcludes;
//   isError;

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