[ViewVC] Diff of: group/trunk/mdtools/interface

Comparing trunk/mdtools/interface_implementation/SimSetup.cpp (file contents):
Revision 11 by mmeineke, Tue Jul 9 18:40:59 2002 UTC vs.
Revision 215 by chuckv, Thu Dec 19 21:59:51 2002 UTC

#	Line 6 \| Line 6
6		#include "parse_me.h"
7		#include "LRI.hpp"
8		#include "Integrator.hpp"
9	+	#include "simError.h"
10
11	+	#ifdef IS_MPI
12	+	#include "mpiBASS.h"
13	+	#include "mpiSimulation.hpp"
14	+	#include "bassDiag.hpp"
15	+	#endif
16	+
17		SimSetup::SimSetup(){
18		stamps = new MakeStamps();
19		globals = new Globals();
20	+
21	+	#ifdef IS_MPI
22	+	strcpy( checkPointMsg, "SimSetup creation successful" );
23	+	MPIcheckPoint();
24	+	#endif // IS_MPI
25		}
26
27		SimSetup::~SimSetup(){
#	Line 19 \| Line 31 \| *void SimSetup::parseFile( char fileName ){**
31
32		void SimSetup::parseFile( char* fileName ){
33
34	<	inFileName = fileName;
35	<	set_interface_stamps( stamps, globals );
36	<	yacc_BASS( fileName );
34	>	#ifdef IS_MPI
35	>	if( worldRank == 0 ){
36	>	#endif // is_mpi
37	>
38	>	inFileName = fileName;
39	>	set_interface_stamps( stamps, globals );
40	>
41	>	#ifdef IS_MPI
42	>	mpiEventInit();
43	>	#endif
44	>
45	>	yacc_BASS( fileName );
46	>
47	>	#ifdef IS_MPI
48	>	throwMPIEvent(NULL);
49	>	}
50	>	else receiveParse();
51	>	#endif
52	>
53		}
54
55	+	#ifdef IS_MPI
56	+	void SimSetup::receiveParse(void){
57	+
58	+	set_interface_stamps( stamps, globals );
59	+	mpiEventInit();
60	+	MPIcheckPoint();
61	+	mpiEventLoop();
62	+
63	+	}
64	+
65	+
66	+	void SimSetup::testMe(void){
67	+	bassDiag* dumpMe = new bassDiag(globals,stamps);
68	+	dumpMe->dumpStamps();
69	+	delete dumpMe;
70	+	}
71	+	#endif
72	+
73		void SimSetup::createSim( void ){
74
75		MakeStamps *the_stamps;
76		Globals* the_globals;
77	<	int i;
77	>	int i, j;
78
79		// get the stamps and globals;
80		the_stamps = stamps;
#	Line 43 \| Line 89 \| void SimSetup::createSim( void ){
89		n_components = the_globals->getNComponents();
90		strcpy( force_field, the_globals->getForceField() );
91		strcpy( ensemble, the_globals->getEnsemble() );
92	<
92	>
93		if( !strcmp( force_field, "TraPPE" ) ) the_ff = new TraPPEFF();
94		else if( !strcmp( force_field, "DipoleTest" ) ) the_ff = new DipoleTestFF();
95		else if( !strcmp( force_field, "TraPPE_Ex" ) ) the_ff = new TraPPE_ExFF();
96	+	else if( !strcmp( force_field, "LJ" ) ) the_ff = new LJ_FF();
97		else{
98	<	std::cerr<< "SimSetup Error. Unrecognized force field -> "
99	<	<< force_field << "\n";
100	<	exit(8);
98	>	sprintf( painCave.errMsg,
99	>	"SimSetup Error. Unrecognized force field -> %s\n",
100	>	force_field );
101	>	painCave.isFatal = 1;
102	>	simError();
103		}
104
105	+	#ifdef IS_MPI
106	+	strcpy( checkPointMsg, "ForceField creation successful" );
107	+	MPIcheckPoint();
108	+	#endif // is_mpi
109	+
110		// get the components and calculate the tot_nMol and indvidual n_mol
111		the_components = the_globals->getComponents();
112		components_nmol = new int[n_components];
113		comp_stamps = new MoleculeStamp*[n_components];
114	<
114	>
115		if( !the_globals->haveNMol() ){
116	<	// we don't have the total number of molecules, so we assume it is
116	>	// we don't have the total number of molecules, so we assume it is
117		// given in each component
118
119		tot_nmol = 0;
120		for( i=0; i<n_components; i++ ){
121	<
121	>
122		if( !the_components[i]->haveNMol() ){
123		// we have a problem
124	<	std::cerr << "SimSetup Error. No global NMol or component NMol"
125	<	<< " given. Cannot calculate the number of atoms.\n";
126	<	exit( 8 );
124	>	sprintf( painCave.errMsg,
125	>	"SimSetup Error. No global NMol or component NMol"
126	>	" given. Cannot calculate the number of atoms.\n" );
127	>	painCave.isFatal = 1;
128	>	simError();
129		}
130
131		tot_nmol += the_components[i]->getNMol();
#	Line 77 \| Line 133 \| void SimSetup::createSim( void ){
133		}
134		}
135		else{
136	<	std::cerr << "NOT A SUPPORTED FEATURE\n";
136	>	sprintf( painCave.errMsg,
137	>	"SimSetup error.\n"
138	>	"\tSorry, the ability to specify total"
139	>	" nMols and then give molfractions in the components\n"
140	>	"\tis not currently supported."
141	>	" Please give nMol in the components.\n" );
142	>	painCave.isFatal = 1;
143	>	simError();
144
82	–	// tot_nmol = the_globals->getNMol();
145
146	<	// //we have the total number of molecules, now we check for molfractions
147	<	// for( i=0; i<n_components; i++ ){
148	<
149	<	// if( !the_components[i]->haveMolFraction() ){
150	<
151	<	// if( !the_components[i]->haveNMol() ){
152	<	// //we have a problem
153	<	// std::cerr << "SimSetup error. Neither molFraction nor "
154	<	// << " nMol was given in component
155	<
146	>	// tot_nmol = the_globals->getNMol();
147	>
148	>	// //we have the total number of molecules, now we check for molfractions
149	>	// for( i=0; i<n_components; i++ ){
150	>
151	>	// if( !the_components[i]->haveMolFraction() ){
152	>
153	>	// if( !the_components[i]->haveNMol() ){
154	>	// //we have a problem
155	>	// std::cerr << "SimSetup error. Neither molFraction nor "
156	>	// << " nMol was given in component
157	>
158		}
159
160	+	#ifdef IS_MPI
161	+	strcpy( checkPointMsg, "Have the number of components" );
162	+	MPIcheckPoint();
163	+	#endif // is_mpi
164	+
165		// make an array of molecule stamps that match the components used.
166	+	// also extract the used stamps out into a separate linked list
167
168	+	simnfo->nComponents = n_components;
169	+	simnfo->componentsNmol = components_nmol;
170	+	simnfo->compStamps = comp_stamps;
171	+	simnfo->headStamp = new LinkedMolStamp();
172	+
173	+	char* id;
174	+	LinkedMolStamp* headStamp = simnfo->headStamp;
175	+	LinkedMolStamp* currentStamp = NULL;
176		for( i=0; i<n_components; i++ ){
177
178	<	comp_stamps[i] =
179	<	the_stamps->getMolecule( the_components[i]->getType() );
178	>	id = the_components[i]->getType();
179	>	comp_stamps[i] = NULL;
180	>
181	>	// check to make sure the component isn't already in the list
182	>
183	>	comp_stamps[i] = headStamp->match( id );
184	>	if( comp_stamps[i] == NULL ){
185	>
186	>	// extract the component from the list;
187	>
188	>	currentStamp = the_stamps->extractMolStamp( id );
189	>	if( currentStamp == NULL ){
190	>	sprintf( painCave.errMsg,
191	>	"SimSetup error: Component \"%s\" was not found in the "
192	>	"list of declared molecules\n"
193	>	id );
194	>	painCave.isFatal = 1;
195	>	simError();
196	>	}
197	>
198	>	headStamp->add( currentStamp );
199	>	comp_stamps[i] = headStamp->match( id );
200	>	}
201		}
202
203	+	#ifdef IS_MPI
204	+	strcpy( checkPointMsg, "Component stamps successfully extracted\n" );
205	+	MPIcheckPoint();
206	+	#endif // is_mpi
207
208
209	+
210	+
211		// caclulate the number of atoms, bonds, bends and torsions
212
213		tot_atoms = 0;
#	Line 111 \| Line 216 \| void SimSetup::createSim( void ){
216		tot_torsions = 0;
217		for( i=0; i<n_components; i++ ){
218
219	<	tot_atoms += components_nmol[i] * comp_stamps[i]->getNAtoms();
220	<	tot_bonds += components_nmol[i] * comp_stamps[i]->getNBonds();
221	<	tot_bends += components_nmol[i] * comp_stamps[i]->getNBends();
219	>	tot_atoms += components_nmol[i] * comp_stamps[i]->getNAtoms();
220	>	tot_bonds += components_nmol[i] * comp_stamps[i]->getNBonds();
221	>	tot_bends += components_nmol[i] * comp_stamps[i]->getNBends();
222		tot_torsions += components_nmol[i] * comp_stamps[i]->getNTorsions();
223		}
224	<
224	>
225		tot_SRI = tot_bonds + tot_bends + tot_torsions;
226	<
226	>
227		simnfo->n_atoms = tot_atoms;
228		simnfo->n_bonds = tot_bonds;
229		simnfo->n_bends = tot_bends;
230		simnfo->n_torsions = tot_torsions;
231		simnfo->n_SRI = tot_SRI;
232	+	simnfo->n_mol = tot_nmol;
233
128	–	// create the atom and short range interaction arrays
234
235	<	the_atoms = new Atom*[tot_atoms];
236	<	// the_molecules = new Molecule[tot_nmol];
235	>	#ifdef IS_MPI
236	>
237	>	// divide the molecules among processors here.
238
239	+	mpiSim = new mpiSimulation( simnfo );
240
241	<	if( tot_SRI ){
242	<	the_sris = new SRI*[tot_SRI];
243	<	the_excludes = new ex_pair[tot_SRI];
241	>	mpiSim->divideLabor();
242	>
243	>	// set up the local variables
244	>
245	>	int localMol, allMol;
246	>	int local_atoms, local_bonds, local_bends, local_torsions, local_SRI;
247	>
248	>	allMol = 0;
249	>	localMol = 0;
250	>	local_atoms = 0;
251	>	local_bonds = 0;
252	>	local_bends = 0;
253	>	local_torsions = 0;
254	>	for( i=0; i<n_components; i++ ){
255	>
256	>	for( j=0; j<components_nmol[i]; j++ ){
257	>
258	>	if( mpiSim->getMyMolStart() <= allMol &&
259	>	allMol <= mpiSim->getMyMolEnd() ){
260	>
261	>	local_atoms += comp_stamps[i]->getNAtoms();
262	>	local_bonds += comp_stamps[i]->getNBonds();
263	>	local_bends += comp_stamps[i]->getNBends();
264	>	local_torsions += comp_stamps[i]->getNTorsions();
265	>	localMol++;
266	>	}
267	>	allMol++;
268	>	}
269		}
270	+	local_SRI = local_bonds + local_bends + local_torsions;
271	+
272
273	+	simnfo->n_atoms = mpiSim->getMyNlocal();
274	+
275	+	if( local_atoms != simnfo->n_atoms ){
276	+	sprintf( painCave.errMsg,
277	+	"SimSetup error: mpiSim's localAtom (%d) and SimSetup's"
278	+	" localAtom (%d) are note equal.\n",
279	+	simnfo->n_atoms,
280	+	local_atoms );
281	+	painCave.isFatal = 1;
282	+	simError();
283	+	}
284	+
285	+	simnfo->n_bonds = local_bonds;
286	+	simnfo->n_bends = local_bends;
287	+	simnfo->n_torsions = local_torsions;
288	+	simnfo->n_SRI = local_SRI;
289	+	simnfo->n_mol = localMol;
290	+
291	+	strcpy( checkPointMsg, "Passed nlocal consistency check." );
292	+	MPIcheckPoint();
293	+
294	+
295	+	#endif // is_mpi
296	+
297	+
298	+	// create the atom and short range interaction arrays
299	+
300	+	Atom::createArrays(simnfo->n_atoms);
301	+	the_atoms = new Atom*[simnfo->n_atoms];
302	+	the_molecules = new Molecule[simnfo->n_mol];
303	+
304	+
305	+	if( simnfo->n_SRI ){
306	+	the_sris = new SRI*[simnfo->n_SRI];
307	+	the_excludes = new ex_pair[simnfo->n_SRI];
308	+	}
309	+
310		// set the arrays into the SimInfo object
311
312		simnfo->atoms = the_atoms;
#	Line 143 \| Line 314 \| void SimSetup::createSim( void ){
314		simnfo->n_exclude = tot_SRI;
315		simnfo->excludes = the_excludes;
316
317	+
318		// initialize the arrays
319	<
319	>
320		the_ff->setSimInfo( simnfo );
321	<
321	>
322		makeAtoms();
323
324		if( tot_bonds ){
#	Line 161 \| Line 333 \| void SimSetup::createSim( void ){
333		makeTorsions();
334		}
335
164	–	// makeMolecules();
336
337		// get some of the tricky things that may still be in the globals
338
339		if( simnfo->n_dipoles ){
340
341		if( !the_globals->haveRRF() ){
342	<	std::cerr << "SimSetup Error, system has dipoles, but no rRF was set.\n";
343	<	exit(8);
342	>	sprintf( painCave.errMsg,
343	>	"SimSetup Error, system has dipoles, but no rRF was set.\n");
344	>	painCave.isFatal = 1;
345	>	simError();
346		}
347		if( !the_globals->haveDielectric() ){
348	<	std::cerr << "SimSetup Error, system has dipoles, but no"
349	<	<< " dielectric was set.\n";
350	<	exit(8);
348	>	sprintf( painCave.errMsg,
349	>	"SimSetup Error, system has dipoles, but no"
350	>	" dielectric was set.\n" );
351	>	painCave.isFatal = 1;
352	>	simError();
353		}
354
355		simnfo->rRF = the_globals->getRRF();
356		simnfo->dielectric = the_globals->getDielectric();
357		}
358
359	+	#ifdef IS_MPI
360	+	strcpy( checkPointMsg, "rRf and dielectric check out" );
361	+	MPIcheckPoint();
362	+	#endif // is_mpi
363	+
364		if( the_globals->haveBox() ){
365		simnfo->box_x = the_globals->getBox();
366		simnfo->box_y = the_globals->getBox();
367		simnfo->box_z = the_globals->getBox();
368		}
369		else if( the_globals->haveDensity() ){
370	<
370	>
371		double vol;
372		vol = (double)tot_nmol / the_globals->getDensity();
373		simnfo->box_x = pow( vol, ( 1.0 / 3.0 ) );
#	Line 196 \| Line 376 \| void SimSetup::createSim( void ){
376		}
377		else{
378		if( !the_globals->haveBoxX() ){
379	<	std::cerr << "SimSetup error, no periodic BoxX size given.\n";
380	<	exit(8);
379	>	sprintf( painCave.errMsg,
380	>	"SimSetup error, no periodic BoxX size given.\n" );
381	>	painCave.isFatal = 1;
382	>	simError();
383		}
384		simnfo->box_x = the_globals->getBoxX();
385
386		if( !the_globals->haveBoxY() ){
387	<	std::cerr << "SimSetup error, no periodic BoxY size given.\n";
388	<	exit(8);
387	>	sprintf( painCave.errMsg,
388	>	"SimSetup error, no periodic BoxY size given.\n" );
389	>	painCave.isFatal = 1;
390	>	simError();
391		}
392		simnfo->box_y = the_globals->getBoxY();
393
394		if( !the_globals->haveBoxZ() ){
395	<	std::cerr << "SimSetup error, no periodic BoxZ size given.\n";
396	<	exit(8);
395	>	sprintf( painCave.errMsg,
396	>	"SimSetup error, no periodic BoxZ size given.\n" );
397	>	painCave.isFatal = 1;
398	>	simError();
399		}
400		simnfo->box_z = the_globals->getBoxZ();
401		}
216	–
217	–	if( the_globals->haveInitialConfig() ){
218	–	InitializeFromFile* fileInit;
219	–	fileInit = new InitializeFromFile( the_globals->getInitialConfig() );
220	–
221	–	fileInit->read_xyz( simnfo ); // default velocities on
402
403	<	delete fileInit;
404	<	}
405	<	else{
406	<	initFromBass();
227	<	}
403	>	#ifdef IS_MPI
404	>	strcpy( checkPointMsg, "Box size set up" );
405	>	MPIcheckPoint();
406	>	#endif // is_mpi
407
408	<	if( the_globals->haveFinalConfig() ){
409	<	strcpy( simnfo->finalName, the_globals->getFinalConfig() );
410	<	}
411	<	else{
412	<	strcpy( simnfo->finalName, inFileName );
413	<	char* endTest;
414	<	int nameLength = strlen( simnfo->finalName );
415	<	endTest = &(simnfo->finalName[nameLength - 5]);
416	<	if( !strcmp( endTest, ".bass" ) ){
417	<	strcpy( endTest, ".eor" );
408	>
409	>
410	>	if( the_globals->haveInitialConfig() ){
411	>
412	>	InitializeFromFile* fileInit;
413	>	#ifdef IS_MPI // is_mpi
414	>	if( worldRank == 0 ){
415	>	#endif //is_mpi
416	>	fileInit = new InitializeFromFile( the_globals->getInitialConfig() );
417	>	#ifdef IS_MPI
418	>	}else fileInit = new InitializeFromFile( NULL );
419	>	#endif
420	>	fileInit->read_xyz( simnfo ); // default velocities on
421	>
422	>	delete fileInit;
423	>	}
424	>	else{
425	>
426	>	#ifdef IS_MPI
427	>
428	>	// no init from bass
429	>
430	>	sprintf( painCave.errMsg,
431	>	"Cannot intialize a parallel simulation without an initial configuration file.\n" );
432	>	painCave.isFatal;
433	>	simError();
434	>
435	>	#else
436	>
437	>	initFromBass();
438	>
439	>
440	>	#endif
441	>	}
442	>
443	>	#ifdef IS_MPI
444	>	strcpy( checkPointMsg, "Successfully read in the initial configuration" );
445	>	MPIcheckPoint();
446	>	#endif // is_mpi
447	>
448	>
449	>
450	>
451	>
452	>
453	>
454	>	#ifdef IS_MPI
455	>	if( worldRank == 0 ){
456	>	#endif // is_mpi
457	>
458	>	if( the_globals->haveFinalConfig() ){
459	>	strcpy( simnfo->finalName, the_globals->getFinalConfig() );
460		}
240	–	else if( !strcmp( endTest, ".BASS" ) ){
241	–	strcpy( endTest, ".eor" );
242	–	}
461		else{
462	<	endTest = &(simnfo->finalName[nameLength - 4]);
463	<	if( !strcmp( endTest, ".bss" ) ){
462	>	strcpy( simnfo->finalName, inFileName );
463	>	char* endTest;
464	>	int nameLength = strlen( simnfo->finalName );
465	>	endTest = &(simnfo->finalName[nameLength - 5]);
466	>	if( !strcmp( endTest, ".bass" ) ){
467		strcpy( endTest, ".eor" );
468		}
469	<	else if( !strcmp( endTest, ".mdl" ) ){
469	>	else if( !strcmp( endTest, ".BASS" ) ){
470		strcpy( endTest, ".eor" );
471		}
472		else{
473	<	strcat( simnfo->finalName, ".eor" );
473	>	endTest = &(simnfo->finalName[nameLength - 4]);
474	>	if( !strcmp( endTest, ".bss" ) ){
475	>	strcpy( endTest, ".eor" );
476	>	}
477	>	else if( !strcmp( endTest, ".mdl" ) ){
478	>	strcpy( endTest, ".eor" );
479	>	}
480	>	else{
481	>	strcat( simnfo->finalName, ".eor" );
482	>	}
483		}
484		}
255	–	}
485
486	<	// make the sample and status out names
487	<
488	<	strcpy( simnfo->sampleName, inFileName );
489	<	char* endTest;
490	<	int nameLength = strlen( simnfo->sampleName );
491	<	endTest = &(simnfo->sampleName[nameLength - 5]);
492	<	if( !strcmp( endTest, ".bass" ) ){
264	<	strcpy( endTest, ".dump" );
265	<	}
266	<	else if( !strcmp( endTest, ".BASS" ) ){
267	<	strcpy( endTest, ".dump" );
268	<	}
269	<	else{
270	<	endTest = &(simnfo->sampleName[nameLength - 4]);
271	<	if( !strcmp( endTest, ".bss" ) ){
486	>	// make the sample and status out names
487	>
488	>	strcpy( simnfo->sampleName, inFileName );
489	>	char* endTest;
490	>	int nameLength = strlen( simnfo->sampleName );
491	>	endTest = &(simnfo->sampleName[nameLength - 5]);
492	>	if( !strcmp( endTest, ".bass" ) ){
493		strcpy( endTest, ".dump" );
494		}
495	<	else if( !strcmp( endTest, ".mdl" ) ){
495	>	else if( !strcmp( endTest, ".BASS" ) ){
496		strcpy( endTest, ".dump" );
497		}
498		else{
499	<	strcat( simnfo->sampleName, ".dump" );
499	>	endTest = &(simnfo->sampleName[nameLength - 4]);
500	>	if( !strcmp( endTest, ".bss" ) ){
501	>	strcpy( endTest, ".dump" );
502	>	}
503	>	else if( !strcmp( endTest, ".mdl" ) ){
504	>	strcpy( endTest, ".dump" );
505	>	}
506	>	else{
507	>	strcat( simnfo->sampleName, ".dump" );
508	>	}
509		}
510	<	}
511	<
512	<	strcpy( simnfo->statusName, inFileName );
513	<	nameLength = strlen( simnfo->statusName );
514	<	endTest = &(simnfo->statusName[nameLength - 5]);
285	<	if( !strcmp( endTest, ".bass" ) ){
286	<	strcpy( endTest, ".stat" );
287	<	}
288	<	else if( !strcmp( endTest, ".BASS" ) ){
289	<	strcpy( endTest, ".stat" );
290	<	}
291	<	else{
292	<	endTest = &(simnfo->statusName[nameLength - 4]);
293	<	if( !strcmp( endTest, ".bss" ) ){
510	>
511	>	strcpy( simnfo->statusName, inFileName );
512	>	nameLength = strlen( simnfo->statusName );
513	>	endTest = &(simnfo->statusName[nameLength - 5]);
514	>	if( !strcmp( endTest, ".bass" ) ){
515		strcpy( endTest, ".stat" );
516		}
517	<	else if( !strcmp( endTest, ".mdl" ) ){
517	>	else if( !strcmp( endTest, ".BASS" ) ){
518		strcpy( endTest, ".stat" );
519		}
520		else{
521	<	strcat( simnfo->statusName, ".stat" );
521	>	endTest = &(simnfo->statusName[nameLength - 4]);
522	>	if( !strcmp( endTest, ".bss" ) ){
523	>	strcpy( endTest, ".stat" );
524	>	}
525	>	else if( !strcmp( endTest, ".mdl" ) ){
526	>	strcpy( endTest, ".stat" );
527	>	}
528	>	else{
529	>	strcat( simnfo->statusName, ".stat" );
530	>	}
531		}
532	+
533	+	#ifdef IS_MPI
534		}
535	+	#endif // is_mpi
536
537		// set the status, sample, and themal kick times
538	<
538	>
539		if( the_globals->haveSampleTime() ){
540	<	simnfo->sampleTime = the_globals->getSampleTime();
540	>	simnfo->sampleTime = the_globals->getSampleTime();
541		simnfo->statusTime = simnfo->sampleTime;
542		simnfo->thermalTime = simnfo->sampleTime;
543		}
544		else{
545	<	simnfo->sampleTime = the_globals->getRunTime();
545	>	simnfo->sampleTime = the_globals->getRunTime();
546		simnfo->statusTime = simnfo->sampleTime;
547		simnfo->thermalTime = simnfo->sampleTime;
548		}
#	Line 325 \| Line 558 \| void SimSetup::createSim( void ){
558		// check for the temperature set flag
559
560		if( the_globals->haveTempSet() ) simnfo->setTemp = the_globals->getTempSet();
561	<
562	<
561	>
562	>
563		// make the longe range forces and the integrator
564	<
564	>
565		new AllLong( simnfo );
566	<
566	>
567		if( !strcmp( force_field, "TraPPE" ) ) new Verlet( *simnfo );
568		if( !strcmp( force_field, "DipoleTest" ) ) new Symplectic( simnfo );
569		if( !strcmp( force_field, "TraPPE_Ex" ) ) new Symplectic( simnfo );
570		}
571
572		void SimSetup::makeAtoms( void ){
573	<
573	>
574		int i, j, k, index;
575		double ux, uy, uz, uSqr, u;
576		AtomStamp* current_atom;
577		DirectionalAtom* dAtom;
578	+	int molIndex, molStart, molEnd, nMemb, lMolIndex;
579
580	+	lMolIndex = 0;
581	+	molIndex = 0;
582		index = 0;
583		for( i=0; i<n_components; i++ ){
584	<
584	>
585		for( j=0; j<components_nmol[i]; j++ ){
350	–
351	–	for( k=0; k<comp_stamps[i]->getNAtoms(); k++ ){
352	–
353	–	current_atom = comp_stamps[i]->getAtom( k );
354	–	if( current_atom->haveOrientation() ){
586
587	<	dAtom = new DirectionalAtom;
588	<	simnfo->n_oriented++;
589	<	the_atoms[index] = dAtom;
590	<
591	<	ux = current_atom->getOrntX();
592	<	uy = current_atom->getOrntY();
593	<	uz = current_atom->getOrntZ();
587	>	#ifdef IS_MPI
588	>	if( mpiSim->getMyMolStart() <= molIndex &&
589	>	molIndex <= mpiSim->getMyMolEnd() ){
590	>	#endif // is_mpi
591	>
592	>	molStart = index;
593	>	nMemb = comp_stamps[i]->getNAtoms();
594	>	for( k=0; k<comp_stamps[i]->getNAtoms(); k++ ){
595
596	<	uSqr = (ux * ux) + (uy * uy) + (uz * uz);
596	>	current_atom = comp_stamps[i]->getAtom( k );
597	>	if( current_atom->haveOrientation() ){
598	>
599	>	dAtom = new DirectionalAtom(index);
600	>	simnfo->n_oriented++;
601	>	the_atoms[index] = dAtom;
602	>
603	>	ux = current_atom->getOrntX();
604	>	uy = current_atom->getOrntY();
605	>	uz = current_atom->getOrntZ();
606	>
607	>	uSqr = (ux * ux) + (uy * uy) + (uz * uz);
608	>
609	>	u = sqrt( uSqr );
610	>	ux = ux / u;
611	>	uy = uy / u;
612	>	uz = uz / u;
613	>
614	>	dAtom->setSUx( ux );
615	>	dAtom->setSUy( uy );
616	>	dAtom->setSUz( uz );
617	>	}
618	>	else{
619	>	the_atoms[index] = new GeneralAtom(index);
620	>	}
621	>	the_atoms[index]->setType( current_atom->getType() );
622	>	the_atoms[index]->setIndex( index );
623
624	<	u = sqrt( uSqr );
625	<	ux = ux / u;
368	<	uy = uy / u;
369	<	uz = uz / u;
370	<
371	<	dAtom->setSUx( ux );
372	<	dAtom->setSUy( uy );
373	<	dAtom->setSUz( uz );
624	>	// increment the index and repeat;
625	>	index++;
626		}
375	–	else{
376	–	the_atoms[index] = new GeneralAtom;
377	–	}
378	–	the_atoms[index]->setType( current_atom->getType() );
379	–	the_atoms[index]->setIndex( index );
627
628	<	// increment the index and repeat;
629	<	index++;
628	>	molEnd = index -1;
629	>	the_molecules[lMolIndex].setNMembers( nMemb );
630	>	the_molecules[lMolIndex].setStartAtom( molStart );
631	>	the_molecules[lMolIndex].setEndAtom( molEnd );
632	>	the_molecules[lMolIndex].setStampID( i );
633	>	lMolIndex++;
634	>
635	>	#ifdef IS_MPI
636		}
637	+	#endif //is_mpi
638	+
639	+	molIndex++;
640		}
641		}
642	<
642	>
643		the_ff->initializeAtoms();
644		}
645
646		void SimSetup::makeBonds( void ){
647
648	<	int i, j, k, index, offset;
648	>	int i, j, k, index, offset, molIndex;
649		bond_pair* the_bonds;
650		BondStamp* current_bond;
651
652		the_bonds = new bond_pair[tot_bonds];
653		index = 0;
654		offset = 0;
655	+	molIndex = 0;g1
656	+
657		for( i=0; i<n_components; i++ ){
658	<
658	>
659		for( j=0; j<components_nmol[i]; j++ ){
402	–
403	–	for( k=0; k<comp_stamps[i]->getNBonds(); k++ ){
404	–
405	–	current_bond = comp_stamps[i]->getBond( k );
406	–	the_bonds[index].a = current_bond->getA() + offset;
407	–	the_bonds[index].b = current_bond->getB() + offset;
660
661	<	the_excludes[index].i = the_bonds[index].a;
662	<	the_excludes[index].j = the_bonds[index].b;
663	<
664	<	// increment the index and repeat;
665	<	index++;
666	<	}
667	<	offset += comp_stamps[i]->getNAtoms();
668	<	}
661	>	#ifdef IS_MPI
662	>	if( mpiSim->getMyMolStart() <= molIndex &&
663	>	molIndex <= mpiSim->getMyMolEnd() ){
664	>	#endif // is_mpi
665	>
666	>	for( k=0; k<comp_stamps[i]->getNBonds(); k++ ){
667	>
668	>	current_bond = comp_stamps[i]->getBond( k );
669	>	the_bonds[index].a = current_bond->getA() + offset;
670	>	the_bonds[index].b = current_bond->getB() + offset;
671	>
672	>	the_excludes[index].i = the_bonds[index].a;
673	>	the_excludes[index].j = the_bonds[index].b;
674	>
675	>	// increment the index and repeat;
676	>	index++;
677	>	}
678	>	offset += comp_stamps[i]->getNAtoms();
679	>
680	>	#ifdef IS_MPI
681	>	}
682	>	#endif is_mpi
683	>
684	>	molIndex++;
685	>	}
686		}
687	<
687	>
688		the_ff->initializeBonds( the_bonds );
689		}
690
691		void SimSetup::makeBends( void ){
692
693	<	int i, j, k, index, offset;
693	>	int i, j, k, index, offset, molIndex;
694		bend_set* the_bends;
695		BendStamp* current_bend;
696
697		the_bends = new bend_set[tot_bends];
698		index = 0;
699		offset = 0;
700	+	molIndex = 0;
701		for( i=0; i<n_components; i++ ){
702	<
702	>
703		for( j=0; j<components_nmol[i]; j++ ){
434	–
435	–	for( k=0; k<comp_stamps[i]->getNBends(); k++ ){
436	–
437	–	current_bend = comp_stamps[i]->getBend( k );
438	–	the_bends[index].a = current_bend->getA() + offset;
439	–	the_bends[index].b = current_bend->getB() + offset;
440	–	the_bends[index].c = current_bend->getC() + offset;
704
705	<	the_excludes[index + tot_bonds].i = the_bends[index].a;
706	<	the_excludes[index + tot_bonds].j = the_bends[index].c;
705	>	#ifdef IS_MPI
706	>	if( mpiSim->getMyMolStart() <= molIndex &&
707	>	molIndex <= mpiSim->getMyMolEnd() ){
708	>	#endif // is_mpi
709
710	<	// increment the index and repeat;
711	<	index++;
710	>	for( k=0; k<comp_stamps[i]->getNBends(); k++ ){
711	>
712	>	current_bend = comp_stamps[i]->getBend( k );
713	>	the_bends[index].a = current_bend->getA() + offset;
714	>	the_bends[index].b = current_bend->getB() + offset;
715	>	the_bends[index].c = current_bend->getC() + offset;
716	>
717	>	the_excludes[index + tot_bonds].i = the_bends[index].a;
718	>	the_excludes[index + tot_bonds].j = the_bends[index].c;
719	>
720	>	// increment the index and repeat;
721	>	index++;
722	>	}
723	>	offset += comp_stamps[i]->getNAtoms();
724	>
725	>	#ifdef IS_MPI
726		}
727	<	offset += comp_stamps[i]->getNAtoms();
727	>	#endif //is_mpi
728	>
729	>	molIndex++;
730		}
731		}
732	<
732	>
733		the_ff->initializeBends( the_bends );
734		}
735
736		void SimSetup::makeTorsions( void ){
737
738	<	int i, j, k, index, offset;
738	>	int i, j, k, index, offset, molIndex;
739		torsion_set* the_torsions;
740		TorsionStamp* current_torsion;
741
742		the_torsions = new torsion_set[tot_torsions];
743		index = 0;
744		offset = 0;
745	+	molIndex = 0;
746		for( i=0; i<n_components; i++ ){
747	<
747	>
748		for( j=0; j<components_nmol[i]; j++ ){
749	<
749	>
750	>	#ifdef IS_MPI
751	>	if( mpiSim->getMyMolStart() <= molIndex &&
752	>	molIndex <= mpiSim->getMyMolEnd() ){
753	>	#endif // is_mpi
754	>
755		for( k=0; k<comp_stamps[i]->getNTorsions(); k++ ){
756	<
756	>
757		current_torsion = comp_stamps[i]->getTorsion( k );
758		the_torsions[index].a = current_torsion->getA() + offset;
759		the_torsions[index].b = current_torsion->getB() + offset;
#	Line 480 \| Line 767 \| void SimSetup::makeTorsions( void ){
767		index++;
768		}
769		offset += comp_stamps[i]->getNAtoms();
770	+
771	+	#ifdef IS_MPI
772	+	}
773	+	#endif //is_mpi
774	+
775	+	molIndex++;
776		}
777		}
778	<
778	>
779		the_ff->initializeTorsions( the_torsions );
780		}
781
489	–	void SimSetup::makeMolecules( void ){
490	–
491	–	//empy for now
492	–	}
493	–
782		void SimSetup::initFromBass( void ){
783
784		int i, j, k;
#	Line 518 \| Line 806 \| void SimSetup::initFromBass( void ){
806		n_per_extra = (int)ceil( temp1 );
807
808		if( n_per_extra > 4){
809	<	std::cerr << "THere has been an error in constructing the non-complete lattice.\n";
810	<	exit(8);
809	>	sprintf( painCave.errMsg,
810	>	"SimSetup error. There has been an error in constructing"
811	>	" the non-complete lattice.\n" );
812	>	painCave.isFatal = 1;
813	>	simError();
814		}
815		}
816		else{
#	Line 528 \| Line 819 \| void SimSetup::initFromBass( void ){
819		celly = simnfo->box_y / temp3;
820		cellz = simnfo->box_z / temp3;
821		}
822	<
822	>
823		current_mol = 0;
824		current_comp_mol = 0;
825		current_comp = 0;
826		current_atom_ndx = 0;
827	<
827	>
828		for( i=0; i < n_cells ; i++ ){
829		for( j=0; j < n_cells; j++ ){
830		for( k=0; k < n_cells; k++ ){
831	<
831	>
832		makeElement( i * cellx,
833		j * celly,
834		k * cellz );
835	<
835	>
836		makeElement( i * cellx + 0.5 * cellx,
837		j * celly + 0.5 * celly,
838		k * cellz );
839	<
839	>
840		makeElement( i * cellx,
841		j * celly + 0.5 * celly,
842		k * cellz + 0.5 * cellz );
843	<
843	>
844		makeElement( i * cellx + 0.5 * cellx,
845		j * celly,
846		k * cellz + 0.5 * cellz );
#	Line 559 \| Line 850 \| void SimSetup::initFromBass( void ){
850
851		if( have_extra ){
852		done = 0;
853	<
853	>
854		int start_ndx;
855		for( i=0; i < (n_cells+1) && !done; i++ ){
856		for( j=0; j < (n_cells+1) && !done; j++ ){
857	<
857	>
858		if( i < n_cells ){
859	<
859	>
860		if( j < n_cells ){
861		start_ndx = n_cells;
862		}
863		else start_ndx = 0;
864		}
865		else start_ndx = 0;
866	<
866	>
867		for( k=start_ndx; k < (n_cells+1) && !done; k++ ){
868	<
868	>
869		makeElement( i * cellx,
870		j * celly,
871		k * cellz );
872		done = ( current_mol >= tot_nmol );
873	<
873	>
874		if( !done && n_per_extra > 1 ){
875		makeElement( i * cellx + 0.5 * cellx,
876		j * celly + 0.5 * celly,
877		k * cellz );
878		done = ( current_mol >= tot_nmol );
879		}
880	<
880	>
881		if( !done && n_per_extra > 2){
882		makeElement( i * cellx,
883		j * celly + 0.5 * celly,
884		k * cellz + 0.5 * cellz );
885		done = ( current_mol >= tot_nmol );
886		}
887	<
887	>
888		if( !done && n_per_extra > 3){
889		makeElement( i * cellx + 0.5 * cellx,
890		j * celly,
#	Line 604 \| Line 895 \| void SimSetup::initFromBass( void ){
895		}
896		}
897		}
898	<
899	<
898	>
899	>
900		for( i=0; i<simnfo->n_atoms; i++ ){
901		simnfo->atoms[i]->set_vx( 0.0 );
902		simnfo->atoms[i]->set_vy( 0.0 );
#	Line 621 \| Line 912 \| void SimSetup::makeElement( double x, double y, double
912		double rotMat[3][3];
913
914		for( k=0; k<comp_stamps[current_comp]->getNAtoms(); k++ ){
915	<
915	>
916		current_atom = comp_stamps[current_comp]->getAtom( k );
917		if( !current_atom->havePosition() ){
918	<	std::cerr << "Component " << comp_stamps[current_comp]->getID()
919	<	<< ", atom " << current_atom->getType()
920	<	<< " does not have a position specified.\n"
921	<	<< "The initialization routine is unable to give a start"
922	<	<< " position.\n";
923	<	exit(8);
918	>	sprintf( painCave.errMsg,
919	>	"SimSetup:initFromBass error.\n"
920	>	"\tComponent %s, atom %s does not have a position specified.\n"
921	>	"\tThe initialization routine is unable to give a start"
922	>	" position.\n",
923	>	comp_stamps[current_comp]->getID(),
924	>	current_atom->getType() );
925	>	painCave.isFatal = 1;
926	>	simError();
927		}
928	<
928	>
929		the_atoms[current_atom_ndx]->setX( x + current_atom->getPosX() );
930		the_atoms[current_atom_ndx]->setY( y + current_atom->getPosY() );
931		the_atoms[current_atom_ndx]->setZ( z + current_atom->getPosZ() );
932	<
932	>
933		if( the_atoms[current_atom_ndx]->isDirectional() ){
934	<
934	>
935		dAtom = (DirectionalAtom *)the_atoms[current_atom_ndx];
936	<
936	>
937		rotMat[0][0] = 1.0;
938		rotMat[0][1] = 0.0;
939		rotMat[0][2] = 0.0;
#	Line 657 \| Line 951 \| void SimSetup::makeElement( double x, double y, double
951
952		current_atom_ndx++;
953		}
954	<
954	>
955		current_mol++;
956		current_comp_mol++;
957
958		if( current_comp_mol >= components_nmol[current_comp] ){
959	<
959	>
960		current_comp_mol = 0;
961		current_comp++;
962		}

Diff Legend

-–
+Removed lines
-+
+Added lines
-<
+Changed lines
->
+Changed lines

Comparing trunk/mdtools/interface_implementation/SimSetup.cpp (file contents): Revision 11 by mmeineke, Tue Jul 9 18:40:59 2002 UTC vs. Revision 215 by chuckv, Thu Dec 19 21:59:51 2002 UTC

Diff Legend

Comparing trunk/mdtools/interface_implementation/SimSetup.cpp (file contents):
Revision 11 by mmeineke, Tue Jul 9 18:40:59 2002 UTC vs.
Revision 215 by chuckv, Thu Dec 19 21:59:51 2002 UTC