Ignore:
File:
1 edited

Legend:

Unmodified
Added
Removed

  • src/config.cpp

    ra67d19 r68f03d  
    88#include <cstring>
    99
     10#include "World.hpp"
    1011#include "atom.hpp"
    1112#include "bond.hpp"
    … …  
    500501//        case 'j': // BoxLength
    501502//          Log() << Verbose(0) << "enter lower triadiagonalo form of basis matrix" << endl << endl;
    502 //          double * const cell_size = World::get()->cell_size;
     503//          double * const cell_size = World::getInstance().getDomain();
    503504//          for (int i=0;i<6;i++) {
    504505//            Log() << Verbose(0) << "Cell size" << i << ": ";
    … …  
    676677{
    677678  int MaxTypes = 0;
    678   element *elementhash[MAX_ELEMENTS];
     679  const element *elementhash[MAX_ELEMENTS];
    679680  char name[MAX_ELEMENTS];
    680681  char keyword[MAX_ELEMENTS];
    … …  
    734735            sprintf(keyword,"%s_%i",name, j+1);
    735736            if (repetition == 0) {
    736               neues = new atom();
     737              neues = World::getInstance().createAtom();
    737738              AtomList[i][j] = neues;
    738739              LinearList[ FileBuffer->LineMapping[FileBuffer->CurrentLine] ] = neues;
    … …  
    741742              neues = AtomList[i][j];
    742743            status = (status &&
    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) &&
     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) &&
    746747                    ParseForParameter(verbose,FileBuffer, keyword, 0, 4, 1, int_type, &neues->FixedIon, 1, (repetition == 0) ? critical : optional));
    747748            if (!status) break;
    … …  
    757758            // put into trajectories list
    758759            for (int d=0;d<NDIM;d++)
    759               neues->Trajectory.R.at(repetition).x[d] = neues->x.x[d];
     760              neues->Trajectory.R.at(repetition)[d] = neues->x[d];
    760761
    761762            // parse velocities if present
    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.;
     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.;
    768769            for (int d=0;d<NDIM;d++)
    769               neues->Trajectory.U.at(repetition).x[d] = neues->v.x[d];
     770              neues->Trajectory.U.at(repetition)[d] = neues->v[d];
    770771
    771772            // parse forces if present
    … …  
    777778              value[2] = 0.;
    778779            for (int d=0;d<NDIM;d++)
    779               neues->Trajectory.F.at(repetition).x[d] = value[d];
     780              neues->Trajectory.F.at(repetition)[d] = value[d];
    780781
    781782  //            Log() << Verbose(0) << "Parsed position of step " << (repetition) << ": (";
    … …  
    813814          sprintf(keyword,"%s_%i",name, j+1);
    814815          if (repetition == 0) {
    815             neues = new atom();
     816            neues = World::getInstance().createAtom();
    816817            AtomList[i][j] = neues;
    817818            LinearList[ FileBuffer->LineMapping[FileBuffer->CurrentLine] ] = neues;
    … …  
    820821            neues = AtomList[i][j];
    821822          // then parse for each atom the coordinates as often as present
    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);
     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);
    825826          ParseForParameter(verbose,FileBuffer, keyword, 0, 4, 1, int_type, &neues->FixedIon, repetition,critical);
    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.;
     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.;
    832833          // here we don't care if forces are present (last in trajectories is always equal to current position)
    833834          neues->type = elementhash[i]; // find element type
    … …  
    852853void config::Load(const char * const filename, const string &BondGraphFileName, const periodentafel * const periode, MoleculeListClass * const &MolList)
    853854{
    854   molecule *mol = new molecule(periode);
     855  molecule *mol = World::getInstance().createMolecule();
    855856  ifstream *file = new ifstream(filename);
    856857  if (file == NULL) {
    … …  
    967968  // Unit cell and magnetic field
    968969  ParseForParameter(verbose,FileBuffer, "BoxLength", 0, 3, 3, lower_trigrid, BoxLength, 1, critical); /* Lattice->RealBasis */
    969   double * const cell_size = World::get()->cell_size;
     970  double * const cell_size = World::getInstance().getDomain();
    970971  cell_size[0] = BoxLength[0];
    971972  cell_size[1] = BoxLength[3];
    … …  
    10911092void config::LoadOld(const char * const filename, const string &BondGraphFileName, const periodentafel * const periode, MoleculeListClass * const &MolList)
    10921093{
    1093   molecule *mol = new molecule(periode);
     1094  molecule *mol = World::getInstance().createMolecule();
    10941095  ifstream *file = new ifstream(filename);
    10951096  if (file == NULL) {
    … …  
    11091110  string zeile;
    11101111  string dummy;
    1111   element *elementhash[128];
     1112  const element *elementhash[128];
    11121113  int Z = -1;
    11131114  int No = -1;
    … …  
    11721173
    11731174  ParseForParameter(verbose,file, "BoxLength", 0, 3, 3, lower_trigrid, BoxLength, 1, critical); /* Lattice->RealBasis */
    1174   double * const cell_size = World::get()->cell_size;
     1175  double * const cell_size = World::getInstance().getDomain();
    11751176  cell_size[0] = BoxLength[0];
    11761177  cell_size[1] = BoxLength[3];
    … …  
    12911292        }
    12921293        istringstream input2(zeile);
    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
     1294        atom *neues = World::getInstance().createAtom();
     1295        input2 >> neues->x[0]; // x
     1296        input2 >> neues->x[1]; // y
     1297        input2 >> neues->x[2]; // z
    12971298        input2 >> l;
    12981299        neues->type = elementhash[No]; // find element type
    … …  
    13161317  // bring MaxTypes up to date
    13171318  mol->CountElements();
    1318   const double * const cell_size = World::get()->cell_size;
     1319  const double * const cell_size = World::getInstance().getDomain();
    13191320  ofstream * const output = new ofstream(filename, ios::out);
    13201321  if (output != NULL) {
    … …  
    15741575             'a'+(unsigned char)(AtomNo % 26),           /* letter for chain */
    15751576             MolNo,         /* residue sequence number */
    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 */
     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 */
    15791580             (double)Walker->type->Valence,         /* occupancy */
    15801581             (double)Walker->type->NoValenceOrbitals,          /* temperature factor */
    … …  
    16281629           'a'+(unsigned char)(AtomNo % 26),           /* letter for chain */
    16291630           0,         /* residue sequence number */
    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 */
     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 */
    16331634           (double)Walker->type->Valence,         /* occupancy */
    16341635           (double)Walker->type->NoValenceOrbitals,          /* temperature factor */
    … …  
    16791680    Walker = Walker->next;
    16801681    *output << Walker->nr << "\t";
    1681     *output << Walker->Name << "\t";
     1682    *output << Walker->getName() << "\t";
    16821683    *output << mol->name << "\t";
    16831684    *output << 0 << "\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";
     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";
    16861687    *output << Walker->type->symbol << "\t";
    16871688    for (BondList::iterator runner = Walker->ListOfBonds.begin(); runner != Walker->ListOfBonds.end(); runner++)
    … …  
    17551756        Walker = Walker->next;
    17561757        *output << AtomNo+1 << "\t";
    1757         *output << Walker->Name << "\t";
     1758        *output << Walker->getName() << "\t";
    17581759        *output << (*MolWalker)->name << "\t";
    17591760        *output << MolCounter+1 << "\t";
    1760         *output << Walker->node->x[0] << "\t" << Walker->node->x[1] << "\t" << Walker->node->x[2] << "\t";
     1761        *output << Walker->node->at(0) << "\t" << Walker->node->at(1) << "\t" << Walker->node->at(2) << "\t";
    17611762        *output << (double)Walker->type->Valence << "\t";
    17621763        *output << Walker->type->symbol << "\t";
    … …  
    17821783
    17831784  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 */
     1793void 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);
    17841909};
    17851910
Note: See TracChangeset for help on using the changeset viewer.