ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/branches/new_design/OOPSE-3.0/src/brains/SimCreator.cpp
(Generate patch)

Comparing branches/new_design/OOPSE-3.0/src/brains/SimCreator.cpp (file contents):
Revision 1712 by tim, Thu Nov 4 20:55:01 2004 UTC vs.
Revision 1725 by tim, Wed Nov 10 22:01:06 2004 UTC

# Line 22 | Line 22
22   * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
23   *
24   */
25 –
26 – /**
27 –  * @file SimCreator.cpp
28 –  * @author tlin
29 –  * @date 11/03/2004
30 –  * @time 13:51am
31 –  * @version 1.0
32 –  */
25  
26 < #include "brains/SimCreator.hpp"  
26 > /**
27 > * @file SimCreator.cpp
28 > * @author tlin
29 > * @date 11/03/2004
30 > * @time 13:51am
31 > * @version 1.0
32 > */
33 >
34 > #include "brains/SimCreator.hpp"
35   namespace oopse {
36  
37   void SimSetup::parseFile(char* mdfile,  MakeStamps* stamps, Globals* globals){
# Line 63 | Line 63 | SimModel* SimCreator::createSim(const std::string& mdf
63   }
64  
65  
66 < SimModel* SimCreator::createSim(const std::string& mdfile) {
67 <    MakeStamps* stamps;
66 > SimInfo* SimCreator::createSim(const std::string& mdfile) {
67 >    MakeStamps* stamps = new MakeStamps();
68      
69      Globals* globals = new Globals();
70  
# Line 72 | Line 72 | SimModel* SimCreator::createSim(const std::string& mdf
72      parseFile(mdfile, stamps, globals);
73  
74      //create the force field
75 <    ForceFiled* ff = ForceFiledFactory::getInstance()->createForceField(globals->getForceField());
75 >    ForceFiled* ff = ForceFieldFactory::getInstance()->createObject(globals->getForceField());
76      
77 <    //create SimModel
78 <    model = new SimInfo();
77 >    //create SimInfo
78 >    SimInfo* info = new SimInfo();
79 >    info->setGlobals(globals);
80 >    info->setForceField(ff);
81 >
82 >    //extract the molecule stamps and add them into SimInfo
83 >    compList(stamps, info);
84      
85 <    //gather parameters
86 <    gatherParameters(model);
85 >    //gather parameters (SimCreator only retrieves the parameters which will be used to create
86 >    // the simulation)
87 >    gatherParameters(info);
88      
83 –
84 –
89      //divide the molecules and determine the global index of molecules
90 <
90 >    
91      //create the molecules
92 +    createMolecules(info);
93  
94 <    //create atoms, bonds, bend, torsions, rigidbodies
95 <
94 >    //set the global index of atoms, rigidbodies and cutoffgroups (only need to be set once, the
95 >    //global index will never change again). Local indices of atoms and rigidbodies are already set by
96 >    //MoleculeCreator class which actually delegates the responsibility to LocalIndexManager.
97 >    setGlobalIndices(info);
98 >    
99      //allocate memory for DataStorage(circular reference, need to break it)
100 <    model->setSnapshotManager(new SimSnapshotManager(model);
100 >    info->setSnapshotManager(new SimSnapshotManager(info);
101      
102 <    //load initial coordinates
103 <    DumpReader reader();
104 <    reader.readFrame(model);
102 >    //load initial coordinates, some extra information are pushed into SimInfo's property map ( such as
103 >    //eta, chi for NPT integrator)
104 >    DumpReader reader(info->getInitFilename());
105 >    int nframes = reader->getNframes();
106 >
107 >    if (nframes > 0) {
108 >        reader.readFrame(info, nframes - 1);
109 >    } else {
110 >        //invalid initial coordinate file
111 >        
112 >    }
113 >
114      
115      //initialize fortran    
116 <    return model;
116 >    
117 >    return info;
118   }    
119  
120 < void SimCreator::gatherParameters(SimModel* model) {
121 <    model->addProperty(new StringGenericData("Ensemble", globals->getForceFiled()));
122 <    model->addProperty(new DoubleGenericData("dt"), globals->getDt());
120 > void SimCreator::gatherParameters(SimInfo* info) {
121 >    //model->addProperty(new StringGenericData("Ensemble", globals->getForceFiled()));
122 >    //model->addProperty(new DoubleGenericData("dt"), globals->getDt());
123 >
124 >    //setup seed for random number generator
125 >    int seedValue;
126 >    Globals* globals = info->getGlobals();
127 >    
128 >    if (globals->haveSeed()){
129 >        seedValue = globals->getSeed();
130 >
131 >        if(seedValue / 1000000000 == 0){
132 >            sprintf(painCave.errMsg,
133 >            "Seed for sprng library should contain at least 9 digits\n"
134 >            "OOPSE will generate a seed for user\n");
135 >            painCave.isFatal = 0;
136 >            simError();
137 >
138 >    //using seed generated by system instead of invalid seed set by user
139 > #ifndef IS_MPI
140 >            seedValue = make_sprng_seed();
141 > #else
142 >            if (worldRank == 0){
143 >                seedValue = make_sprng_seed();
144 >            }
145 >            MPI_Bcast(&seedValue, 1, MPI_INT, 0, MPI_COMM_WORLD);  
146 > #endif      
147 >        } //end if (seedValue /1000000000 == 0)
148 >    } else{
149 >    
150 > #ifndef IS_MPI
151 >        seedValue = make_sprng_seed();
152 > #else
153 >        if (worldRank == 0){
154 >            seedValue = make_sprng_seed();
155 >        }
156 >        MPI_Bcast(&seedValue, 1, MPI_INT, 0, MPI_COMM_WORLD);  
157 > #endif
158 >    }//end of globals->haveSeed()
159 >    info->setSeed(seedValue);
160 >
161 >    //
162 >    std::string prefix;
163 > #ifdef IS_MPI
164 >    if (worldRank == 0){
165 > #endif // is_mpi
166 >      
167 >        if(globals->haveFinalConfig()) {
168 >            prefix = getPrefix(globals->getFinalConfig());  
169 >        } else {
170 >            prefix = getPrefix(mdfile_);
171 >        }
172 >        
173 >        info->setFinalConfFileName(prefix + ".eor");
174 >        info->setDumpFileName(prefix + ".dump");
175 >        info->setStatFileName(prefix + ".stat");
176 >
177 > #ifdef IS_MPI
178 >    }
179 > #endif    
180 >
181 >
182 >
183   }
184  
185   #ifdef IS_MPI
186 < void SimCreator::mpiMolDivide() {
186 > void SimCreator::mpiMolDivide(){
187 >
188 >    mpiSim = new mpiSimulation(info);
189  
190 +    mpiSim->divideLabor();
191  
192   }
193 + #endif
194  
195 + Molecule* SimCreator::createMolecules(SimInfo* info) {
196 +    MoleculeCreator molCreator;
197 +    int stampId;
198 +
199 +    
200 +    for (int i = 0; i < info->getNGlobalMolecules(); i++){
201 +        if (info->getMolToProc(i) == worldRank){
202 +
203 +            stampId = info->getMoleculeStampId(i);
204 +            Molecule* mol = molCreator.createMolecule(info->getForceField(),
205 +                                                                                info->getMoleculeStamp(stampId), stampId, i);
206 +                        
207 +            info->addMolecule(mol);
208 +        }
209 +    }        
210 +      
211 + }
212 +
213 +
214 + void SimSetup::compList(MakeStamps* stamps,SimInfo* info) {
215 +    int i;
216 +    char* id;
217 +    MoleculeStamp* currentStamp;
218 +    Component* the_components = info->getGlobals()->getComponents();
219 +    int n_components = info->getGlobals()->getNComponents();
220 +
221 +    if (!globals->haveNMol()){
222 +        // we don't have the total number of molecules, so we assume it is
223 +        // given in each component
224 +
225 +        for (i = 0; i < n_components; i++){
226 +            if (!the_components[i]->haveNMol()){
227 +                // we have a problem
228 +                sprintf(painCave.errMsg,
229 +                    "SimSetup Error. No global NMol or component NMol given.\n"
230 +                    "\tCannot calculate the number of atoms.\n");
231 +                painCave.isFatal = 1;
232 +                simError();
233 +            }
234 +
235 +            id = the_components[i]->getType();
236 +            currentStamp = stamps->extractMolStamp(id);
237 +            if (currentStamp == NULL){
238 +                sprintf(painCave.errMsg,
239 +                    "SimSetup error: Component \"%s\" was not found in the "
240 +                    "list of declared molecules\n",
241 +                    id);
242 +                painCave.isFatal = 1;
243 +                simError();
244 +            }      
245 +
246 +            info->addMoleculeStamp(currentStamp, the_components[i]->getNMol);
247 +          
248 +        } //end for (i = 0; i < n_components; i++)
249 +        
250 +    } else{
251 +        sprintf(painCave.errMsg,
252 +                "SimSetup error.\n"
253 +                "\tSorry, the ability to specify total"
254 +                " nMols and then give molfractions in the components\n"
255 +                "\tis not currently supported."
256 +                " Please give nMol in the components.\n");
257 +        painCave.isFatal = 1;
258 +        simError();
259 +    }
260 +  
261 + #ifdef IS_MPI
262 +  strcpy(checkPointMsg, "Component stamps successfully extracted\n");
263 +  MPIcheckPoint();
264 + #endif // is_mpi
265 + }
266 +
267 + void SimCreator::setGlobalIndices(SimInfo* info) {
268 +    typename SimInfo::MoleculeIterator mi;
269 +    typename Molecule::AtomIterator ai;
270 +    typename Molecule::RigidBodyIterator ri;
271 +    typename Molecule::CutoffGroupIterator ci;
272 +    Molecule* mol;
273 +    Atom* atom;
274 +    RigidBody* rb;
275 +    CutoffGroup* cg;
276 +    int beginAtomIndex;
277 +    int beginRigidBodyIndex;
278 +    int beginCutoffGroupIndex;
279 +
280 + #ifdef IS_MPI
281 +    beginAtomIndex = 0;
282 +    beginRigidBodyIndex = 0;
283 +    beginCutoffGroupIndex = 0;
284 + #else
285 +    int nproc;
286 +    int myNode;
287 +
288 +    myNode = worldRank;
289 +    MPI_Comm_size(MPI_COMM_WORLD, &nproc);
290 +    
291 +    std::vector<int> tmpAtomsInProc(nproc, 0);
292 +    std::vector<int> tmpRigidBodiesInProc(nproc, 0);
293 +    std::vector<int> tmpCutoffGroupsInProc(nproc, 0);    
294 +    std::vector<int> NumAtomsInProc(nproc, 0);
295 +    std::vector<int> NumRigidBodiesInProc(nproc, 0);
296 +    std::vector<int> NumCutoffGroupsInProc(nproc, 0);
297 +
298 +    tmpAtomInProc[myNode] = info->getNAtoms();
299 +    tmpRigidBodiesInProc[myNode] = info->getNRigidBodiess();
300 +    tmpCutoffGroupsInProc[myNode] = info->getNCutoffGroupss();
301 +    
302 +    //do MPI_ALLREDUCE to exchange the total number of atoms, rigidbodies and cutoff groups
303 +    MPI_Allreduce(&tmpAtomsInProc[0], &NumAtomsInProc[0], nproc, MPI_INT,MPI_SUM, MPI_COMM_WORLD);
304 +    MPI_Allreduce(&tmpRigidBodiesInProc[0], &NumRigidBodiesInProc[0], nproc, MPI_INT,MPI_SUM, MPI_COMM_WORLD);
305 +    MPI_Allreduce(&tmpCutoffGroupsInProc[0], &NumCutoffGroupsInProc[0], nproc, MPI_INT,MPI_SUM, MPI_COMM_WORLD);
306 +
307 +    beginAtomIndex = 0;
308 +    beginRigidBodyIndex = 0;
309 +    beginCutoffGroupIndex = 0;
310 +
311 +    for (int i = 0; i < nproc; i++) {
312 +        beginAtomIndex += NumAtomsInProc[i];
313 +        beginRigidBodyIndex += NumRigidBodiesInProc[i];
314 +        beginCutoffGroupIndex += NumCutoffGroupsInProc[i];        
315 +    }
316 +    
317   #endif
318 +    
319 +    for (mol = info->beginMolecule(mi); mol != NULL; mol  = info->nextMolecule(mi)) {
320 +
321 +        //local index(index in DataStorge) of atom is important
322 +        for (atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) {
323 +            atom->setGlobalIndex(beginAtomIndex++);
324 +        }
325 +
326 +        for (rb = mol->beginRigidBody(ri); rb != NULL; rb = mol->nextRigidBody(ri)) {
327 +            rb->setGlobalIndex(beginRigidBodyIndex++);
328 +        }
329 +        
330 +        //local index of cutoff group is trivial, it only depends on the order of travesing
331 +        for (cg = mol->beginCutoffGroup(ci); cg != NULL; cg = mol->nextCutoffGroup(ci)) {
332 +            cg->setGlobalIndex(beginCutoffGroupIndex++);
333 +        }        
334 +        
335 +    }
336 + }
337 +
338   } //end namespace oopse

Diff Legend

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