C++ Species::reservoir方法代码示例

   /*
   * Evaluate Rosenbluth weight, and add to accumulator.
   */
   void McChemicalPotential::sample(long iStep)
   {
      if (isAtInterval(iStep))  {

         Species* speciesPtr;
         Molecule* molPtr;
         Molecule::BondIterator bondIter;
         Atom* endPtr;
         double w;
         double rosenbluth = 1;
         double de;
         double e = 0;

         speciesPtr = &(simulation().species(speciesId_));

         // Pop a new molecule off the species reservoir
         molPtr = &(speciesPtr->reservoir().pop());
         system().addMolecule(*molPtr);

         // Loop over molecule growth trials
         for (int i = 0; i < nMoleculeTrial_; i++) {

            // Pick a random position for the first atom
            endPtr = &molPtr->atom(0);
            boundary().randomPosition(random(), endPtr->position());

            e = system().pairPotential().atomEnergy(*endPtr);
            rosenbluth = boltzmann(e);
            system().pairPotential().addAtom(*endPtr);

            for (molPtr->begin(bondIter); bondIter.notEnd(); ++bondIter) {
                addEndAtom(&(bondIter->atom(1)), &(bondIter->atom(0)), bondIter->typeId(), w, de);
                e += de;
                rosenbluth *= w;
                system().pairPotential().addAtom(bondIter->atom(1));
            }

            rosenbluth = rosenbluth / pow(nTrial_,molPtr->nAtom()-1);
            accumulator_.sample(rosenbluth, outputFile_);

            system().pairPotential().deleteAtom(*endPtr);
            for (molPtr->begin(bondIter); bondIter.notEnd(); ++bondIter) {
                system().pairPotential().deleteAtom(bondIter->atom(1));
            }
         }

         // Return additional molecule to reservoir
         system().removeMolecule(*molPtr);
         speciesPtr->reservoir().push(*molPtr);
      }

   }

   /*
   * Initialize all Molecule and Atom objects for one Species (private).
   *
   * This function creates associations between Species, Molecule, and
   * Atom objects for all molecules of one species, and sets atom typeIds.
   *
   * For each molecule, it sets the id, species pointer, nAtom, and the 
   * firstAtom pointer. The molecule id is only unique within each species.
   *
   * For each atom, it sets the molecule pointer and an integer typeId.
   *
   * This method also pushes all molecules of the species onto the
   * reservoir, pushing them in order of decreasing molecule id.
   */
   void Simulation::initializeSpecies(int iSpecies)
   {

      Species*  speciesPtr;
      Molecule* moleculePtr;
      Atom*     atomPtr;
      int       iMol, iAtom;
      int       capacity, nAtom;

      speciesPtr = &species(iSpecies);
      capacity   = speciesPtr->capacity();
      nAtom      = speciesPtr->nAtom();

      // Initialize pointers before loop
      moleculePtr = &molecules_[firstMoleculeIds_[iSpecies]];
      atomPtr     = &atoms_[firstAtomIds_[iSpecies]];

      // Loop over all molecules in Species
      for (iMol = 0; iMol < capacity; ++iMol) {

         // Initialize a Molecule
         moleculePtr->setId(iMol);
         moleculePtr->setSpecies(*speciesPtr);
         moleculePtr->setNAtom(nAtom);
         moleculePtr->setFirstAtom(*atomPtr);

         // Loop over atoms in a molecule, set molecule and atom TypeId
         for (iAtom = 0; iAtom < nAtom; ++iAtom) {
            atomPtr->setMolecule(*moleculePtr);
            atomPtr->setTypeId(speciesPtr->atomTypeId(iAtom));
            ++atomPtr;
         }

         ++moleculePtr;
      }

      // Push all molecules of this species onto the reservoir stack
      // Push on in reverse order, so that they pop off in sequence
      moleculePtr = &molecules_[firstMoleculeIds_[iSpecies] + capacity - 1];
      for (iMol = 0; iMol < capacity; ++iMol) {
         speciesPtr->reservoir().push(*moleculePtr);
         --moleculePtr;
      }

   }