Ignore:
File:
1 edited

Legend:

Unmodified
Added
Removed

  • src/builder.cpp

    rf65e1f r9907e7  
    5353
    5454#include <cstring>
     55#include <cstdlib>
    5556
    5657#include "analysis_bonds.hpp"
    … …  
    864865
    865866        mol->CountAtoms(); // recount atoms
    866         if (mol->AtomCount != 0) {  // if there is more than none
    867           count = mol->AtomCount;  // is changed becausing of adding, thus has to be stored away beforehand
     867        if (mol->getAtomCount() != 0) {  // if there is more than none
     868          count = mol->getAtomCount();  // is changed becausing of adding, thus has to be stored away beforehand
    868869          Elements = new element *[count];
    869870          vectors = new Vector *[count];
    … …  
    12951296  // generate some KeySets
    12961297  DoLog(0) && (Log() << Verbose(0) << "Generating KeySets." << endl);
    1297   KeySet TestSets[mol->AtomCount+1];
     1298  KeySet TestSets[mol->getAtomCount()+1];
    12981299  i=1;
    12991300  while (Walker->next != mol->end) {
    … …  
    13061307  DoLog(0) && (Log() << Verbose(0) << "Testing insertion of already present item in KeySets." << endl);
    13071308  KeySetTestPair test;
    1308   test = TestSets[mol->AtomCount-1].insert(Walker->nr);
     1309  test = TestSets[mol->getAtomCount()-1].insert(Walker->nr);
    13091310  if (test.second) {
    13101311    DoLog(1) && (Log() << Verbose(1) << "Insertion worked?!" << endl);
    … …  
    13121313    DoLog(1) && (Log() << Verbose(1) << "Insertion rejected: Present object is " << (*test.first) << "." << endl);
    13131314  }
    1314   TestSets[mol->AtomCount].insert(mol->end->previous->nr);
    1315   TestSets[mol->AtomCount].insert(mol->end->previous->previous->previous->nr);
     1315  TestSets[mol->getAtomCount()].insert(mol->end->previous->nr);
     1316  TestSets[mol->getAtomCount()].insert(mol->end->previous->previous->previous->nr);
    13161317
    13171318  // constructing Graph structure
    … …  
    13211322  // insert KeySets into Subgraphs
    13221323  DoLog(0) && (Log() << Verbose(0) << "Inserting KeySets into Subgraph class." << endl);
    1323   for (int j=0;j<mol->AtomCount;j++) {
     1324  for (int j=0;j<mol->getAtomCount();j++) {
    13241325    Subgraphs.insert(GraphPair (TestSets[j],pair<int, double>(counter++, 1.)));
    13251326  }
    13261327  DoLog(0) && (Log() << Verbose(0) << "Testing insertion of already present item in Subgraph." << endl);
    13271328  GraphTestPair test2;
    1328   test2 = Subgraphs.insert(GraphPair (TestSets[mol->AtomCount],pair<int, double>(counter++, 1.)));
     1329  test2 = Subgraphs.insert(GraphPair (TestSets[mol->getAtomCount()],pair<int, double>(counter++, 1.)));
    13291330  if (test2.second) {
    13301331    DoLog(1) && (Log() << Verbose(1) << "Insertion worked?!" << endl);
    … …  
    14941495 */
    14951496static int ParseCommandLineOptions(int argc, char **argv, MoleculeListClass *&molecules, periodentafel *&periode,
    1496                                    config& configuration, char *&ConfigFileName, set<int> &ArgcList)
     1497                                   config& configuration, char **ConfigFileName, set<int> &ArgcList)
    14971498{
    14981499  Vector x,y,z,n;  // coordinates for absolute point in cell volume
    … …  
    15121513  molecule *mol = NULL;
    15131514  string BondGraphFileName("\n");
    1514   strncpy(configuration.databasepath, LocalPath, MAXSTRINGSIZE-1);
     1515  bool DatabasePathGiven = false;
    15151516
    15161517  if (argc > 1) { // config file specified as option
    … …  
    15291530            break;
    15301531          case 'v':
    1531             ArgcList.insert(argptr-1);
    1532             return(1);
    1533             break;
    1534           case 'V':
     1532            setVerbosity(atoi(argv[argptr]));
    15351533            ArgcList.insert(argptr-1);
    15361534            ArgcList.insert(argptr);
    15371535            argptr++;
    15381536            break;
     1537          case 'V':
     1538            ArgcList.insert(argptr-1);
     1539            return(1);
     1540            break;
    15391541          case 'B':
     1542            ArgcList.insert(argptr-1);
     1543            ArgcList.insert(argptr);
     1544            ArgcList.insert(argptr+1);
     1545            ArgcList.insert(argptr+2);
     1546            ArgcList.insert(argptr+3);
     1547            ArgcList.insert(argptr+4);
     1548            ArgcList.insert(argptr+5);
     1549            argptr+=6;
    15401550            if (ExitFlag == 0) ExitFlag = 1;
    1541             if ((argptr+5 >= argc) || (argv[argptr][0] == '-') || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1])) || (!IsValidNumber(argv[argptr+2])) || (!IsValidNumber(argv[argptr+3])) || (!IsValidNumber(argv[argptr+4])) || (!IsValidNumber(argv[argptr+5])) ) {
    1542               ExitFlag = 255;
    1543               DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for bounding in box: -B <xx> <xy> <xz> <yy> <yz> <zz>" << endl);
    1544               performCriticalExit();
    1545             } else {
    1546               SaveFlag = true;
    1547               j = -1;
    1548               DoLog(1) && (Log() << Verbose(1) << "Centering atoms in config file within given simulation box." << endl);
    1549               double * const cell_size = World::getInstance().getDomain();
    1550               for (int i=0;i<6;i++) {
    1551                 cell_size[i] = atof(argv[argptr+i]);
    1552               }
    1553               argptr+=6;
    1554             }
    15551551            break;
    15561552          case 'e':
    … …  
    15611557              DoLog(0) && (Log() << Verbose(0) << "Using " << argv[argptr] << " as elements database." << endl);
    15621558              strncpy (configuration.databasepath, argv[argptr], MAXSTRINGSIZE-1);
     1559              DatabasePathGiven = true;
    15631560              argptr+=1;
    15641561            }
    … …  
    15741571            }
    15751572            break;
     1573          case 'M':
     1574            if ((argptr >= argc) || (argv[argptr][0] == '-')) {
     1575              ExitFlag = 255;
     1576              DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for setting MPQC basis: -M <basis name>" << endl);
     1577              performCriticalExit();
     1578            } else {
     1579              configuration.basis = argv[argptr];
     1580              DoLog(1) && (Log() << Verbose(1) << "Setting MPQC basis to " << configuration.basis << "." << endl);
     1581              argptr+=1;
     1582            }
     1583            break;
    15761584          case 'n':
    15771585            DoLog(0) && (Log() << Verbose(0) << "I won't parse trajectories." << endl);
    … …  
    15811589            {
    15821590              World::getInstance().setDefaultName(argv[argptr]);
    1583               DoLog(0) && (Log() << Verbose(0) << "Default name of new molecules set to " << *World::getInstance().getDefaultName() << "." << endl);
     1591              DoLog(0) && (Log() << Verbose(0) << "Default name of new molecules set to " << World::getInstance().getDefaultName() << "." << endl);
    15841592            }
    15851593            break;
    … …  
    15931601
    15941602    // 3a. Parse the element database
    1595     if (periode->LoadPeriodentafel(configuration.databasepath)) {
    1596       DoLog(0) && (Log() << Verbose(0) << "Element list loaded successfully." << endl);
    1597       //periode->Output();
    1598     } else {
    1599       DoLog(0) && (Log() << Verbose(0) << "Element list loading failed." << endl);
    1600       return 1;
    1601     }
     1603    if (DatabasePathGiven)
     1604      if (periode->LoadPeriodentafel(configuration.databasepath)) {
     1605        DoLog(0) && (Log() << Verbose(0) << "Element list loaded successfully." << endl);
     1606        //periode->Output();
     1607      } else {
     1608        DoLog(0) && (Log() << Verbose(0) << "Element list loading failed." << endl);
     1609        return 1;
     1610      }
    16021611    // 3b. Find config file name and parse if possible, also BondGraphFileName
    16031612    if (argv[1][0] != '-') {
    … …  
    16131622        } else {
    16141623          DoLog(0) && (Log() << Verbose(0) << "Empty configuration file." << endl);
    1615           ConfigFileName = argv[1];
     1624          strcpy(*ConfigFileName, argv[1]);
    16161625          configPresent = empty;
    16171626          output.close();
    … …  
    16191628      } else {
    16201629        test.close();
    1621         ConfigFileName = argv[1];
     1630        strcpy(*ConfigFileName, argv[1]);
    16221631        DoLog(1) && (Log() << Verbose(1) << "Specified config file found, parsing ... ");
    1623         switch (configuration.TestSyntax(ConfigFileName, periode)) {
     1632        switch (configuration.TestSyntax(*ConfigFileName, periode)) {
    16241633          case 1:
    16251634            DoLog(0) && (Log() << Verbose(0) << "new syntax." << endl);
    1626             configuration.Load(ConfigFileName, BondGraphFileName, periode, molecules);
     1635            configuration.Load(*ConfigFileName, BondGraphFileName, periode, molecules);
    16271636            configPresent = present;
    16281637            break;
    16291638          case 0:
    16301639            DoLog(0) && (Log() << Verbose(0) << "old syntax." << endl);
    1631             configuration.LoadOld(ConfigFileName, BondGraphFileName, periode, molecules);
     1640            configuration.LoadOld(*ConfigFileName, BondGraphFileName, periode, molecules);
    16321641            configPresent = present;
    16331642            break;
    … …  
    16501659       mol = World::getInstance().createMolecule();
    16511660       mol->ActiveFlag = true;
    1652        if (ConfigFileName != NULL)
    1653          mol->SetNameFromFilename(ConfigFileName);
     1661       if (*ConfigFileName != NULL)
     1662         mol->SetNameFromFilename(*ConfigFileName);
    16541663       molecules->insert(mol);
    16551664     }
    … …  
    17051714                if (first->type != NULL) {
    17061715                  mol->AddAtom(first);  // add to molecule
    1707                   if ((configPresent == empty) && (mol->AtomCount != 0))
     1716                  if ((configPresent == empty) && (mol->getAtomCount() != 0))
    17081717                    configPresent = present;
    17091718                } else
    … …  
    17181727        if (configPresent == present) {
    17191728          switch(argv[argptr-1][1]) {
    1720             case 'M':
    1721               if ((argptr >= argc) || (argv[argptr][0] == '-')) {
    1722                 ExitFlag = 255;
    1723                 DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for setting MPQC basis: -B <basis name>" << endl);
    1724                 performCriticalExit();
    1725               } else {
    1726                 configuration.basis = argv[argptr];
    1727                 DoLog(1) && (Log() << Verbose(1) << "Setting MPQC basis to " << configuration.basis << "." << endl);
    1728                 argptr+=1;
    1729               }
    1730               break;
    17311729            case 'D':
    17321730              if (ExitFlag == 0) ExitFlag = 1;
    … …  
    17341732                DoLog(1) && (Log() << Verbose(1) << "Depth-First-Search Analysis." << endl);
    17351733                MoleculeLeafClass *Subgraphs = NULL;      // list of subgraphs from DFS analysis
    1736                 int *MinimumRingSize = new int[mol->AtomCount];
     1734                int *MinimumRingSize = new int[mol->getAtomCount()];
    17371735                atom ***ListOfLocalAtoms = NULL;
    17381736                class StackClass<bond *> *BackEdgeStack = NULL;
    … …  
    17541752                  delete(Subgraphs);
    17551753                  for (int i=0;i<FragmentCounter;i++)
    1756                     Free(&ListOfLocalAtoms[i]);
    1757                   Free(&ListOfLocalAtoms);
     1754                    delete[](ListOfLocalAtoms[i]);
     1755                  delete[](ListOfLocalAtoms);
    17581756                }
    17591757                delete(BackEdgeStack);
    … …  
    17951793                        if ((argptr+6 >= argc) || (!IsValidNumber(argv[argptr+1])) || (!IsValidNumber(argv[argptr+5])) || (!IsValidNumber(argv[argptr+6])) || (!IsValidNumber(argv[argptr+2])) || (argv[argptr+1][0] == '-') || (argv[argptr+2][0] == '-') || (argv[argptr+3][0] == '-') || (argv[argptr+4][0] == '-')) {
    17961794                          ExitFlag = 255;
    1797                           DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for pair correlation analysis: -C E <Z1> <Z2> <output> <bin output>" << endl);
     1795                          DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for pair correlation analysis: -C E <Z1> <Z2> <output> <bin output> <binstart> <binend>" << endl);
    17981796                          performCriticalExit();
    17991797                        } else {
    … …  
    18101808                          else
    18111809                            correlationmap = PairCorrelation(molecules, elemental, elemental2);
    1812                           //OutputCorrelationToSurface(&output, correlationmap);
     1810                          OutputPairCorrelation(&output, correlationmap);
    18131811                          BinPairMap *binmap = BinData( correlationmap, 0.5, BinStart, BinEnd );
    18141812                          OutputCorrelation ( &binoutput, binmap );
    … …  
    18261824                        if ((argptr+8 >= argc) || (!IsValidNumber(argv[argptr+1])) ||  (!IsValidNumber(argv[argptr+2])) || (!IsValidNumber(argv[argptr+3])) || (!IsValidNumber(argv[argptr+4])) || (!IsValidNumber(argv[argptr+7])) || (!IsValidNumber(argv[argptr+8])) || (argv[argptr+1][0] == '-') || (argv[argptr+2][0] == '-') || (argv[argptr+3][0] == '-') || (argv[argptr+4][0] == '-') || (argv[argptr+5][0] == '-') || (argv[argptr+6][0] == '-')) {
    18271825                          ExitFlag = 255;
    1828                           DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for pair correlation analysis: -C P <Z1> <x> <y> <z> <output> <bin output>" << endl);
     1826                          DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for pair correlation analysis: -C P <Z1> <x> <y> <z> <output> <bin output> <binstart> <binend>" << endl);
    18291827                          performCriticalExit();
    18301828                        } else {
    … …  
    18411839                          else
    18421840                            correlationmap = CorrelationToPoint(molecules, elemental, Point);
    1843                           //OutputCorrelationToSurface(&output, correlationmap);
     1841                          OutputCorrelationToPoint(&output, correlationmap);
    18441842                          BinPairMap *binmap = BinData( correlationmap, 0.5, BinStart, BinEnd );
    18451843                          OutputCorrelation ( &binoutput, binmap );
    … …  
    18821880                          int counter  = 0;
    18831881                          for (MoleculeList::iterator BigFinder = molecules->ListOfMolecules.begin(); BigFinder != molecules->ListOfMolecules.end(); BigFinder++) {
    1884                             if ((Boundary == NULL) || (Boundary->AtomCount < (*BigFinder)->AtomCount)) {
     1882                            if ((Boundary == NULL) || (Boundary->getAtomCount() < (*BigFinder)->getAtomCount())) {
    18851883                              Boundary = *BigFinder;
    18861884                            }
    18871885                            counter++;
    18881886                          }
    1889                           bool *Actives = Malloc<bool>(counter, "ParseCommandLineOptions() - case C -- *Actives");
     1887                          bool *Actives = new bool[counter];
    18901888                          counter = 0;
    18911889                          for (MoleculeList::iterator BigFinder = molecules->ListOfMolecules.begin(); BigFinder != molecules->ListOfMolecules.end(); BigFinder++) {
    … …  
    19131911                          output.close();
    19141912                          binoutput.close();
     1913                          counter = 0;
    19151914                          for (MoleculeList::iterator BigFinder = molecules->ListOfMolecules.begin(); BigFinder != molecules->ListOfMolecules.end(); BigFinder++)
    19161915                            (*BigFinder)->ActiveFlag = Actives[counter++];
    1917                           Free(&Actives);
     1916                          delete[](Actives);
    19181917                          delete(LCList);
    19191918                          delete(TesselStruct);
    … …  
    19411940                performCriticalExit();
    19421941              } else {
     1942                mol->getAtomCount();
    19431943                SaveFlag = true;
    19441944                DoLog(1) && (Log() << Verbose(1) << "Changing atom " << argv[argptr] << " to element " << argv[argptr+1] << "." << endl);
    … …  
    19501950            case 'F':
    19511951              if (ExitFlag == 0) ExitFlag = 1;
    1952               MaxDistance = -1;
    1953               if (argv[argptr-1][2] == 'F') { // option is -FF?
    1954                 // fetch first argument as max distance to surface
    1955                 MaxDistance = atof(argv[argptr++]);
    1956                 DoLog(0) && (Log() << Verbose(0) << "Filling with maximum layer distance of " << MaxDistance << "." << endl);
    1957               }
    1958               if ((argptr+7 >=argc) || (argv[argptr][0] == '-') || (!IsValidNumber(argv[argptr+1])) || (!IsValidNumber(argv[argptr+2])) || (!IsValidNumber(argv[argptr+3])) || (!IsValidNumber(argv[argptr+4])) || (!IsValidNumber(argv[argptr+5])) || (!IsValidNumber(argv[argptr+6])) || (!IsValidNumber(argv[argptr+7]))) {
    1959                 ExitFlag = 255;
    1960                 DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for filling box with water: -F <xyz of filler> <dist_x> <dist_y> <dist_z> <boundary> <randatom> <randmol> <DoRotate>" << endl);
    1961                 performCriticalExit();
    1962               } else {
    1963                 SaveFlag = true;
    1964                 DoLog(1) && (Log() << Verbose(1) << "Filling Box with water molecules." << endl);
    1965                 // construct water molecule
    1966                 molecule *filler = World::getInstance().createMolecule();
    1967                 if (!filler->AddXYZFile(argv[argptr])) {
    1968                   DoeLog(0) && (eLog()<< Verbose(0) << "Could not parse filler molecule from " << argv[argptr] << "." << endl);
    1969                 }
    1970                 filler->SetNameFromFilename(argv[argptr]);
    1971                 configuration.BG->ConstructBondGraph(filler);
    1972                 molecule *Filling = NULL;
    1973                 atom *second = NULL, *third = NULL;
    1974                 first = World::getInstance().createAtom();
    1975                 first->type = periode->FindElement(1);
    1976                 first->x = Vector(0.441, -0.143, 0.);
    1977                 filler->AddAtom(first);
    1978                 second = World::getInstance().createAtom();
    1979                 second->type = periode->FindElement(1);
    1980                 second->x = Vector(-0.464, 1.137, 0.0);
    1981                 filler->AddAtom(second);
    1982                 third = World::getInstance().createAtom();
    1983                 third->type = periode->FindElement(8);
    1984                 third->x = Vector(-0.464, 0.177, 0.);
    1985                 filler->AddAtom(third);
    1986                 filler->AddBond(first, third, 1);
    1987                 filler->AddBond(second, third, 1);
    1988                 // call routine
    1989                 double distance[NDIM];
    1990                 for (int i=0;i<NDIM;i++)
    1991                   distance[i] = atof(argv[argptr+i+1]);
    1992                 Filling = FillBoxWithMolecule(molecules, filler, configuration, MaxDistance, distance, atof(argv[argptr+4]), atof(argv[argptr+5]), atof(argv[argptr+6]), atoi(argv[argptr+7]));
    1993                 if (Filling != NULL) {
    1994                   Filling->ActiveFlag = false;
    1995                   molecules->insert(Filling);
    1996                 }
    1997                 World::getInstance().destroyMolecule(filler);
    1998                 argptr+=6;
    1999               }
     1952              ArgcList.insert(argptr-1);
     1953              ArgcList.insert(argptr);
     1954              ArgcList.insert(argptr+1);
     1955              ArgcList.insert(argptr+2);
     1956              ArgcList.insert(argptr+3);
     1957              ArgcList.insert(argptr+4);
     1958              ArgcList.insert(argptr+5);
     1959              ArgcList.insert(argptr+6);
     1960              ArgcList.insert(argptr+7);
     1961              ArgcList.insert(argptr+8);
     1962              ArgcList.insert(argptr+9);
     1963              ArgcList.insert(argptr+10);
     1964              ArgcList.insert(argptr+11);
     1965              ArgcList.insert(argptr+12);
     1966              argptr+=13;
    20001967              break;
    20011968            case 'A':
    … …  
    20071974              } else {
    20081975                DoLog(0) && (Log() << Verbose(0) << "Parsing bonds from " << argv[argptr] << "." << endl);
    2009                 ifstream *input = new ifstream(argv[argptr]);
    2010                 mol->CreateAdjacencyListFromDbondFile(input);
    2011                 input->close();
     1976                ifstream input(argv[argptr]);
     1977                mol->CreateAdjacencyListFromDbondFile(&input);
     1978                input.close();
    20121979                argptr+=1;
    20131980              }
    … …  
    20592026                int counter  = 0;
    20602027                for (MoleculeList::iterator BigFinder = molecules->ListOfMolecules.begin(); BigFinder != molecules->ListOfMolecules.end(); BigFinder++) {
    2061                   (*BigFinder)->CountAtoms();
    2062                   if ((Boundary == NULL) || (Boundary->AtomCount < (*BigFinder)->AtomCount)) {
     2028                  if ((Boundary == NULL) || (Boundary->getAtomCount() < (*BigFinder)->getAtomCount())) {
    20632029                    Boundary = *BigFinder;
    20642030                  }
    20652031                  counter++;
    20662032                }
    2067                 DoLog(1) && (Log() << Verbose(1) << "Biggest molecule has " << Boundary->AtomCount << " atoms." << endl);
     2033                DoLog(1) && (Log() << Verbose(1) << "Biggest molecule has " << Boundary->getAtomCount() << " atoms." << endl);
    20682034                start = clock();
    20692035                LCList = new LinkedCell(Boundary, atof(argv[argptr])*2.);
    … …  
    21002066              if ((argptr >= argc) || (argv[argptr][0] == '-')) {
    21012067                ExitFlag = 255;
    2102                 DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for storing tempature: -L <step0> <step1> <prefix> <identity mapping?>" << endl);
     2068                DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for linear interpolation: -L <step0> <step1> <prefix> <identity mapping?>" << endl);
    21032069                performCriticalExit();
    21042070              } else {
    … …  
    21322098            case 'R':
    21332099              if (ExitFlag == 0) ExitFlag = 1;
    2134               if ((argptr+1 >= argc) || (argv[argptr][0] == '-') || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1])))  {
     2100              if ((argptr+3 >= argc) || (argv[argptr][0] == '-') || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1])) || (!IsValidNumber(argv[argptr+2])) || (!IsValidNumber(argv[argptr+3])))  {
    21352101                ExitFlag = 255;
    2136                 DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for removing atoms: -R <id> <distance>" << endl);
     2102                DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for removing atoms: -R <x> <y> <z> <distance>" << endl);
    21372103                performCriticalExit();
    21382104              } else {
    21392105                SaveFlag = true;
    2140                 DoLog(1) && (Log() << Verbose(1) << "Removing atoms around " << argv[argptr] << " with radius " << argv[argptr+1] << "." << endl);
    2141                 double tmp1 = atof(argv[argptr+1]);
    2142                 atom *third = mol->FindAtom(atoi(argv[argptr]));
    2143                 atom *first = mol->start;
    2144                 if ((third != NULL) && (first != mol->end)) {
    2145                   atom *second = first->next;
    2146                   while(second != mol->end) {
    2147                     first = second;
    2148                     second = first->next;
    2149                     if (first->x.DistanceSquared(third->x) > tmp1*tmp1) // distance to first above radius ...
    2150                       mol->RemoveAtom(first);
     2106                const double radius = atof(argv[argptr+3]);
     2107                Vector point(atof(argv[argptr]),atof(argv[argptr+1]),atof(argv[argptr+2]));
     2108                DoLog(1) && (Log() << Verbose(1) << "Removing atoms around " << point << " with radius " << radius << "." << endl);
     2109                atom *Walker = NULL;
     2110                molecule::iterator advancer = mol->begin();
     2111                for(molecule::iterator iter = advancer; advancer != mol->end();) {
     2112                  iter = advancer++;
     2113                  if ((*iter)->x.DistanceSquared(point) > radius*radius){ // distance to first above radius ...
     2114                    Walker = (*iter);
     2115                    DoLog(1) && (Log() << Verbose(1) << "Removing atom " << *Walker << "." << endl);
     2116                    mol->RemoveAtom(*(iter));
     2117                    World::getInstance().destroyAtom(Walker);
    21512118                  }
    2152                 } else {
    2153                   DoeLog(1) && (eLog()<< Verbose(1) << "Removal failed due to missing atoms on molecule or wrong id." << endl);
    21542119                }
    2155                 argptr+=2;
     2120                argptr+=4;
    21562121              }
    21572122              break;
    … …  
    21902155            case 's':
    21912156              if (ExitFlag == 0) ExitFlag = 1;
    2192               if ((argptr >= argc) || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1])) || (!IsValidNumber(argv[argptr+2])) ) {
     2157              if ((argptr+2 >= argc) || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1])) || (!IsValidNumber(argv[argptr+2])) ) {
    21932158                ExitFlag = 255;
    21942159                DoeLog(0) && (eLog()<< Verbose(0) << "Not enough or invalid arguments given for scaling: -s <factor_x> [factor_y] [factor_z]" << endl);
    21952160                performCriticalExit();
    21962161              } else {
    2197                 SaveFlag = true;
    2198                 j = -1;
    2199                 DoLog(1) && (Log() << Verbose(1) << "Scaling all ion positions by factor." << endl);
    2200                 factor = new double[NDIM];
    2201                 factor[0] = atof(argv[argptr]);
    2202                 factor[1] = atof(argv[argptr+1]);
    2203                 factor[2] = atof(argv[argptr+2]);
    2204                 mol->Scale((const double ** const)&factor);
    2205                 double * const cell_size = World::getInstance().getDomain();
    2206                 for (int i=0;i<NDIM;i++) {
    2207                   j += i+1;
    2208                   x[i] = atof(argv[NDIM+i]);
    2209                   cell_size[j]*=factor[i];
    2210                 }
    2211                 delete[](factor);
     2162                ArgcList.insert(argptr-1);
     2163                ArgcList.insert(argptr);
     2164                ArgcList.insert(argptr+1);
     2165                ArgcList.insert(argptr+2);
    22122166                argptr+=3;
    22132167              }
    … …  
    22202174                performCriticalExit();
    22212175              } else {
    2222                 SaveFlag = true;
    2223                 j = -1;
    2224                 DoLog(1) && (Log() << Verbose(1) << "Centering atoms in config file within given simulation box." << endl);
    2225                 double * const cell_size = World::getInstance().getDomain();
    2226                 for (int i=0;i<6;i++) {
    2227                   cell_size[i] = atof(argv[argptr+i]);
    2228                 }
    2229                 // center
    2230                 mol->CenterInBox();
     2176                ArgcList.insert(argptr-1);
     2177                ArgcList.insert(argptr);
     2178                ArgcList.insert(argptr+1);
     2179                ArgcList.insert(argptr+2);
     2180                ArgcList.insert(argptr+3);
     2181                ArgcList.insert(argptr+4);
     2182                ArgcList.insert(argptr+5);
    22312183                argptr+=6;
    22322184              }
    … …  
    22392191                performCriticalExit();
    22402192              } else {
    2241                 SaveFlag = true;
    2242                 j = -1;
    2243                 DoLog(1) && (Log() << Verbose(1) << "Centering atoms in config file within given simulation box." << endl);
    2244                 double * const cell_size = World::getInstance().getDomain();
    2245                 for (int i=0;i<6;i++) {
    2246                   cell_size[i] = atof(argv[argptr+i]);
    2247                 }
    2248                 // center
    2249                 mol->BoundInBox();
     2193                ArgcList.insert(argptr-1);
     2194                ArgcList.insert(argptr);
     2195                ArgcList.insert(argptr+1);
     2196                ArgcList.insert(argptr+2);
     2197                ArgcList.insert(argptr+3);
     2198                ArgcList.insert(argptr+4);
     2199                ArgcList.insert(argptr+5);
    22502200                argptr+=6;
    22512201              }
    … …  
    22582208                performCriticalExit();
    22592209              } else {
    2260                 SaveFlag = true;
    2261                 j = -1;
    2262                 DoLog(1) && (Log() << Verbose(1) << "Centering atoms in config file within given additional boundary." << endl);
    2263                 // make every coordinate positive
    2264                 mol->CenterEdge(&x);
    2265                 // update Box of atoms by boundary
    2266                 mol->SetBoxDimension(&x);
    2267                 // translate each coordinate by boundary
    2268                 double * const cell_size = World::getInstance().getDomain();
    2269                 j=-1;
    2270                 for (int i=0;i<NDIM;i++) {
    2271                   j += i+1;
    2272                   x[i] = atof(argv[argptr+i]);
    2273                   cell_size[j] += x[i]*2.;
    2274                 }
    2275                 mol->Translate((const Vector *)&x);
     2210                ArgcList.insert(argptr-1);
     2211                ArgcList.insert(argptr);
     2212                ArgcList.insert(argptr+1);
     2213                ArgcList.insert(argptr+2);
    22762214                argptr+=3;
    22772215              }
    … …  
    22792217            case 'O':
    22802218              if (ExitFlag == 0) ExitFlag = 1;
    2281               SaveFlag = true;
    2282               DoLog(1) && (Log() << Verbose(1) << "Centering atoms on edge and setting box dimensions." << endl);
    2283               x.Zero();
    2284               mol->CenterEdge(&x);
    2285               mol->SetBoxDimension(&x);
     2219              ArgcList.insert(argptr-1);
    22862220              argptr+=0;
    22872221              break;
    … …  
    22932227                performCriticalExit();
    22942228              } else {
     2229                mol->getAtomCount();
    22952230                SaveFlag = true;
    22962231                DoLog(1) && (Log() << Verbose(1) << "Removing atom " << argv[argptr] << "." << endl);
    … …  
    23122247                mol->CreateAdjacencyList(atof(argv[argptr]), configuration.GetIsAngstroem(), &BondGraph::CovalentMinMaxDistance, NULL);
    23132248                DoLog(0) && (Log() << Verbose(0) << "Fragmenting molecule with current connection matrix ..." << endl);
    2314                 if (mol->first->next != mol->last) {
     2249                if (mol->hasBondStructure()) {
    23152250                  ExitFlag = mol->FragmentMolecule(atoi(argv[argptr+1]), &configuration);
    23162251                }
    … …  
    24012336                performCriticalExit();
    24022337              } else {
    2403                 SaveFlag = true;
    2404                 double * const cell_size = World::getInstance().getDomain();
    2405                 for (int axis = 1; axis <= NDIM; axis++) {
    2406                   int faktor = atoi(argv[argptr++]);
    2407                   int count;
    2408                   const element ** Elements;
    2409                   Vector ** vectors;
    2410                   if (faktor < 1) {
    2411                     DoeLog(1) && (eLog()<< Verbose(1) << "Repetition factor mus be greater than 1!" << endl);
    2412                     faktor = 1;
    2413                   }
    2414                   mol->CountAtoms();  // recount atoms
    2415                   if (mol->AtomCount != 0) {  // if there is more than none
    2416                     count = mol->AtomCount;   // is changed becausing of adding, thus has to be stored away beforehand
    2417                     Elements = new const element *[count];
    2418                     vectors = new Vector *[count];
    2419                     j = 0;
    2420                     first = mol->start;
    2421                     while (first->next != mol->end) {  // make a list of all atoms with coordinates and element
    2422                       first = first->next;
    2423                       Elements[j] = first->type;
    2424                       vectors[j] = &first->x;
    2425                       j++;
    2426                     }
    2427                     if (count != j)
    2428                       DoeLog(1) && (eLog()<< Verbose(1) << "AtomCount " << count << " is not equal to number of atoms in molecule " << j << "!" << endl);
    2429                     x.Zero();
    2430                     y.Zero();
    2431                     y[abs(axis)-1] = cell_size[(abs(axis) == 2) ? 2 : ((abs(axis) == 3) ? 5 : 0)] * abs(axis)/axis; // last term is for sign, first is for magnitude
    2432                     for (int i=1;i<faktor;i++) {  // then add this list with respective translation factor times
    2433                       x += y; // per factor one cell width further
    2434                       for (int k=count;k--;) { // go through every atom of the original cell
    2435                         first = World::getInstance().createAtom(); // create a new body
    2436                         first->x = (*vectors[k]) + x;
    2437                         first->type = Elements[k];  // insert original element
    2438                         mol->AddAtom(first);        // and add to the molecule (which increments ElementsInMolecule, AtomCount, ...)
    2439                       }
    2440                     }
    2441                     // free memory
    2442                     delete[](Elements);
    2443                     delete[](vectors);
    2444                     // correct cell size
    2445                     if (axis < 0) { // if sign was negative, we have to translate everything
    2446                       x =(-(faktor-1)) * y;
    2447                       mol->Translate(&x);
    2448                     }
    2449                     cell_size[(abs(axis) == 2) ? 2 : ((abs(axis) == 3) ? 5 : 0)] *= faktor;
    2450                   }
    2451                 }
     2338                ArgcList.insert(argptr-1);
     2339                ArgcList.insert(argptr);
     2340                ArgcList.insert(argptr+1);
     2341                ArgcList.insert(argptr+2);
     2342                argptr+=3;
    24522343              }
    24532344              break;
    … …  
    24612352    } while (argptr < argc);
    24622353    if (SaveFlag)
    2463       configuration.SaveAll(ConfigFileName, periode, molecules);
     2354      configuration.SaveAll(*ConfigFileName, periode, molecules);
    24642355  } else {  // no arguments, hence scan the elements db
    24652356    if (periode->LoadPeriodentafel(configuration.databasepath))
    … …  
    24762367void cleanUp(){
    24772368  World::purgeInstance();
    2478   Log() << Verbose(0) <<  "Maximum of allocated memory: "
    2479     << MemoryUsageObserver::getInstance()->getMaximumUsedMemory() << endl;
    2480   Log() << Verbose(0) <<  "Remaining non-freed memory: "
    2481     << MemoryUsageObserver::getInstance()->getUsedMemorySize() << endl;
    2482   MemoryUsageObserver::purgeInstance();
    24832369  logger::purgeInstance();
    24842370  errorLogger::purgeInstance();
    … …  
    24882374  ActionRegistry::purgeInstance();
    24892375  ActionHistory::purgeInstance();
     2376  Memory::getState();
    24902377}
    24912378
    … …  
    24932380{
    24942381    config *configuration = World::getInstance().getConfig();
     2382    // while we are non interactive, we want to abort from asserts
     2383    //ASSERT_DO(Assert::Abort);
     2384    molecule *mol = NULL;
    24952385    Vector x, y, z, n;
    24962386    ifstream test;
    … …  
    24992389    char **Arguments = NULL;
    25002390    int ArgcSize = 0;
     2391    int ExitFlag = 0;
    25012392    bool ArgumentsCopied = false;
    2502 
     2393    char *ConfigFileName = new char[MAXSTRINGSIZE];
     2394
     2395    // print version check whether arguments are present at all
    25032396    cout << ESPACKVersion << endl;
     2397    if (argc < 2) {
     2398      cout << "Obtain help with " << argv[0] << " -h." << endl;
     2399      cleanUp();
     2400      Memory::getState();
     2401      return(1);
     2402    }
     2403
    25042404
    25052405    setVerbosity(0);
    25062406    // need to init the history before any action is created
    25072407    ActionHistory::init();
     2408
     2409    // In the interactive mode, we can leave the user the choice in case of error
     2410    ASSERT_DO(Assert::Ask);
     2411
     2412    // from this moment on, we need to be sure to deeinitialize in the correct order
     2413    // this is handled by the cleanup function
     2414    atexit(cleanUp);
    25082415
    25092416    // Parse command line options and if present create respective UI
    … …  
    25122419      ArgcList.insert(0); // push back program!
    25132420      ArgcList.insert(1); // push back config file name
    2514       char ConfigFileName[MAXSTRINGSIZE];
    25152421      // handle arguments by ParseCommandLineOptions()
    2516       ParseCommandLineOptions(argc,argv,World::getInstance().getMolecules(),World::getInstance().getPeriode(),*World::getInstance().getConfig(), (char *&)ConfigFileName, ArgcList);
     2422      ExitFlag = ParseCommandLineOptions(argc,argv,World::getInstance().getMolecules(),World::getInstance().getPeriode(),*World::getInstance().getConfig(), &ConfigFileName, ArgcList);
    25172423      // copy all remaining arguments to a new argv
    2518       Arguments = Malloc<char *>(ArgcList.size(), "main - **Arguments");
     2424      Arguments = new char *[ArgcList.size()];
    25192425      cout << "The following arguments are handled by CommandLineParser: ";
    25202426      for (set<int>::iterator ArgcRunner = ArgcList.begin(); ArgcRunner != ArgcList.end(); ++ArgcRunner) {
    2521         Arguments[ArgcSize] = Malloc<char>(strlen(argv[*ArgcRunner])+2, "main - *Arguments[]");
     2427        Arguments[ArgcSize] = new char[strlen(argv[*ArgcRunner])+2];
    25222428        strcpy(Arguments[ArgcSize], argv[*ArgcRunner]);
    25232429        cout << " " << argv[*ArgcRunner];
    … …  
    25442450    }
    25452451
    2546     if(World::getInstance().getPeriode()->StorePeriodentafel(configuration->databasepath))
    2547         Log() << Verbose(0) << "Saving of elements.db successful." << endl;
    2548 
    2549     else
    2550         Log() << Verbose(0) << "Saving of elements.db failed." << endl;
     2452    Log() << Verbose(0) << "Saving to " << ConfigFileName << "." << endl;
     2453    World::getInstance().getConfig()->SaveAll(ConfigFileName, World::getInstance().getPeriode(), World::getInstance().getMolecules());
    25512454
    25522455  // free the new argv
    25532456  if (ArgumentsCopied) {
    25542457    for (int i=0; i<ArgcSize;i++)
    2555       Free(&Arguments[i]);
    2556     Free(&Arguments);
     2458      delete[](Arguments[i]);
     2459    delete[](Arguments);
    25572460  }
    2558 
    2559   cleanUp();
    2560   Memory::getState();
    2561   return (0);
     2461  delete[](ConfigFileName);
     2462
     2463  return (ExitFlag == 1 ? 0 : ExitFlag);
    25622464}
    25632465
Note: See TracChangeset for help on using the changeset viewer.