--- trunk/OOPSE/libmdtools/mpiSimulation.cpp 2003/09/25 19:27:15 787 +++ trunk/OOPSE/libmdtools/mpiSimulation.cpp 2004/06/01 21:45:22 1217 @@ -1,8 +1,8 @@ #ifdef IS_MPI #include -#include -#include -#include +#include +#include +#include #include #include "mpiSimulation.hpp" @@ -15,14 +15,15 @@ mpiSimulation::mpiSimulation(SimInfo* the_entryPlug) mpiSimulation::mpiSimulation(SimInfo* the_entryPlug) { entryPlug = the_entryPlug; - mpiPlug = new mpiSimData; + parallelData = new mpiSimData; - MPI_Comm_size(MPI_COMM_WORLD, &(mpiPlug->numberProcessors) ); - mpiPlug->myNode = worldRank; + MPI_Comm_size(MPI_COMM_WORLD, &(parallelData->nProcessors) ); + parallelData->myNode = worldRank; MolToProcMap = new int[entryPlug->n_mol]; MolComponentType = new int[entryPlug->n_mol]; AtomToProcMap = new int[entryPlug->n_atoms]; + GroupToProcMap = new int[entryPlug->ngroup]; mpiSim = this; wrapMeSimParallel( this ); @@ -34,70 +35,81 @@ mpiSimulation::~mpiSimulation(){ delete[] MolToProcMap; delete[] MolComponentType; delete[] AtomToProcMap; + delete[] GroupToProcMap; - delete mpiPlug; + delete parallelData; // perhaps we should let fortran know the party is over. } -int* mpiSimulation::divideLabor( void ){ +void mpiSimulation::divideLabor( ){ - int* globalIndex; - int nComponents; MoleculeStamp** compStamps; randomSPRNG *myRandom; int* componentsNmol; int* AtomsPerProc; + int* GroupsPerProc; double numerator; double denominator; double precast; double x, y, a; int old_atoms, add_atoms, new_atoms; + int old_groups, add_groups, new_groups; int nTarget; - int molIndex, atomIndex; + int molIndex, atomIndex, groupIndex; int done; - int i, j, loops, which_proc, nmol_local, natoms_local; - int nmol_global, natoms_global; + int i, j, loops, which_proc; + int nmol_global, nmol_local; + int ngroups_global, ngroups_local; + int natoms_global, natoms_local; + int ncutoff_groups, nAtomsInGroups; int local_index; int baseSeed = entryPlug->getSeed(); + CutoffGroupStamp* cg; nComponents = entryPlug->nComponents; compStamps = entryPlug->compStamps; componentsNmol = entryPlug->componentsNmol; - AtomsPerProc = new int[mpiPlug->numberProcessors]; + AtomsPerProc = new int[parallelData->nProcessors]; + GroupsPerProc = new int[parallelData->nProcessors]; - mpiPlug->nAtomsGlobal = entryPlug->n_atoms; - mpiPlug->nBondsGlobal = entryPlug->n_bonds; - mpiPlug->nBendsGlobal = entryPlug->n_bends; - mpiPlug->nTorsionsGlobal = entryPlug->n_torsions; - mpiPlug->nSRIGlobal = entryPlug->n_SRI; - mpiPlug->nMolGlobal = entryPlug->n_mol; + parallelData->nAtomsGlobal = entryPlug->n_atoms; + parallelData->nBondsGlobal = entryPlug->n_bonds; + parallelData->nBendsGlobal = entryPlug->n_bends; + parallelData->nTorsionsGlobal = entryPlug->n_torsions; + parallelData->nSRIGlobal = entryPlug->n_SRI; + parallelData->nGroupsGlobal = entryPlug->ngroup; + parallelData->nMolGlobal = entryPlug->n_mol; - myRandom = new randomSPRNG( baseSeed ); - a = 3.0 * (double)mpiPlug->nMolGlobal / (double)mpiPlug->nAtomsGlobal; + a = 3.0 * (double)parallelData->nMolGlobal / (double)parallelData->nAtomsGlobal; // Initialize things that we'll send out later: - for (i = 0; i < mpiPlug->numberProcessors; i++ ) { + for (i = 0; i < parallelData->nProcessors; i++ ) { AtomsPerProc[i] = 0; + GroupsPerProc[i] = 0; } - for (i = 0; i < mpiPlug->nMolGlobal; i++ ) { + for (i = 0; i < parallelData->nMolGlobal; i++ ) { // default to an error condition: MolToProcMap[i] = -1; MolComponentType[i] = -1; } - for (i = 0; i < mpiPlug->nAtomsGlobal; i++ ) { + for (i = 0; i < parallelData->nAtomsGlobal; i++ ) { // default to an error condition: AtomToProcMap[i] = -1; } + for (i = 0; i < parallelData->nGroupsGlobal; i++ ) { + // default to an error condition: + GroupToProcMap[i] = -1; + } - if (mpiPlug->myNode == 0) { + if (parallelData->myNode == 0) { numerator = (double) entryPlug->n_atoms; - denominator = (double) mpiPlug->numberProcessors; + denominator = (double) parallelData->nProcessors; precast = numerator / denominator; nTarget = (int)( precast + 0.5 ); @@ -111,6 +123,7 @@ int* mpiSimulation::divideLabor( void ){ } atomIndex = 0; + groupIndex = 0; for (i = 0; i < molIndex; i++ ) { @@ -122,7 +135,7 @@ int* mpiSimulation::divideLabor( void ){ // Pick a processor at random - which_proc = (int) (myRandom->getRandom() * mpiPlug->numberProcessors); + which_proc = (int) (myRandom->getRandom() * parallelData->nProcessors); // How many atoms does this processor have? @@ -130,6 +143,16 @@ int* mpiSimulation::divideLabor( void ){ add_atoms = compStamps[MolComponentType[i]]->getNAtoms(); new_atoms = old_atoms + add_atoms; + old_groups = GroupsPerProc[which_proc]; + ncutoff_groups = compStamps[MolComponentType[i]]->getNCutoffGroups(); + nAtomsInGroups = 0; + for (j = 0; j < ncutoff_groups; j++) { + cg = compStamps[MolComponentType[i]]->getCutoffGroup(j); + nAtomsInGroups += cg->getNMembers(); + } + add_groups = add_atoms - nAtomsInGroups + ncutoff_groups; + new_groups = old_groups + add_groups; + // If we've been through this loop too many times, we need // to just give up and assign the molecule to this processor // and be done with it. @@ -149,6 +172,11 @@ int* mpiSimulation::divideLabor( void ){ AtomToProcMap[atomIndex] = which_proc; atomIndex++; } + GroupsPerProc[which_proc] += add_groups; + for (j=0; j < add_groups; j++) { + GroupToProcMap[groupIndex] = which_proc; + groupIndex++; + } done = 1; continue; } @@ -163,6 +191,11 @@ int* mpiSimulation::divideLabor( void ){ AtomToProcMap[atomIndex] = which_proc; atomIndex++; } + GroupsPerProc[which_proc] += add_groups; + for (j=0; j < add_groups; j++) { + GroupToProcMap[groupIndex] = which_proc; + groupIndex++; + } done = 1; continue; } @@ -186,6 +219,11 @@ int* mpiSimulation::divideLabor( void ){ AtomToProcMap[atomIndex] = which_proc; atomIndex++; } + GroupsPerProc[which_proc] += add_groups; + for (j=0; j < add_groups; j++) { + GroupToProcMap[groupIndex] = which_proc; + groupIndex++; + } done = 1; continue; } else { @@ -195,59 +233,86 @@ int* mpiSimulation::divideLabor( void ){ } } + // Spray out this nonsense to all other processors: - MPI_Bcast(MolToProcMap, mpiPlug->nMolGlobal, + //std::cerr << "node 0 mol2proc = \n"; + //for (i = 0; i < parallelData->nMolGlobal; i++) + // std::cerr << i << "\t" << MolToProcMap[i] << "\n"; + + MPI_Bcast(MolToProcMap, parallelData->nMolGlobal, MPI_INT, 0, MPI_COMM_WORLD); - MPI_Bcast(AtomToProcMap, mpiPlug->nAtomsGlobal, + MPI_Bcast(AtomToProcMap, parallelData->nAtomsGlobal, MPI_INT, 0, MPI_COMM_WORLD); - MPI_Bcast(MolComponentType, mpiPlug->nMolGlobal, + MPI_Bcast(GroupToProcMap, parallelData->nGroupsGlobal, MPI_INT, 0, MPI_COMM_WORLD); - MPI_Bcast(AtomsPerProc, mpiPlug->numberProcessors, + MPI_Bcast(MolComponentType, parallelData->nMolGlobal, + MPI_INT, 0, MPI_COMM_WORLD); + + MPI_Bcast(AtomsPerProc, parallelData->nProcessors, MPI_INT, 0, MPI_COMM_WORLD); + + MPI_Bcast(GroupsPerProc, parallelData->nProcessors, + MPI_INT, 0, MPI_COMM_WORLD); } else { // Listen to your marching orders from processor 0: - MPI_Bcast(MolToProcMap, mpiPlug->nMolGlobal, + MPI_Bcast(MolToProcMap, parallelData->nMolGlobal, MPI_INT, 0, MPI_COMM_WORLD); - MPI_Bcast(AtomToProcMap, mpiPlug->nAtomsGlobal, + MPI_Bcast(AtomToProcMap, parallelData->nAtomsGlobal, MPI_INT, 0, MPI_COMM_WORLD); - MPI_Bcast(MolComponentType, mpiPlug->nMolGlobal, + MPI_Bcast(GroupToProcMap, parallelData->nGroupsGlobal, MPI_INT, 0, MPI_COMM_WORLD); + + MPI_Bcast(MolComponentType, parallelData->nMolGlobal, + MPI_INT, 0, MPI_COMM_WORLD); - MPI_Bcast(AtomsPerProc, mpiPlug->numberProcessors, + MPI_Bcast(AtomsPerProc, parallelData->nProcessors, MPI_INT, 0, MPI_COMM_WORLD); + MPI_Bcast(GroupsPerProc, parallelData->nProcessors, + MPI_INT, 0, MPI_COMM_WORLD); + } - // Let's all check for sanity: nmol_local = 0; - for (i = 0 ; i < mpiPlug->nMolGlobal; i++ ) { - if (MolToProcMap[i] == mpiPlug->myNode) { + for (i = 0 ; i < parallelData->nMolGlobal; i++ ) { + if (MolToProcMap[i] == parallelData->myNode) { nmol_local++; } } natoms_local = 0; - for (i = 0; i < mpiPlug->nAtomsGlobal; i++) { - if (AtomToProcMap[i] == mpiPlug->myNode) { + for (i = 0; i < parallelData->nAtomsGlobal; i++) { + if (AtomToProcMap[i] == parallelData->myNode) { natoms_local++; } } + ngroups_local = 0; + for (i = 0; i < parallelData->nGroupsGlobal; i++) { + if (GroupToProcMap[i] == parallelData->myNode) { + ngroups_local++; + } + } + MPI_Allreduce(&nmol_local,&nmol_global,1,MPI_INT,MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&natoms_local,&natoms_global,1,MPI_INT, MPI_SUM, MPI_COMM_WORLD); + + MPI_Allreduce(&ngroups_local,&ngroups_global,1,MPI_INT, + MPI_SUM, MPI_COMM_WORLD); if( nmol_global != entryPlug->n_mol ){ sprintf( painCave.errMsg, @@ -267,39 +332,91 @@ int* mpiSimulation::divideLabor( void ){ simError(); } + if( ngroups_global != entryPlug->ngroup ){ + sprintf( painCave.errMsg, + "The sum of all ngroups_local, %d, did not equal the " + "total number of cutoffGroups, %d.\n", + ngroups_global, entryPlug->ngroup ); + painCave.isFatal = 1; + simError(); + } + sprintf( checkPointMsg, "Successfully divided the molecules among the processors.\n" ); MPIcheckPoint(); - mpiPlug->myNMol = nmol_local; - mpiPlug->myNlocal = natoms_local; + parallelData->nMolLocal = nmol_local; + parallelData->nAtomsLocal = natoms_local; + parallelData->nGroupsLocal = ngroups_local; - globalIndex = new int[mpiPlug->myNlocal]; + globalAtomIndex.resize(parallelData->nAtomsLocal); + globalToLocalAtom.resize(parallelData->nAtomsGlobal); local_index = 0; - for (i = 0; i < mpiPlug->nAtomsGlobal; i++) { - if (AtomToProcMap[i] == mpiPlug->myNode) { - globalIndex[local_index] = i; + for (i = 0; i < parallelData->nAtomsGlobal; i++) { + if (AtomToProcMap[i] == parallelData->myNode) { + globalAtomIndex[local_index] = i; + + globalToLocalAtom[i] = local_index; local_index++; + } + else + globalToLocalAtom[i] = -1; } + + globalGroupIndex.resize(parallelData->nGroupsLocal); + globalToLocalGroup.resize(parallelData->nGroupsGlobal); + local_index = 0; + for (i = 0; i < parallelData->nGroupsGlobal; i++) { + if (GroupToProcMap[i] == parallelData->myNode) { + globalGroupIndex[local_index] = i; + + globalToLocalGroup[i] = local_index; + local_index++; + + } + else + globalToLocalGroup[i] = -1; + } + + globalMolIndex.resize(parallelData->nMolLocal); + globalToLocalMol.resize(parallelData->nMolGlobal); + local_index = 0; + for (i = 0; i < parallelData->nMolGlobal; i++) { + if (MolToProcMap[i] == parallelData->myNode) { + globalMolIndex[local_index] = i; + globalToLocalMol[i] = local_index; + local_index++; + } + else + globalToLocalMol[i] = -1; + } - return globalIndex; } void mpiSimulation::mpiRefresh( void ){ int isError, i; - int *globalIndex = new int[mpiPlug->myNlocal]; + int *localToGlobalAtomIndex = new int[parallelData->nAtomsLocal]; + int *localToGlobalGroupIndex = new int[parallelData->nGroupsLocal]; - // Fortran indexing needs to be increased by 1 in order to get the 2 languages to - // not barf + // Fortran indexing needs to be increased by 1 in order to get the 2 + // languages to not barf - for(i=0; imyNlocal; i++) globalIndex[i] = entryPlug->atoms[i]->getGlobalIndex()+1; + for(i = 0; i < parallelData->nAtomsLocal; i++) + localToGlobalAtomIndex[i] = globalAtomIndex[i] + 1; + for(i = 0; i < parallelData->nGroupsLocal; i++) + localToGlobalGroupIndex[i] = globalGroupIndex[i] + 1; isError = 0; - setFsimParallel( mpiPlug, &(entryPlug->n_atoms), globalIndex, &isError ); + + setFsimParallel( parallelData, + &(parallelData->nAtomsLocal), localToGlobalAtomIndex, + &(parallelData->nGroupsLocal), localToGlobalGroupIndex, + &isError ); + if( isError ){ sprintf( painCave.errMsg, @@ -308,8 +425,10 @@ void mpiSimulation::mpiRefresh( void ){ simError(); } - delete[] globalIndex; + delete[] localToGlobalGroupIndex; + delete[] localToGlobalAtomIndex; + sprintf( checkPointMsg, " mpiRefresh successful.\n" ); MPIcheckPoint();