Ignore:
File:
1 edited

Legend:

Unmodified
Added
Removed
  • src/config.cpp

    r68f03d ra67d19  
    88#include <cstring>
    99
    10 #include "World.hpp"
    1110#include "atom.hpp"
    1211#include "bond.hpp"
     
    501500//        case 'j': // BoxLength
    502501//          Log() << Verbose(0) << "enter lower triadiagonalo form of basis matrix" << endl << endl;
    503 //          double * const cell_size = World::getInstance().getDomain();
     502//          double * const cell_size = World::get()->cell_size;
    504503//          for (int i=0;i<6;i++) {
    505504//            Log() << Verbose(0) << "Cell size" << i << ": ";
     
    677676{
    678677  int MaxTypes = 0;
    679   const element *elementhash[MAX_ELEMENTS];
     678  element *elementhash[MAX_ELEMENTS];
    680679  char name[MAX_ELEMENTS];
    681680  char keyword[MAX_ELEMENTS];
     
    735734            sprintf(keyword,"%s_%i",name, j+1);
    736735            if (repetition == 0) {
    737               neues = World::getInstance().createAtom();
     736              neues = new atom();
    738737              AtomList[i][j] = neues;
    739738              LinearList[ FileBuffer->LineMapping[FileBuffer->CurrentLine] ] = neues;
     
    742741              neues = AtomList[i][j];
    743742            status = (status &&
    744                     ParseForParameter(verbose,FileBuffer, keyword, 0, 1, 1, double_type, &neues->x[0], 1, (repetition == 0) ? critical : optional) &&
    745                     ParseForParameter(verbose,FileBuffer, keyword, 0, 2, 1, double_type, &neues->x[1], 1, (repetition == 0) ? critical : optional) &&
    746                     ParseForParameter(verbose,FileBuffer, keyword, 0, 3, 1, double_type, &neues->x[2], 1, (repetition == 0) ? critical : optional) &&
     743                    ParseForParameter(verbose,FileBuffer, keyword, 0, 1, 1, double_type, &neues->x.x[0], 1, (repetition == 0) ? critical : optional) &&
     744                    ParseForParameter(verbose,FileBuffer, keyword, 0, 2, 1, double_type, &neues->x.x[1], 1, (repetition == 0) ? critical : optional) &&
     745                    ParseForParameter(verbose,FileBuffer, keyword, 0, 3, 1, double_type, &neues->x.x[2], 1, (repetition == 0) ? critical : optional) &&
    747746                    ParseForParameter(verbose,FileBuffer, keyword, 0, 4, 1, int_type, &neues->FixedIon, 1, (repetition == 0) ? critical : optional));
    748747            if (!status) break;
     
    758757            // put into trajectories list
    759758            for (int d=0;d<NDIM;d++)
    760               neues->Trajectory.R.at(repetition)[d] = neues->x[d];
     759              neues->Trajectory.R.at(repetition).x[d] = neues->x.x[d];
    761760
    762761            // parse velocities if present
    763             if(!ParseForParameter(verbose,FileBuffer, keyword, 0, 5, 1, double_type, &neues->v[0], 1,optional))
    764               neues->v[0] = 0.;
    765             if(!ParseForParameter(verbose,FileBuffer, keyword, 0, 6, 1, double_type, &neues->v[1], 1,optional))
    766               neues->v[1] = 0.;
    767             if(!ParseForParameter(verbose,FileBuffer, keyword, 0, 7, 1, double_type, &neues->v[2], 1,optional))
    768               neues->v[2] = 0.;
     762            if(!ParseForParameter(verbose,FileBuffer, keyword, 0, 5, 1, double_type, &neues->v.x[0], 1,optional))
     763              neues->v.x[0] = 0.;
     764            if(!ParseForParameter(verbose,FileBuffer, keyword, 0, 6, 1, double_type, &neues->v.x[1], 1,optional))
     765              neues->v.x[1] = 0.;
     766            if(!ParseForParameter(verbose,FileBuffer, keyword, 0, 7, 1, double_type, &neues->v.x[2], 1,optional))
     767              neues->v.x[2] = 0.;
    769768            for (int d=0;d<NDIM;d++)
    770               neues->Trajectory.U.at(repetition)[d] = neues->v[d];
     769              neues->Trajectory.U.at(repetition).x[d] = neues->v.x[d];
    771770
    772771            // parse forces if present
     
    778777              value[2] = 0.;
    779778            for (int d=0;d<NDIM;d++)
    780               neues->Trajectory.F.at(repetition)[d] = value[d];
     779              neues->Trajectory.F.at(repetition).x[d] = value[d];
    781780
    782781  //            Log() << Verbose(0) << "Parsed position of step " << (repetition) << ": (";
     
    814813          sprintf(keyword,"%s_%i",name, j+1);
    815814          if (repetition == 0) {
    816             neues = World::getInstance().createAtom();
     815            neues = new atom();
    817816            AtomList[i][j] = neues;
    818817            LinearList[ FileBuffer->LineMapping[FileBuffer->CurrentLine] ] = neues;
     
    821820            neues = AtomList[i][j];
    822821          // then parse for each atom the coordinates as often as present
    823           ParseForParameter(verbose,FileBuffer, keyword, 0, 1, 1, double_type, &neues->x[0], repetition,critical);
    824           ParseForParameter(verbose,FileBuffer, keyword, 0, 2, 1, double_type, &neues->x[1], repetition,critical);
    825           ParseForParameter(verbose,FileBuffer, keyword, 0, 3, 1, double_type, &neues->x[2], repetition,critical);
     822          ParseForParameter(verbose,FileBuffer, keyword, 0, 1, 1, double_type, &neues->x.x[0], repetition,critical);
     823          ParseForParameter(verbose,FileBuffer, keyword, 0, 2, 1, double_type, &neues->x.x[1], repetition,critical);
     824          ParseForParameter(verbose,FileBuffer, keyword, 0, 3, 1, double_type, &neues->x.x[2], repetition,critical);
    826825          ParseForParameter(verbose,FileBuffer, keyword, 0, 4, 1, int_type, &neues->FixedIon, repetition,critical);
    827           if(!ParseForParameter(verbose,FileBuffer, keyword, 0, 5, 1, double_type, &neues->v[0], repetition,optional))
    828             neues->v[0] = 0.;
    829           if(!ParseForParameter(verbose,FileBuffer, keyword, 0, 6, 1, double_type, &neues->v[1], repetition,optional))
    830             neues->v[1] = 0.;
    831           if(!ParseForParameter(verbose,FileBuffer, keyword, 0, 7, 1, double_type, &neues->v[2], repetition,optional))
    832             neues->v[2] = 0.;
     826          if(!ParseForParameter(verbose,FileBuffer, keyword, 0, 5, 1, double_type, &neues->v.x[0], repetition,optional))
     827            neues->v.x[0] = 0.;
     828          if(!ParseForParameter(verbose,FileBuffer, keyword, 0, 6, 1, double_type, &neues->v.x[1], repetition,optional))
     829            neues->v.x[1] = 0.;
     830          if(!ParseForParameter(verbose,FileBuffer, keyword, 0, 7, 1, double_type, &neues->v.x[2], repetition,optional))
     831            neues->v.x[2] = 0.;
    833832          // here we don't care if forces are present (last in trajectories is always equal to current position)
    834833          neues->type = elementhash[i]; // find element type
     
    853852void config::Load(const char * const filename, const string &BondGraphFileName, const periodentafel * const periode, MoleculeListClass * const &MolList)
    854853{
    855   molecule *mol = World::getInstance().createMolecule();
     854  molecule *mol = new molecule(periode);
    856855  ifstream *file = new ifstream(filename);
    857856  if (file == NULL) {
     
    968967  // Unit cell and magnetic field
    969968  ParseForParameter(verbose,FileBuffer, "BoxLength", 0, 3, 3, lower_trigrid, BoxLength, 1, critical); /* Lattice->RealBasis */
    970   double * const cell_size = World::getInstance().getDomain();
     969  double * const cell_size = World::get()->cell_size;
    971970  cell_size[0] = BoxLength[0];
    972971  cell_size[1] = BoxLength[3];
     
    10921091void config::LoadOld(const char * const filename, const string &BondGraphFileName, const periodentafel * const periode, MoleculeListClass * const &MolList)
    10931092{
    1094   molecule *mol = World::getInstance().createMolecule();
     1093  molecule *mol = new molecule(periode);
    10951094  ifstream *file = new ifstream(filename);
    10961095  if (file == NULL) {
     
    11101109  string zeile;
    11111110  string dummy;
    1112   const element *elementhash[128];
     1111  element *elementhash[128];
    11131112  int Z = -1;
    11141113  int No = -1;
     
    11731172
    11741173  ParseForParameter(verbose,file, "BoxLength", 0, 3, 3, lower_trigrid, BoxLength, 1, critical); /* Lattice->RealBasis */
    1175   double * const cell_size = World::getInstance().getDomain();
     1174  double * const cell_size = World::get()->cell_size;
    11761175  cell_size[0] = BoxLength[0];
    11771176  cell_size[1] = BoxLength[3];
     
    12921291        }
    12931292        istringstream input2(zeile);
    1294         atom *neues = World::getInstance().createAtom();
    1295         input2 >> neues->x[0]; // x
    1296         input2 >> neues->x[1]; // y
    1297         input2 >> neues->x[2]; // z
     1293        atom *neues = new atom();
     1294        input2 >> neues->x.x[0]; // x
     1295        input2 >> neues->x.x[1]; // y
     1296        input2 >> neues->x.x[2]; // z
    12981297        input2 >> l;
    12991298        neues->type = elementhash[No]; // find element type
     
    13171316  // bring MaxTypes up to date
    13181317  mol->CountElements();
    1319   const double * const cell_size = World::getInstance().getDomain();
     1318  const double * const cell_size = World::get()->cell_size;
    13201319  ofstream * const output = new ofstream(filename, ios::out);
    13211320  if (output != NULL) {
     
    15751574             'a'+(unsigned char)(AtomNo % 26),           /* letter for chain */
    15761575             MolNo,         /* residue sequence number */
    1577              Walker->node->at(0),                 /* position X in Angstroem */
    1578              Walker->node->at(1),                 /* position Y in Angstroem */
    1579              Walker->node->at(2),                 /* position Z in Angstroem */
     1576             Walker->node->x[0],                 /* position X in Angstroem */
     1577             Walker->node->x[1],                 /* position Y in Angstroem */
     1578             Walker->node->x[2],                 /* position Z in Angstroem */
    15801579             (double)Walker->type->Valence,         /* occupancy */
    15811580             (double)Walker->type->NoValenceOrbitals,          /* temperature factor */
     
    16291628           'a'+(unsigned char)(AtomNo % 26),           /* letter for chain */
    16301629           0,         /* residue sequence number */
    1631            Walker->node->at(0),                 /* position X in Angstroem */
    1632            Walker->node->at(1),                 /* position Y in Angstroem */
    1633            Walker->node->at(2),                 /* position Z in Angstroem */
     1630           Walker->node->x[0],                 /* position X in Angstroem */
     1631           Walker->node->x[1],                 /* position Y in Angstroem */
     1632           Walker->node->x[2],                 /* position Z in Angstroem */
    16341633           (double)Walker->type->Valence,         /* occupancy */
    16351634           (double)Walker->type->NoValenceOrbitals,          /* temperature factor */
     
    16801679    Walker = Walker->next;
    16811680    *output << Walker->nr << "\t";
    1682     *output << Walker->getName() << "\t";
     1681    *output << Walker->Name << "\t";
    16831682    *output << mol->name << "\t";
    16841683    *output << 0 << "\t";
    1685     *output << Walker->node->at(0) << "\t" << Walker->node->at(1) << "\t" << Walker->node->at(2) << "\t";
    1686     *output << static_cast<double>(Walker->type->Valence) << "\t";
     1684    *output << Walker->node->x[0] << "\t" << Walker->node->x[1] << "\t" << Walker->node->x[2] << "\t";
     1685    *output << (double)Walker->type->Valence << "\t";
    16871686    *output << Walker->type->symbol << "\t";
    16881687    for (BondList::iterator runner = Walker->ListOfBonds.begin(); runner != Walker->ListOfBonds.end(); runner++)
     
    17561755        Walker = Walker->next;
    17571756        *output << AtomNo+1 << "\t";
    1758         *output << Walker->getName() << "\t";
     1757        *output << Walker->Name << "\t";
    17591758        *output << (*MolWalker)->name << "\t";
    17601759        *output << MolCounter+1 << "\t";
    1761         *output << Walker->node->at(0) << "\t" << Walker->node->at(1) << "\t" << Walker->node->at(2) << "\t";
     1760        *output << Walker->node->x[0] << "\t" << Walker->node->x[1] << "\t" << Walker->node->x[2] << "\t";
    17621761        *output << (double)Walker->type->Valence << "\t";
    17631762        *output << Walker->type->symbol << "\t";
     
    17831782
    17841783  return true;
    1785 };
    1786 
    1787 
    1788 /** Tries given filename or standard on saving the config file.
    1789  * \param *ConfigFileName name of file
    1790  * \param *periode pointer to periodentafel structure with all the elements
    1791  * \param *molecules list of molecules structure with all the atoms and coordinates
    1792  */
    1793 void config::SaveAll(char *ConfigFileName, periodentafel *periode, MoleculeListClass *molecules)
    1794 {
    1795   char filename[MAXSTRINGSIZE];
    1796   ofstream output;
    1797   molecule *mol = World::getInstance().createMolecule();
    1798   mol->SetNameFromFilename(ConfigFileName);
    1799 
    1800   if (!strcmp(configpath, GetDefaultPath())) {
    1801     eLog() << Verbose(2) << "config is found under different path then stated in config file::defaultpath!" << endl;
    1802   }
    1803 
    1804 
    1805   // first save as PDB data
    1806   if (ConfigFileName != NULL)
    1807     strcpy(filename, ConfigFileName);
    1808   if (output == NULL)
    1809     strcpy(filename,"main_pcp_linux");
    1810   Log() << Verbose(0) << "Saving as pdb input ";
    1811   if (SavePDB(filename, molecules))
    1812     Log() << Verbose(0) << "done." << endl;
    1813   else
    1814     Log() << Verbose(0) << "failed." << endl;
    1815 
    1816   // then save as tremolo data file
    1817   if (ConfigFileName != NULL)
    1818     strcpy(filename, ConfigFileName);
    1819   if (output == NULL)
    1820     strcpy(filename,"main_pcp_linux");
    1821   Log() << Verbose(0) << "Saving as tremolo data input ";
    1822   if (SaveTREMOLO(filename, molecules))
    1823     Log() << Verbose(0) << "done." << endl;
    1824   else
    1825     Log() << Verbose(0) << "failed." << endl;
    1826 
    1827   // translate each to its center and merge all molecules in MoleculeListClass into this molecule
    1828   int N = molecules->ListOfMolecules.size();
    1829   int *src = new int[N];
    1830   N=0;
    1831   for (MoleculeList::iterator ListRunner = molecules->ListOfMolecules.begin(); ListRunner != molecules->ListOfMolecules.end(); ListRunner++) {
    1832     src[N++] = (*ListRunner)->IndexNr;
    1833     (*ListRunner)->Translate(&(*ListRunner)->Center);
    1834   }
    1835   molecules->SimpleMultiAdd(mol, src, N);
    1836   delete[](src);
    1837 
    1838   // ... and translate back
    1839   for (MoleculeList::iterator ListRunner = molecules->ListOfMolecules.begin(); ListRunner != molecules->ListOfMolecules.end(); ListRunner++) {
    1840     (*ListRunner)->Center.Scale(-1.);
    1841     (*ListRunner)->Translate(&(*ListRunner)->Center);
    1842     (*ListRunner)->Center.Scale(-1.);
    1843   }
    1844 
    1845   Log() << Verbose(0) << "Storing configuration ... " << endl;
    1846   // get correct valence orbitals
    1847   mol->CalculateOrbitals(*this);
    1848   InitMaxMinStopStep = MaxMinStopStep = MaxPsiDouble;
    1849   if (ConfigFileName != NULL) { // test the file name
    1850     strcpy(filename, ConfigFileName);
    1851     output.open(filename, ios::trunc);
    1852   } else if (strlen(configname) != 0) {
    1853     strcpy(filename, configname);
    1854     output.open(configname, ios::trunc);
    1855     } else {
    1856       strcpy(filename, DEFAULTCONFIG);
    1857       output.open(DEFAULTCONFIG, ios::trunc);
    1858     }
    1859   output.close();
    1860   output.clear();
    1861   Log() << Verbose(0) << "Saving of config file ";
    1862   if (Save(filename, periode, mol))
    1863     Log() << Verbose(0) << "successful." << endl;
    1864   else
    1865     Log() << Verbose(0) << "failed." << endl;
    1866 
    1867   // and save to xyz file
    1868   if (ConfigFileName != NULL) {
    1869     strcpy(filename, ConfigFileName);
    1870     strcat(filename, ".xyz");
    1871     output.open(filename, ios::trunc);
    1872   }
    1873   if (output == NULL) {
    1874     strcpy(filename,"main_pcp_linux");
    1875     strcat(filename, ".xyz");
    1876     output.open(filename, ios::trunc);
    1877   }
    1878   Log() << Verbose(0) << "Saving of XYZ file ";
    1879   if (mol->MDSteps <= 1) {
    1880     if (mol->OutputXYZ(&output))
    1881       Log() << Verbose(0) << "successful." << endl;
    1882     else
    1883       Log() << Verbose(0) << "failed." << endl;
    1884   } else {
    1885     if (mol->OutputTrajectoriesXYZ(&output))
    1886       Log() << Verbose(0) << "successful." << endl;
    1887     else
    1888       Log() << Verbose(0) << "failed." << endl;
    1889   }
    1890   output.close();
    1891   output.clear();
    1892 
    1893   // and save as MPQC configuration
    1894   if (ConfigFileName != NULL)
    1895     strcpy(filename, ConfigFileName);
    1896   if (output == NULL)
    1897     strcpy(filename,"main_pcp_linux");
    1898   Log() << Verbose(0) << "Saving as mpqc input ";
    1899   if (SaveMPQC(filename, mol))
    1900     Log() << Verbose(0) << "done." << endl;
    1901   else
    1902     Log() << Verbose(0) << "failed." << endl;
    1903 
    1904   if (!strcmp(configpath, GetDefaultPath())) {
    1905     eLog() << Verbose(2) << "config is found under different path then stated in config file::defaultpath!" << endl;
    1906   }
    1907 
    1908   World::getInstance().destroyMolecule(mol);
    19091784};
    19101785
Note: See TracChangeset for help on using the changeset viewer.