Changeset 042f82 for src/builder.cpp


Ignore:
Timestamp:
Jul 23, 2009, 2:21:57 PM (15 years ago)
Author:
Frederik Heber <heber@…>
Branches:
Action_Thermostats, Add_AtomRandomPerturbation, Add_FitFragmentPartialChargesAction, Add_RotateAroundBondAction, Add_SelectAtomByNameAction, Added_ParseSaveFragmentResults, AddingActions_SaveParseParticleParameters, Adding_Graph_to_ChangeBondActions, Adding_MD_integration_tests, Adding_ParticleName_to_Atom, Adding_StructOpt_integration_tests, AtomFragments, Automaking_mpqc_open, AutomationFragmentation_failures, Candidate_v1.5.4, Candidate_v1.6.0, Candidate_v1.6.1, ChangeBugEmailaddress, ChangingTestPorts, ChemicalSpaceEvaluator, CombiningParticlePotentialParsing, Combining_Subpackages, Debian_Package_split, Debian_package_split_molecuildergui_only, Disabling_MemDebug, Docu_Python_wait, EmpiricalPotential_contain_HomologyGraph, EmpiricalPotential_contain_HomologyGraph_documentation, Enable_parallel_make_install, Enhance_userguide, Enhanced_StructuralOptimization, Enhanced_StructuralOptimization_continued, Example_ManyWaysToTranslateAtom, Exclude_Hydrogens_annealWithBondGraph, FitPartialCharges_GlobalError, Fix_BoundInBox_CenterInBox_MoleculeActions, Fix_ChargeSampling_PBC, Fix_ChronosMutex, Fix_FitPartialCharges, Fix_FitPotential_needs_atomicnumbers, Fix_ForceAnnealing, Fix_IndependentFragmentGrids, Fix_ParseParticles, Fix_ParseParticles_split_forward_backward_Actions, Fix_PopActions, Fix_QtFragmentList_sorted_selection, Fix_Restrictedkeyset_FragmentMolecule, Fix_StatusMsg, Fix_StepWorldTime_single_argument, Fix_Verbose_Codepatterns, Fix_fitting_potentials, Fixes, ForceAnnealing_goodresults, ForceAnnealing_oldresults, ForceAnnealing_tocheck, ForceAnnealing_with_BondGraph, ForceAnnealing_with_BondGraph_continued, ForceAnnealing_with_BondGraph_continued_betteresults, ForceAnnealing_with_BondGraph_contraction-expansion, FragmentAction_writes_AtomFragments, FragmentMolecule_checks_bonddegrees, GeometryObjects, Gui_Fixes, Gui_displays_atomic_force_velocity, ImplicitCharges, IndependentFragmentGrids, IndependentFragmentGrids_IndividualZeroInstances, IndependentFragmentGrids_IntegrationTest, IndependentFragmentGrids_Sole_NN_Calculation, JobMarket_RobustOnKillsSegFaults, JobMarket_StableWorkerPool, JobMarket_unresolvable_hostname_fix, MoreRobust_FragmentAutomation, ODR_violation_mpqc_open, PartialCharges_OrthogonalSummation, PdbParser_setsAtomName, PythonUI_with_named_parameters, QtGui_reactivate_TimeChanged_changes, Recreated_GuiChecks, Rewrite_FitPartialCharges, RotateToPrincipalAxisSystem_UndoRedo, SaturateAtoms_findBestMatching, SaturateAtoms_singleDegree, StoppableMakroAction, Subpackage_CodePatterns, Subpackage_JobMarket, Subpackage_LinearAlgebra, Subpackage_levmar, Subpackage_mpqc_open, Subpackage_vmg, Switchable_LogView, ThirdParty_MPQC_rebuilt_buildsystem, TrajectoryDependenant_MaxOrder, TremoloParser_IncreasedPrecision, TremoloParser_MultipleTimesteps, TremoloParser_setsAtomName, Ubuntu_1604_changes, stable
Children:
36ec71
Parents:
205ccd
Message:

fixed indentation from tabs to two spaces.

File:
1 edited

Legend:

Unmodified
Added
Removed

  • src/builder.cpp

    r205ccd r042f82  
    1515 * \section about About the Program
    1616 *
    17  *      Molecuilder is a short program, written in C++, that enables the construction of a coordinate set for the
    18  *      atoms making up an molecule by the successive statement of binding angles and distances and referencing to
    19  *      already constructed atoms.
     17 *  Molecuilder is a short program, written in C++, that enables the construction of a coordinate set for the
     18 *  atoms making up an molecule by the successive statement of binding angles and distances and referencing to
     19 *  already constructed atoms.
    2020 *
    21  *      A configuration file may be written that is compatible to the format used by PCP - a parallel Car-Parrinello
    22  *      molecular dynamics implementation.
     21 *  A configuration file may be written that is compatible to the format used by PCP - a parallel Car-Parrinello
     22 *  molecular dynamics implementation.
    2323 *
    2424 * \section install Installation
    2525 *
    26  *      Installation should without problems succeed as follows:
    27  *      -# ./configure (or: mkdir build;mkdir run;cd build; ../configure --bindir=../run)
    28  *      -# make
    29  *      -# make install
     26 *  Installation should without problems succeed as follows:
     27 *  -# ./configure (or: mkdir build;mkdir run;cd build; ../configure --bindir=../run)
     28 *  -# make
     29 *  -# make install
    3030 *
    31  *      Further useful commands are
    32  *      -# make clean uninstall: deletes .o-files and removes executable from the given binary directory\n
    33  *      -# make doxygen-doc: Creates these html pages out of the documented source
     31 *  Further useful commands are
     32 *  -# make clean uninstall: deletes .o-files and removes executable from the given binary directory\n
     33 *  -# make doxygen-doc: Creates these html pages out of the documented source
    3434 *
    3535 * \section run Running
    3636 *
    37  *      The program can be executed by running: ./molecuilder
     37 *  The program can be executed by running: ./molecuilder
    3838 *
    39  *      Note, that it uses a database, called "elements.db", in the executable's directory. If the file is not found,
    40  *      it is created and any given data on elements of the periodic table will be stored therein and re-used on
    41  *      later re-execution.
     39 *  Note, that it uses a database, called "elements.db", in the executable's directory. If the file is not found,
     40 *  it is created and any given data on elements of the periodic table will be stored therein and re-used on
     41 *  later re-execution.
    4242 *
    4343 * \section ref References
    4444 *
    45  *      For the special configuration file format, see the documentation of pcp.
     45 *  For the special configuration file format, see the documentation of pcp.
    4646 *
    4747 */
    … …  
    6363static void AddAtoms(periodentafel *periode, molecule *mol)
    6464{
    65         atom *first, *second, *third, *fourth;
    66         Vector **atoms;
    67         Vector x,y,z,n; // coordinates for absolute point in cell volume
    68         double a,b,c;
    69         char choice;    // menu choice char
    70         bool valid;
    71 
    72         cout << Verbose(0) << "===========ADD ATOM============================" << endl;
    73         cout << Verbose(0) << " a - state absolute coordinates of atom" << endl;
    74         cout << Verbose(0) << " b - state relative coordinates of atom wrt to reference point" << endl;
    75         cout << Verbose(0) << " c - state relative coordinates of atom wrt to already placed atom" << endl;
    76         cout << Verbose(0) << " d - state two atoms, two angles and a distance" << endl;
    77         cout << Verbose(0) << " e - least square distance position to a set of atoms" << endl;
    78         cout << Verbose(0) << "all else - go back" << endl;
    79         cout << Verbose(0) << "===============================================" << endl;
    80         cout << Verbose(0) << "Note: Specifiy angles in degrees not multiples of Pi!" << endl;
    81         cout << Verbose(0) << "INPUT: ";
    82         cin >> choice;
    83 
    84         switch (choice) {
     65  atom *first, *second, *third, *fourth;
     66  Vector **atoms;
     67  Vector x,y,z,n;  // coordinates for absolute point in cell volume
     68  double a,b,c;
     69  char choice;  // menu choice char
     70  bool valid;
     71
     72  cout << Verbose(0) << "===========ADD ATOM============================" << endl;
     73  cout << Verbose(0) << " a - state absolute coordinates of atom" << endl;
     74  cout << Verbose(0) << " b - state relative coordinates of atom wrt to reference point" << endl;
     75  cout << Verbose(0) << " c - state relative coordinates of atom wrt to already placed atom" << endl;
     76  cout << Verbose(0) << " d - state two atoms, two angles and a distance" << endl;
     77  cout << Verbose(0) << " e - least square distance position to a set of atoms" << endl;
     78  cout << Verbose(0) << "all else - go back" << endl;
     79  cout << Verbose(0) << "===============================================" << endl;
     80  cout << Verbose(0) << "Note: Specifiy angles in degrees not multiples of Pi!" << endl;
     81  cout << Verbose(0) << "INPUT: ";
     82  cin >> choice;
     83
     84  switch (choice) {
    8585    default:
    8686      cout << Verbose(0) << "Not a valid choice." << endl;
    8787      break;
    88                         case 'a': // absolute coordinates of atom
     88      case 'a': // absolute coordinates of atom
    8989        cout << Verbose(0) << "Enter absolute coordinates." << endl;
    9090        first = new atom;
    9191        first->x.AskPosition(mol->cell_size, false);
    92         first->type = periode->AskElement();    // give type
    93         mol->AddAtom(first);    // add to molecule
    94                                 break;
    95 
    96                         case 'b': // relative coordinates of atom wrt to reference point
     92        first->type = periode->AskElement();  // give type
     93        mol->AddAtom(first);  // add to molecule
     94        break;
     95
     96      case 'b': // relative coordinates of atom wrt to reference point
    9797        first = new atom;
    9898        valid = true;
    … …  
    106106          cout << Verbose(0) << "\n";
    107107        } while (!(valid = mol->CheckBounds((const Vector *)&first->x)));
    108         first->type = periode->AskElement();    // give type
    109         mol->AddAtom(first);    // add to molecule
    110                                 break;
    111 
    112                         case 'c': // relative coordinates of atom wrt to already placed atom
     108        first->type = periode->AskElement();  // give type
     109        mol->AddAtom(first);  // add to molecule
     110        break;
     111
     112      case 'c': // relative coordinates of atom wrt to already placed atom
    113113        first = new atom;
    114114        valid = true;
    … …  
    122122          }
    123123        } while (!(valid = mol->CheckBounds((const Vector *)&first->x)));
    124         first->type = periode->AskElement();    // give type
    125         mol->AddAtom(first);    // add to molecule
     124        first->type = periode->AskElement();  // give type
     125        mol->AddAtom(first);  // add to molecule
    126126        break;
    127127
    … …  
    224224          cout << Verbose(0) << endl;
    225225        } while (!(valid = mol->CheckBounds((const Vector *)&first->x)));
    226         first->type = periode->AskElement();    // give type
    227         mol->AddAtom(first);    // add to molecule
    228                                 break;
    229 
    230                         case 'e': // least square distance position to a set of atoms
     226        first->type = periode->AskElement();  // give type
     227        mol->AddAtom(first);  // add to molecule
     228        break;
     229
     230      case 'e': // least square distance position to a set of atoms
    231231        first = new atom;
    232232        atoms = new (Vector*[128]);
    … …  
    248248
    249249          first->x.Output((ofstream *)&cout);
    250           first->type = periode->AskElement();  // give type
    251           mol->AddAtom(first);  // add to molecule
     250          first->type = periode->AskElement();  // give type
     251          mol->AddAtom(first);  // add to molecule
    252252        } else {
    253253          delete first;
    254254          cout << Verbose(0) << "Please enter at least two vectors!\n";
    255255        }
    256                                 break;
    257         };
     256        break;
     257  };
    258258};
    259259
    … …  
    263263static void CenterAtoms(molecule *mol)
    264264{
    265         Vector x, y, helper;
    266         char choice;    // menu choice char
    267 
    268         cout << Verbose(0) << "===========CENTER ATOMS=========================" << endl;
    269         cout << Verbose(0) << " a - on origin" << endl;
    270         cout << Verbose(0) << " b - on center of gravity" << endl;
    271         cout << Verbose(0) << " c - within box with additional boundary" << endl;
    272         cout << Verbose(0) << " d - within given simulation box" << endl;
    273         cout << Verbose(0) << "all else - go back" << endl;
    274         cout << Verbose(0) << "===============================================" << endl;
    275         cout << Verbose(0) << "INPUT: ";
    276         cin >> choice;
    277 
    278         switch (choice) {
    279                 default:
    280                         cout << Verbose(0) << "Not a valid choice." << endl;
    281                         break;
    282                 case 'a':
    283                         cout << Verbose(0) << "Centering atoms in config file on origin." << endl;
    284                         mol->CenterOrigin((ofstream *)&cout, &x);
    285                         break;
    286                 case 'b':
    287                         cout << Verbose(0) << "Centering atoms in config file on center of gravity." << endl;
    288                         mol->CenterGravity((ofstream *)&cout, &x);
    289                         break;
    290                 case 'c':
    291                         cout << Verbose(0) << "Centering atoms in config file within given additional boundary." << endl;
    292                         for (int i=0;i<NDIM;i++) {
    293                                 cout << Verbose(0) << "Enter axis " << i << " boundary: ";
    294                                 cin >> y.x[i];
    295                         }
    296                         mol->CenterEdge((ofstream *)&cout, &x); // make every coordinate positive
    297                         mol->Translate(&y); // translate by boundary
    298                         helper.CopyVector(&y);
    299                         helper.Scale(2.);
    300                         helper.AddVector(&x);
    301                         mol->SetBoxDimension(&helper);  // update Box of atoms by boundary
    302                         break;
    303                 case 'd':
    304                         cout << Verbose(1) << "Centering atoms in config file within given simulation box." << endl;
    305                         for (int i=0;i<NDIM;i++) {
    306                                 cout << Verbose(0) << "Enter axis " << i << " boundary: ";
    307                                 cin >> x.x[i];
    308                         }
    309                         // center
    310                         mol->CenterInBox((ofstream *)&cout);
    311                         // update Box of atoms by boundary
    312                         mol->SetBoxDimension(&x);
    313                         break;
    314         }
     265  Vector x, y, helper;
     266  char choice;  // menu choice char
     267
     268  cout << Verbose(0) << "===========CENTER ATOMS=========================" << endl;
     269  cout << Verbose(0) << " a - on origin" << endl;
     270  cout << Verbose(0) << " b - on center of gravity" << endl;
     271  cout << Verbose(0) << " c - within box with additional boundary" << endl;
     272  cout << Verbose(0) << " d - within given simulation box" << endl;
     273  cout << Verbose(0) << "all else - go back" << endl;
     274  cout << Verbose(0) << "===============================================" << endl;
     275  cout << Verbose(0) << "INPUT: ";
     276  cin >> choice;
     277
     278  switch (choice) {
     279    default:
     280      cout << Verbose(0) << "Not a valid choice." << endl;
     281      break;
     282    case 'a':
     283      cout << Verbose(0) << "Centering atoms in config file on origin." << endl;
     284      mol->CenterOrigin((ofstream *)&cout, &x);
     285      break;
     286    case 'b':
     287      cout << Verbose(0) << "Centering atoms in config file on center of gravity." << endl;
     288      mol->CenterGravity((ofstream *)&cout, &x);
     289      break;
     290    case 'c':
     291      cout << Verbose(0) << "Centering atoms in config file within given additional boundary." << endl;
     292      for (int i=0;i<NDIM;i++) {
     293        cout << Verbose(0) << "Enter axis " << i << " boundary: ";
     294        cin >> y.x[i];
     295      }
     296      mol->CenterEdge((ofstream *)&cout, &x);  // make every coordinate positive
     297      mol->Translate(&y); // translate by boundary
     298      helper.CopyVector(&y);
     299      helper.Scale(2.);
     300      helper.AddVector(&x);
     301      mol->SetBoxDimension(&helper);  // update Box of atoms by boundary
     302      break;
     303    case 'd':
     304      cout << Verbose(1) << "Centering atoms in config file within given simulation box." << endl;
     305      for (int i=0;i<NDIM;i++) {
     306        cout << Verbose(0) << "Enter axis " << i << " boundary: ";
     307        cin >> x.x[i];
     308      }
     309      // center
     310      mol->CenterInBox((ofstream *)&cout);
     311      // update Box of atoms by boundary
     312      mol->SetBoxDimension(&x);
     313      break;
     314  }
    315315};
    316316
    … …  
    321321static void AlignAtoms(periodentafel *periode, molecule *mol)
    322322{
    323         atom *first, *second, *third;
    324         Vector x,n;
    325         char choice;    // menu choice char
    326 
    327         cout << Verbose(0) << "===========ALIGN ATOMS=========================" << endl;
    328         cout << Verbose(0) << " a - state three atoms defining align plane" << endl;
    329         cout << Verbose(0) << " b - state alignment vector" << endl;
    330         cout << Verbose(0) << " c - state two atoms in alignment direction" << endl;
    331         cout << Verbose(0) << " d - align automatically by least square fit" << endl;
    332         cout << Verbose(0) << "all else - go back" << endl;
    333         cout << Verbose(0) << "===============================================" << endl;
    334         cout << Verbose(0) << "INPUT: ";
    335         cin >> choice;
    336 
    337         switch (choice) {
    338                 default:
    339                 case 'a': // three atoms defining mirror plane
    340                         first = mol->AskAtom("Enter first atom: ");
    341                         second = mol->AskAtom("Enter second atom: ");
    342                         third = mol->AskAtom("Enter third atom: ");
    343 
    344                         n.MakeNormalVector((const Vector *)&first->x,(const Vector *)&second->x,(const Vector *)&third->x);
    345                         break;
    346                 case 'b': // normal vector of mirror plane
    347                         cout << Verbose(0) << "Enter normal vector of mirror plane." << endl;
    348                         n.AskPosition(mol->cell_size,0);
    349                         n.Normalize();
    350                         break;
    351                 case 'c': // three atoms defining mirror plane
    352                         first = mol->AskAtom("Enter first atom: ");
    353                         second = mol->AskAtom("Enter second atom: ");
    354 
    355                         n.CopyVector((const Vector *)&first->x);
    356                         n.SubtractVector((const Vector *)&second->x);
    357                         n.Normalize();
    358                         break;
    359                 case 'd':
    360                         char shorthand[4];
    361                         Vector a;
    362                         struct lsq_params param;
    363                         do {
    364                                 fprintf(stdout, "Enter the element of atoms to be chosen: ");
    365                                 fscanf(stdin, "%3s", shorthand);
    366                         } while ((param.type = periode->FindElement(shorthand)) == NULL);
    367                         cout << Verbose(0) << "Element is " << param.type->name << endl;
    368                         mol->GetAlignvector(&param);
    369                         for (int i=NDIM;i--;) {
    370                                 x.x[i] = gsl_vector_get(param.x,i);
    371                                 n.x[i] = gsl_vector_get(param.x,i+NDIM);
    372                         }
    373                         gsl_vector_free(param.x);
    374                         cout << Verbose(0) << "Offset vector: ";
    375                         x.Output((ofstream *)&cout);
    376                         cout << Verbose(0) << endl;
    377                         n.Normalize();
    378                         break;
    379         };
    380         cout << Verbose(0) << "Alignment vector: ";
    381         n.Output((ofstream *)&cout);
    382         cout << Verbose(0) << endl;
    383         mol->Align(&n);
     323  atom *first, *second, *third;
     324  Vector x,n;
     325  char choice;  // menu choice char
     326
     327  cout << Verbose(0) << "===========ALIGN ATOMS=========================" << endl;
     328  cout << Verbose(0) << " a - state three atoms defining align plane" << endl;
     329  cout << Verbose(0) << " b - state alignment vector" << endl;
     330  cout << Verbose(0) << " c - state two atoms in alignment direction" << endl;
     331  cout << Verbose(0) << " d - align automatically by least square fit" << endl;
     332  cout << Verbose(0) << "all else - go back" << endl;
     333  cout << Verbose(0) << "===============================================" << endl;
     334  cout << Verbose(0) << "INPUT: ";
     335  cin >> choice;
     336
     337  switch (choice) {
     338    default:
     339    case 'a': // three atoms defining mirror plane
     340      first = mol->AskAtom("Enter first atom: ");
     341      second = mol->AskAtom("Enter second atom: ");
     342      third = mol->AskAtom("Enter third atom: ");
     343
     344      n.MakeNormalVector((const Vector *)&first->x,(const Vector *)&second->x,(const Vector *)&third->x);
     345      break;
     346    case 'b': // normal vector of mirror plane
     347      cout << Verbose(0) << "Enter normal vector of mirror plane." << endl;
     348      n.AskPosition(mol->cell_size,0);
     349      n.Normalize();
     350      break;
     351    case 'c': // three atoms defining mirror plane
     352      first = mol->AskAtom("Enter first atom: ");
     353      second = mol->AskAtom("Enter second atom: ");
     354
     355      n.CopyVector((const Vector *)&first->x);
     356      n.SubtractVector((const Vector *)&second->x);
     357      n.Normalize();
     358      break;
     359    case 'd':
     360      char shorthand[4];
     361      Vector a;
     362      struct lsq_params param;
     363      do {
     364        fprintf(stdout, "Enter the element of atoms to be chosen: ");
     365        fscanf(stdin, "%3s", shorthand);
     366      } while ((param.type = periode->FindElement(shorthand)) == NULL);
     367      cout << Verbose(0) << "Element is " << param.type->name << endl;
     368      mol->GetAlignvector(&param);
     369      for (int i=NDIM;i--;) {
     370        x.x[i] = gsl_vector_get(param.x,i);
     371        n.x[i] = gsl_vector_get(param.x,i+NDIM);
     372      }
     373      gsl_vector_free(param.x);
     374      cout << Verbose(0) << "Offset vector: ";
     375      x.Output((ofstream *)&cout);
     376      cout << Verbose(0) << endl;
     377      n.Normalize();
     378      break;
     379  };
     380  cout << Verbose(0) << "Alignment vector: ";
     381  n.Output((ofstream *)&cout);
     382  cout << Verbose(0) << endl;
     383  mol->Align(&n);
    384384};
    385385
    … …  
    389389static void MirrorAtoms(molecule *mol)
    390390{
    391         atom *first, *second, *third;
    392         Vector n;
    393         char choice;    // menu choice char
    394 
    395         cout << Verbose(0) << "===========MIRROR ATOMS=========================" << endl;
    396         cout << Verbose(0) << " a - state three atoms defining mirror plane" << endl;
    397         cout << Verbose(0) << " b - state normal vector of mirror plane" << endl;
    398         cout << Verbose(0) << " c - state two atoms in normal direction" << endl;
    399         cout << Verbose(0) << "all else - go back" << endl;
    400         cout << Verbose(0) << "===============================================" << endl;
    401         cout << Verbose(0) << "INPUT: ";
    402         cin >> choice;
    403 
    404         switch (choice) {
    405                 default:
    406                 case 'a': // three atoms defining mirror plane
    407                         first = mol->AskAtom("Enter first atom: ");
    408                         second = mol->AskAtom("Enter second atom: ");
    409                         third = mol->AskAtom("Enter third atom: ");
    410 
    411                         n.MakeNormalVector((const Vector *)&first->x,(const Vector *)&second->x,(const Vector *)&third->x);
    412                         break;
    413                 case 'b': // normal vector of mirror plane
    414                         cout << Verbose(0) << "Enter normal vector of mirror plane." << endl;
    415                         n.AskPosition(mol->cell_size,0);
    416                         n.Normalize();
    417                         break;
    418                 case 'c': // three atoms defining mirror plane
    419                         first = mol->AskAtom("Enter first atom: ");
    420                         second = mol->AskAtom("Enter second atom: ");
    421 
    422                         n.CopyVector((const Vector *)&first->x);
    423                         n.SubtractVector((const Vector *)&second->x);
    424                         n.Normalize();
    425                         break;
    426         };
    427         cout << Verbose(0) << "Normal vector: ";
    428         n.Output((ofstream *)&cout);
    429         cout << Verbose(0) << endl;
    430         mol->Mirror((const Vector *)&n);
     391  atom *first, *second, *third;
     392  Vector n;
     393  char choice;  // menu choice char
     394
     395  cout << Verbose(0) << "===========MIRROR ATOMS=========================" << endl;
     396  cout << Verbose(0) << " a - state three atoms defining mirror plane" << endl;
     397  cout << Verbose(0) << " b - state normal vector of mirror plane" << endl;
     398  cout << Verbose(0) << " c - state two atoms in normal direction" << endl;
     399  cout << Verbose(0) << "all else - go back" << endl;
     400  cout << Verbose(0) << "===============================================" << endl;
     401  cout << Verbose(0) << "INPUT: ";
     402  cin >> choice;
     403
     404  switch (choice) {
     405    default:
     406    case 'a': // three atoms defining mirror plane
     407      first = mol->AskAtom("Enter first atom: ");
     408      second = mol->AskAtom("Enter second atom: ");
     409      third = mol->AskAtom("Enter third atom: ");
     410
     411      n.MakeNormalVector((const Vector *)&first->x,(const Vector *)&second->x,(const Vector *)&third->x);
     412      break;
     413    case 'b': // normal vector of mirror plane
     414      cout << Verbose(0) << "Enter normal vector of mirror plane." << endl;
     415      n.AskPosition(mol->cell_size,0);
     416      n.Normalize();
     417      break;
     418    case 'c': // three atoms defining mirror plane
     419      first = mol->AskAtom("Enter first atom: ");
     420      second = mol->AskAtom("Enter second atom: ");
     421
     422      n.CopyVector((const Vector *)&first->x);
     423      n.SubtractVector((const Vector *)&second->x);
     424      n.Normalize();
     425      break;
     426  };
     427  cout << Verbose(0) << "Normal vector: ";
     428  n.Output((ofstream *)&cout);
     429  cout << Verbose(0) << endl;
     430  mol->Mirror((const Vector *)&n);
    431431};
    432432
    … …  
    436436static void RemoveAtoms(molecule *mol)
    437437{
    438         atom *first, *second;
    439         int axis;
    440         double tmp1, tmp2;
    441         char choice;    // menu choice char
    442 
    443         cout << Verbose(0) << "===========REMOVE ATOMS=========================" << endl;
    444         cout << Verbose(0) << " a - state atom for removal by number" << endl;
    445         cout << Verbose(0) << " b - keep only in radius around atom" << endl;
    446         cout << Verbose(0) << " c - remove this with one axis greater value" << endl;
    447         cout << Verbose(0) << "all else - go back" << endl;
    448         cout << Verbose(0) << "===============================================" << endl;
    449         cout << Verbose(0) << "INPUT: ";
    450         cin >> choice;
    451 
    452         switch (choice) {
    453                 default:
    454                 case 'a':
    455                         if (mol->RemoveAtom(mol->AskAtom("Enter number of atom within molecule: ")))
    456                                 cout << Verbose(1) << "Atom removed." << endl;
    457                         else
    458                                 cout << Verbose(1) << "Atom not found." << endl;
    459                         break;
    460                 case 'b':
    461                         second = mol->AskAtom("Enter number of atom as reference point: ");
    462                         cout << Verbose(0) << "Enter radius: ";
    463                         cin >> tmp1;
    464                         first = mol->start;
    465                         while(first->next != mol->end) {
    466                                 first = first->next;
    467                                 if (first->x.DistanceSquared((const Vector *)&second->x) > tmp1*tmp1) // distance to first above radius ...
    468                                         mol->RemoveAtom(first);
    469                         }
    470                         break;
    471                 case 'c':
    472                         cout << Verbose(0) << "Which axis is it: ";
    473                         cin >> axis;
    474                         cout << Verbose(0) << "Left inward boundary: ";
    475                         cin >> tmp1;
    476                         cout << Verbose(0) << "Right inward boundary: ";
    477                         cin >> tmp2;
    478                         first = mol->start;
    479                         while(first->next != mol->end) {
    480                                 first = first->next;
    481                                 if ((first->x.x[axis] > tmp2) || (first->x.x[axis] < tmp1)) // out of boundary ...
    482                                         mol->RemoveAtom(first);
    483                         }
    484                         break;
    485         };
    486         //mol->Output((ofstream *)&cout);
    487         choice = 'r';
     438  atom *first, *second;
     439  int axis;
     440  double tmp1, tmp2;
     441  char choice;  // menu choice char
     442
     443  cout << Verbose(0) << "===========REMOVE ATOMS=========================" << endl;
     444  cout << Verbose(0) << " a - state atom for removal by number" << endl;
     445  cout << Verbose(0) << " b - keep only in radius around atom" << endl;
     446  cout << Verbose(0) << " c - remove this with one axis greater value" << endl;
     447  cout << Verbose(0) << "all else - go back" << endl;
     448  cout << Verbose(0) << "===============================================" << endl;
     449  cout << Verbose(0) << "INPUT: ";
     450  cin >> choice;
     451
     452  switch (choice) {
     453    default:
     454    case 'a':
     455      if (mol->RemoveAtom(mol->AskAtom("Enter number of atom within molecule: ")))
     456        cout << Verbose(1) << "Atom removed." << endl;
     457      else
     458        cout << Verbose(1) << "Atom not found." << endl;
     459      break;
     460    case 'b':
     461      second = mol->AskAtom("Enter number of atom as reference point: ");
     462      cout << Verbose(0) << "Enter radius: ";
     463      cin >> tmp1;
     464      first = mol->start;
     465      while(first->next != mol->end) {
     466        first = first->next;
     467        if (first->x.DistanceSquared((const Vector *)&second->x) > tmp1*tmp1) // distance to first above radius ...
     468          mol->RemoveAtom(first);
     469      }
     470      break;
     471    case 'c':
     472      cout << Verbose(0) << "Which axis is it: ";
     473      cin >> axis;
     474      cout << Verbose(0) << "Left inward boundary: ";
     475      cin >> tmp1;
     476      cout << Verbose(0) << "Right inward boundary: ";
     477      cin >> tmp2;
     478      first = mol->start;
     479      while(first->next != mol->end) {
     480        first = first->next;
     481        if ((first->x.x[axis] > tmp2) || (first->x.x[axis] < tmp1)) // out of boundary ...
     482          mol->RemoveAtom(first);
     483      }
     484      break;
     485  };
     486  //mol->Output((ofstream *)&cout);
     487  choice = 'r';
    488488};
    489489
    … …  
    494494static void MeasureAtoms(periodentafel *periode, molecule *mol, config *configuration)
    495495{
    496         atom *first, *second, *third;
    497         Vector x,y;
    498         double min[256], tmp1, tmp2, tmp3;
    499         int Z;
    500         char choice;    // menu choice char
    501 
    502         cout << Verbose(0) << "===========MEASURE ATOMS=========================" << endl;
    503         cout << Verbose(0) << " a - calculate bond length between one atom and all others" << endl;
    504         cout << Verbose(0) << " b - calculate bond length between two atoms" << endl;
    505         cout << Verbose(0) << " c - calculate bond angle" << endl;
    506         cout << Verbose(0) << " d - calculate principal axis of the system" << endl;
    507         cout << Verbose(0) << " e - calculate volume of the convex envelope" << endl;
    508         cout << Verbose(0) << " f - calculate temperature from current velocity" << endl;
    509         cout << Verbose(0) << " g - output all temperatures per step from velocities" << endl;
    510         cout << Verbose(0) << "all else - go back" << endl;
    511         cout << Verbose(0) << "===============================================" << endl;
    512         cout << Verbose(0) << "INPUT: ";
    513         cin >> choice;
    514 
    515         switch(choice) {
    516                 default:
    517                         cout << Verbose(1) << "Not a valid choice." << endl;
    518                         break;
    519                 case 'a':
    520                         first = mol->AskAtom("Enter first atom: ");
    521                         for (int i=MAX_ELEMENTS;i--;)
    522                                 min[i] = 0.;
    523 
    524                         second = mol->start;
    525                         while ((second->next != mol->end)) {
    526                                 second = second->next; // advance
    527                                 Z = second->type->Z;
    528                                 tmp1 = 0.;
    529                                 if (first != second) {
    530                                         x.CopyVector((const Vector *)&first->x);
    531                                         x.SubtractVector((const Vector *)&second->x);
    532                                         tmp1 = x.Norm();
    533                                 }
    534                                 if ((tmp1 != 0.) && ((min[Z] == 0.) || (tmp1 < min[Z]))) min[Z] = tmp1;
    535                                 //cout << Verbose(0) << "Bond length between Atom " << first->nr << " and " << second->nr << ": " << tmp1 << " a.u." << endl;
    536                         }
    537                         for (int i=MAX_ELEMENTS;i--;)
    538                                 if (min[i] != 0.) cout << Verbose(0) << "Minimum Bond length between " << first->type->name << " Atom " << first->nr << " and next Ion of type " << (periode->FindElement(i))->name << ": " << min[i] << " a.u." << endl;
    539                         break;
    540 
    541                 case 'b':
    542                         first = mol->AskAtom("Enter first atom: ");
    543                         second = mol->AskAtom("Enter second atom: ");
    544                         for (int i=NDIM;i--;)
    545                                 min[i] = 0.;
    546                         x.CopyVector((const Vector *)&first->x);
    547                         x.SubtractVector((const Vector *)&second->x);
    548                         tmp1 = x.Norm();
    549                         cout << Verbose(1) << "Distance vector is ";
    550                         x.Output((ofstream *)&cout);
    551                         cout << "." << endl << "Norm of distance is " << tmp1 << "." << endl;
    552                         break;
    553 
    554                 case 'c':
    555                         cout << Verbose(0) << "Evaluating bond angle between three - first, central, last - atoms." << endl;
    556                         first = mol->AskAtom("Enter first atom: ");
    557                         second = mol->AskAtom("Enter central atom: ");
    558                         third   = mol->AskAtom("Enter last atom: ");
    559                         tmp1 = tmp2 = tmp3 = 0.;
    560                         x.CopyVector((const Vector *)&first->x);
    561                         x.SubtractVector((const Vector *)&second->x);
    562                         y.CopyVector((const Vector *)&third->x);
    563                         y.SubtractVector((const Vector *)&second->x);
    564                         cout << Verbose(0) << "Bond angle between first atom Nr." << first->nr << ", central atom Nr." << second->nr << " and last atom Nr." << third->nr << ": ";
    565                         cout << Verbose(0) << (acos(x.ScalarProduct((const Vector *)&y)/(y.Norm()*x.Norm()))/M_PI*180.) << " degrees" << endl;
    566                         break;
    567                 case 'd':
    568                         cout << Verbose(0) << "Evaluating prinicipal axis." << endl;
    569                         cout << Verbose(0) << "Shall we rotate? [0/1]: ";
    570                         cin >> Z;
    571                         if ((Z >=0) && (Z <=1))
    572                                 mol->PrincipalAxisSystem((ofstream *)&cout, (bool)Z);
    573                         else
    574                                 mol->PrincipalAxisSystem((ofstream *)&cout, false);
    575                         break;
    576                 case 'e':
    577                         cout << Verbose(0) << "Evaluating volume of the convex envelope.";
    578                         VolumeOfConvexEnvelope((ofstream *)&cout, NULL, configuration, NULL, mol);
    579                         break;
    580                 case 'f':
    581                         mol->OutputTemperatureFromTrajectories((ofstream *)&cout, mol->MDSteps-1, mol->MDSteps, (ofstream *)&cout);
    582                         break;
    583                 case 'g':
    584                         {
    585                                 char filename[255];
    586                                 cout << "Please enter filename: " << endl;
    587                                 cin >> filename;
    588                                 cout << Verbose(1) << "Storing temperatures in " << filename << "." << endl;
    589                                 ofstream *output = new ofstream(filename, ios::trunc);
    590                                 if (!mol->OutputTemperatureFromTrajectories((ofstream *)&cout, 0, mol->MDSteps, output))
    591                                         cout << Verbose(2) << "File could not be written." << endl;
    592                                 else
    593                                         cout << Verbose(2) << "File stored." << endl;
    594                                 output->close();
    595                                 delete(output);
    596                         }
    597                         break;
    598         }
     496  atom *first, *second, *third;
     497  Vector x,y;
     498  double min[256], tmp1, tmp2, tmp3;
     499  int Z;
     500  char choice;  // menu choice char
     501
     502  cout << Verbose(0) << "===========MEASURE ATOMS=========================" << endl;
     503  cout << Verbose(0) << " a - calculate bond length between one atom and all others" << endl;
     504  cout << Verbose(0) << " b - calculate bond length between two atoms" << endl;
     505  cout << Verbose(0) << " c - calculate bond angle" << endl;
     506  cout << Verbose(0) << " d - calculate principal axis of the system" << endl;
     507  cout << Verbose(0) << " e - calculate volume of the convex envelope" << endl;
     508  cout << Verbose(0) << " f - calculate temperature from current velocity" << endl;
     509  cout << Verbose(0) << " g - output all temperatures per step from velocities" << endl;
     510  cout << Verbose(0) << "all else - go back" << endl;
     511  cout << Verbose(0) << "===============================================" << endl;
     512  cout << Verbose(0) << "INPUT: ";
     513  cin >> choice;
     514
     515  switch(choice) {
     516    default:
     517      cout << Verbose(1) << "Not a valid choice." << endl;
     518      break;
     519    case 'a':
     520      first = mol->AskAtom("Enter first atom: ");
     521      for (int i=MAX_ELEMENTS;i--;)
     522        min[i] = 0.;
     523
     524      second = mol->start;
     525      while ((second->next != mol->end)) {
     526        second = second->next; // advance
     527        Z = second->type->Z;
     528        tmp1 = 0.;
     529        if (first != second) {
     530          x.CopyVector((const Vector *)&first->x);
     531          x.SubtractVector((const Vector *)&second->x);
     532          tmp1 = x.Norm();
     533        }
     534        if ((tmp1 != 0.) && ((min[Z] == 0.) || (tmp1 < min[Z]))) min[Z] = tmp1;
     535        //cout << Verbose(0) << "Bond length between Atom " << first->nr << " and " << second->nr << ": " << tmp1 << " a.u." << endl;
     536      }
     537      for (int i=MAX_ELEMENTS;i--;)
     538        if (min[i] != 0.) cout << Verbose(0) << "Minimum Bond length between " << first->type->name << " Atom " << first->nr << " and next Ion of type " << (periode->FindElement(i))->name << ": " << min[i] << " a.u." << endl;
     539      break;
     540
     541    case 'b':
     542      first = mol->AskAtom("Enter first atom: ");
     543      second = mol->AskAtom("Enter second atom: ");
     544      for (int i=NDIM;i--;)
     545        min[i] = 0.;
     546      x.CopyVector((const Vector *)&first->x);
     547      x.SubtractVector((const Vector *)&second->x);
     548      tmp1 = x.Norm();
     549      cout << Verbose(1) << "Distance vector is ";
     550      x.Output((ofstream *)&cout);
     551      cout << "." << endl << "Norm of distance is " << tmp1 << "." << endl;
     552      break;
     553
     554    case 'c':
     555      cout << Verbose(0) << "Evaluating bond angle between three - first, central, last - atoms." << endl;
     556      first = mol->AskAtom("Enter first atom: ");
     557      second = mol->AskAtom("Enter central atom: ");
     558      third  = mol->AskAtom("Enter last atom: ");
     559      tmp1 = tmp2 = tmp3 = 0.;
     560      x.CopyVector((const Vector *)&first->x);
     561      x.SubtractVector((const Vector *)&second->x);
     562      y.CopyVector((const Vector *)&third->x);
     563      y.SubtractVector((const Vector *)&second->x);
     564      cout << Verbose(0) << "Bond angle between first atom Nr." << first->nr << ", central atom Nr." << second->nr << " and last atom Nr." << third->nr << ": ";
     565      cout << Verbose(0) << (acos(x.ScalarProduct((const Vector *)&y)/(y.Norm()*x.Norm()))/M_PI*180.) << " degrees" << endl;
     566      break;
     567    case 'd':
     568      cout << Verbose(0) << "Evaluating prinicipal axis." << endl;
     569      cout << Verbose(0) << "Shall we rotate? [0/1]: ";
     570      cin >> Z;
     571      if ((Z >=0) && (Z <=1))
     572        mol->PrincipalAxisSystem((ofstream *)&cout, (bool)Z);
     573      else
     574        mol->PrincipalAxisSystem((ofstream *)&cout, false);
     575      break;
     576    case 'e':
     577      cout << Verbose(0) << "Evaluating volume of the convex envelope.";
     578      VolumeOfConvexEnvelope((ofstream *)&cout, NULL, configuration, NULL, mol);
     579      break;
     580    case 'f':
     581      mol->OutputTemperatureFromTrajectories((ofstream *)&cout, mol->MDSteps-1, mol->MDSteps, (ofstream *)&cout);
     582      break;
     583    case 'g':
     584      {
     585        char filename[255];
     586        cout << "Please enter filename: " << endl;
     587        cin >> filename;
     588        cout << Verbose(1) << "Storing temperatures in " << filename << "." << endl;
     589        ofstream *output = new ofstream(filename, ios::trunc);
     590        if (!mol->OutputTemperatureFromTrajectories((ofstream *)&cout, 0, mol->MDSteps, output))
     591          cout << Verbose(2) << "File could not be written." << endl;
     592        else
     593          cout << Verbose(2) << "File stored." << endl;
     594        output->close();
     595        delete(output);
     596      }
     597      break;
     598  }
    599599};
    600600
    … …  
    605605static void FragmentAtoms(molecule *mol, config *configuration)
    606606{
    607         int Order1;
    608         clock_t start, end;
    609 
    610         cout << Verbose(0) << "Fragmenting molecule with current connection matrix ..." << endl;
    611         cout << Verbose(0) << "What's the desired bond order: ";
    612         cin >> Order1;
    613         if (mol->first->next != mol->last) {    // there are bonds
    614                 start = clock();
    615                 mol->FragmentMolecule((ofstream *)&cout, Order1, configuration);
    616                 end = clock();
    617                 cout << Verbose(0) << "Clocks for this operation: " << (end-start) << ", time: " << ((double)(end-start)/CLOCKS_PER_SEC) << "s." << endl;
    618         } else
    619                 cout << Verbose(0) << "Connection matrix has not yet been generated!" << endl;
     607  int Order1;
     608  clock_t start, end;
     609
     610  cout << Verbose(0) << "Fragmenting molecule with current connection matrix ..." << endl;
     611  cout << Verbose(0) << "What's the desired bond order: ";
     612  cin >> Order1;
     613  if (mol->first->next != mol->last) {  // there are bonds
     614    start = clock();
     615    mol->FragmentMolecule((ofstream *)&cout, Order1, configuration);
     616    end = clock();
     617    cout << Verbose(0) << "Clocks for this operation: " << (end-start) << ", time: " << ((double)(end-start)/CLOCKS_PER_SEC) << "s." << endl;
     618  } else
     619    cout << Verbose(0) << "Connection matrix has not yet been generated!" << endl;
    620620};
    621621
    … …  
    10381038static void testroutine(MoleculeListClass *molecules)
    10391039{
    1040         // the current test routine checks the functionality of the KeySet&Graph concept:
    1041         // We want to have a multiindex (the KeySet) describing a unique subgraph
     1040  // the current test routine checks the functionality of the KeySet&Graph concept:
     1041  // We want to have a multiindex (the KeySet) describing a unique subgraph
    10421042  int i, comp, counter=0;
    10431043
    … …  
    10521052  atom *Walker = mol->start;
    10531053
    1054         // generate some KeySets
    1055         cout << "Generating KeySets." << endl;
    1056         KeySet TestSets[mol->AtomCount+1];
    1057         i=1;
    1058         while (Walker->next != mol->end) {
    1059                 Walker = Walker->next;
    1060                 for (int j=0;j<i;j++) {
    1061                         TestSets[j].insert(Walker->nr);
    1062                 }
    1063                 i++;
    1064         }
    1065         cout << "Testing insertion of already present item in KeySets." << endl;
    1066         KeySetTestPair test;
    1067         test = TestSets[mol->AtomCount-1].insert(Walker->nr);
    1068         if (test.second) {
    1069                 cout << Verbose(1) << "Insertion worked?!" << endl;
    1070         } else {
    1071                 cout << Verbose(1) << "Insertion rejected: Present object is " << (*test.first) << "." << endl;
    1072         }
    1073         TestSets[mol->AtomCount].insert(mol->end->previous->nr);
    1074         TestSets[mol->AtomCount].insert(mol->end->previous->previous->previous->nr);
    1075 
    1076         // constructing Graph structure
    1077         cout << "Generating Subgraph class." << endl;
    1078         Graph Subgraphs;
    1079 
    1080         // insert KeySets into Subgraphs
    1081         cout << "Inserting KeySets into Subgraph class." << endl;
    1082         for (int j=0;j<mol->AtomCount;j++) {
    1083                 Subgraphs.insert(GraphPair (TestSets[j],pair<int, double>(counter++, 1.)));
    1084         }
    1085         cout << "Testing insertion of already present item in Subgraph." << endl;
    1086         GraphTestPair test2;
    1087         test2 = Subgraphs.insert(GraphPair (TestSets[mol->AtomCount],pair<int, double>(counter++, 1.)));
    1088         if (test2.second) {
    1089                 cout << Verbose(1) << "Insertion worked?!" << endl;
    1090         } else {
    1091                 cout << Verbose(1) << "Insertion rejected: Present object is " << (*(test2.first)).second.first << "." << endl;
    1092         }
    1093 
    1094         // show graphs
    1095         cout << "Showing Subgraph's contents, checking that it's sorted." << endl;
    1096         Graph::iterator A = Subgraphs.begin();
    1097         while (A !=     Subgraphs.end()) {
    1098                 cout << (*A).second.first << ": ";
    1099                 KeySet::iterator key = (*A).first.begin();
    1100                 comp = -1;
    1101                 while (key != (*A).first.end()) {
    1102                         if ((*key) > comp)
    1103                                 cout << (*key) << " ";
    1104                         else
    1105                                 cout << (*key) << "! ";
    1106                         comp = (*key);
    1107                         key++;
    1108                 }
    1109                 cout << endl;
    1110                 A++;
    1111         }
    1112         delete(mol);
     1054  // generate some KeySets
     1055  cout << "Generating KeySets." << endl;
     1056  KeySet TestSets[mol->AtomCount+1];
     1057  i=1;
     1058  while (Walker->next != mol->end) {
     1059    Walker = Walker->next;
     1060    for (int j=0;j<i;j++) {
     1061      TestSets[j].insert(Walker->nr);
     1062    }
     1063    i++;
     1064  }
     1065  cout << "Testing insertion of already present item in KeySets." << endl;
     1066  KeySetTestPair test;
     1067  test = TestSets[mol->AtomCount-1].insert(Walker->nr);
     1068  if (test.second) {
     1069    cout << Verbose(1) << "Insertion worked?!" << endl;
     1070  } else {
     1071    cout << Verbose(1) << "Insertion rejected: Present object is " << (*test.first) << "." << endl;
     1072  }
     1073  TestSets[mol->AtomCount].insert(mol->end->previous->nr);
     1074  TestSets[mol->AtomCount].insert(mol->end->previous->previous->previous->nr);
     1075
     1076  // constructing Graph structure
     1077  cout << "Generating Subgraph class." << endl;
     1078  Graph Subgraphs;
     1079
     1080  // insert KeySets into Subgraphs
     1081  cout << "Inserting KeySets into Subgraph class." << endl;
     1082  for (int j=0;j<mol->AtomCount;j++) {
     1083    Subgraphs.insert(GraphPair (TestSets[j],pair<int, double>(counter++, 1.)));
     1084  }
     1085  cout << "Testing insertion of already present item in Subgraph." << endl;
     1086  GraphTestPair test2;
     1087  test2 = Subgraphs.insert(GraphPair (TestSets[mol->AtomCount],pair<int, double>(counter++, 1.)));
     1088  if (test2.second) {
     1089    cout << Verbose(1) << "Insertion worked?!" << endl;
     1090  } else {
     1091    cout << Verbose(1) << "Insertion rejected: Present object is " << (*(test2.first)).second.first << "." << endl;
     1092  }
     1093
     1094  // show graphs
     1095  cout << "Showing Subgraph's contents, checking that it's sorted." << endl;
     1096  Graph::iterator A = Subgraphs.begin();
     1097  while (A !=  Subgraphs.end()) {
     1098    cout << (*A).second.first << ": ";
     1099    KeySet::iterator key = (*A).first.begin();
     1100    comp = -1;
     1101    while (key != (*A).first.end()) {
     1102      if ((*key) > comp)
     1103        cout << (*key) << " ";
     1104      else
     1105        cout << (*key) << "! ";
     1106      comp = (*key);
     1107      key++;
     1108    }
     1109    cout << endl;
     1110    A++;
     1111  }
     1112  delete(mol);
    11131113};
    11141114
    … …  
    11211121static void SaveConfig(char *ConfigFileName, config *configuration, periodentafel *periode, MoleculeListClass *molecules)
    11221122{
    1123         char filename[MAXSTRINGSIZE];
    1124         ofstream output;
    1125         molecule *mol = new molecule(periode);
    1126 
    1127         // merge all molecules in MoleculeListClass into this molecule
    1128         int N = molecules->ListOfMolecules.size();
    1129         int *src = new int(N);
    1130         N=0;
    1131         for (MoleculeList::iterator ListRunner = molecules->ListOfMolecules.begin(); ListRunner != molecules->ListOfMolecules.end(); ListRunner++)
    1132           src[N++] = (*ListRunner)->IndexNr;
    1133         molecules->SimpleMultiAdd(mol, src, N);
    1134 
    1135         cout << Verbose(0) << "Storing configuration ... " << endl;
    1136         // get correct valence orbitals
    1137         mol->CalculateOrbitals(*configuration);
    1138         configuration->InitMaxMinStopStep = configuration->MaxMinStopStep = configuration->MaxPsiDouble;
    1139         strcpy(filename, ConfigFileName);
    1140         if (ConfigFileName != NULL) { // test the file name
    1141                 output.open(ConfigFileName, ios::trunc);
    1142         } else if (strlen(configuration->configname) != 0) {
    1143                 strcpy(filename, configuration->configname);
    1144                 output.open(configuration->configname, ios::trunc);
    1145                 } else {
    1146                         strcpy(filename, DEFAULTCONFIG);
    1147                         output.open(DEFAULTCONFIG, ios::trunc);
    1148                 }
    1149         output.close();
    1150         output.clear();
    1151         cout << Verbose(0) << "Saving of config file ";
    1152         if (configuration->Save(filename, periode, mol))
    1153                 cout << "successful." << endl;
    1154         else
    1155                 cout << "failed." << endl;
    1156 
    1157         // and save to xyz file
    1158         if (ConfigFileName != NULL) {
    1159                 strcpy(filename, ConfigFileName);
    1160                 strcat(filename, ".xyz");
    1161                 output.open(filename, ios::trunc);
    1162         }
    1163         if (output == NULL) {
    1164                 strcpy(filename,"main_pcp_linux");
    1165                 strcat(filename, ".xyz");
    1166                 output.open(filename, ios::trunc);
    1167         }
    1168         cout << Verbose(0) << "Saving of XYZ file ";
    1169         if (mol->MDSteps <= 1) {
    1170                 if (mol->OutputXYZ(&output))
    1171                         cout << "successful." << endl;
    1172                 else
    1173                         cout << "failed." << endl;
    1174         } else {
    1175                 if (mol->OutputTrajectoriesXYZ(&output))
    1176                         cout << "successful." << endl;
    1177                 else
    1178                         cout << "failed." << endl;
    1179         }
    1180         output.close();
    1181         output.clear();
    1182 
    1183         // and save as MPQC configuration
    1184         if (ConfigFileName != NULL)
    1185                 strcpy(filename, ConfigFileName);
    1186         if (output == NULL)
    1187                 strcpy(filename,"main_pcp_linux");
    1188         cout << Verbose(0) << "Saving as mpqc input ";
    1189         if (configuration->SaveMPQC(filename, mol))
    1190                 cout << "done." << endl;
    1191         else
    1192                 cout << "failed." << endl;
    1193 
    1194         if (!strcmp(configuration->configpath, configuration->GetDefaultPath())) {
    1195                 cerr << "WARNING: config is found under different path then stated in config file::defaultpath!" << endl;
    1196         }
    1197         delete(mol);
     1123  char filename[MAXSTRINGSIZE];
     1124  ofstream output;
     1125  molecule *mol = new molecule(periode);
     1126
     1127  // merge all molecules in MoleculeListClass into this molecule
     1128  int N = molecules->ListOfMolecules.size();
     1129  int *src = new int(N);
     1130  N=0;
     1131  for (MoleculeList::iterator ListRunner = molecules->ListOfMolecules.begin(); ListRunner != molecules->ListOfMolecules.end(); ListRunner++)
     1132    src[N++] = (*ListRunner)->IndexNr;
     1133  molecules->SimpleMultiAdd(mol, src, N);
     1134
     1135  cout << Verbose(0) << "Storing configuration ... " << endl;
     1136  // get correct valence orbitals
     1137  mol->CalculateOrbitals(*configuration);
     1138  configuration->InitMaxMinStopStep = configuration->MaxMinStopStep = configuration->MaxPsiDouble;
     1139  strcpy(filename, ConfigFileName);
     1140  if (ConfigFileName != NULL) { // test the file name
     1141    output.open(ConfigFileName, ios::trunc);
     1142  } else if (strlen(configuration->configname) != 0) {
     1143    strcpy(filename, configuration->configname);
     1144    output.open(configuration->configname, ios::trunc);
     1145    } else {
     1146      strcpy(filename, DEFAULTCONFIG);
     1147      output.open(DEFAULTCONFIG, ios::trunc);
     1148    }
     1149  output.close();
     1150  output.clear();
     1151  cout << Verbose(0) << "Saving of config file ";
     1152  if (configuration->Save(filename, periode, mol))
     1153    cout << "successful." << endl;
     1154  else
     1155    cout << "failed." << endl;
     1156
     1157  // and save to xyz file
     1158  if (ConfigFileName != NULL) {
     1159    strcpy(filename, ConfigFileName);
     1160    strcat(filename, ".xyz");
     1161    output.open(filename, ios::trunc);
     1162  }
     1163  if (output == NULL) {
     1164    strcpy(filename,"main_pcp_linux");
     1165    strcat(filename, ".xyz");
     1166    output.open(filename, ios::trunc);
     1167  }
     1168  cout << Verbose(0) << "Saving of XYZ file ";
     1169  if (mol->MDSteps <= 1) {
     1170    if (mol->OutputXYZ(&output))
     1171      cout << "successful." << endl;
     1172    else
     1173      cout << "failed." << endl;
     1174  } else {
     1175    if (mol->OutputTrajectoriesXYZ(&output))
     1176      cout << "successful." << endl;
     1177    else
     1178      cout << "failed." << endl;
     1179  }
     1180  output.close();
     1181  output.clear();
     1182
     1183  // and save as MPQC configuration
     1184  if (ConfigFileName != NULL)
     1185    strcpy(filename, ConfigFileName);
     1186  if (output == NULL)
     1187    strcpy(filename,"main_pcp_linux");
     1188  cout << Verbose(0) << "Saving as mpqc input ";
     1189  if (configuration->SaveMPQC(filename, mol))
     1190    cout << "done." << endl;
     1191  else
     1192    cout << "failed." << endl;
     1193
     1194  if (!strcmp(configuration->configpath, configuration->GetDefaultPath())) {
     1195    cerr << "WARNING: config is found under different path then stated in config file::defaultpath!" << endl;
     1196  }
     1197  delete(mol);
    11981198};
    11991199
    … …  
    12101210static int ParseCommandLineOptions(int argc, char **argv, MoleculeListClass *&molecules, periodentafel *&periode, config& configuration, char *&ConfigFileName)
    12111211{
    1212         Vector x,y,z,n; // coordinates for absolute point in cell volume
    1213         double *factor; // unit factor if desired
    1214         ifstream test;
    1215         ofstream output;
    1216         string line;
    1217         atom *first;
    1218         bool SaveFlag = false;
    1219         int ExitFlag = 0;
    1220         int j;
    1221         double volume = 0.;
    1222         enum ConfigStatus config_present = absent;
    1223         clock_t start,end;
    1224         int argptr;
     1212  Vector x,y,z,n;  // coordinates for absolute point in cell volume
     1213  double *factor; // unit factor if desired
     1214  ifstream test;
     1215  ofstream output;
     1216  string line;
     1217  atom *first;
     1218  bool SaveFlag = false;
     1219  int ExitFlag = 0;
     1220  int j;
     1221  double volume = 0.;
     1222  enum ConfigStatus config_present = absent;
     1223  clock_t start,end;
     1224  int argptr;
    12251225  strncpy(configuration.databasepath, LocalPath, MAXSTRINGSIZE-1);
    12261226
    1227         // simply create a new molecule, wherein the config file is loaded and the manipulation takes place
    1228         molecule *mol = new molecule(periode);
    1229         molecules->insert(mol);
    1230 
    1231         if (argc > 1) { // config file specified as option
    1232                 // 1. : Parse options that just set variables or print help
    1233                 argptr = 1;
    1234                 do {
    1235                         if (argv[argptr][0] == '-') {
    1236                                 cout << Verbose(0) << "Recognized command line argument: " << argv[argptr][1] << ".\n";
    1237                                 argptr++;
    1238                                 switch(argv[argptr-1][1]) {
    1239                                         case 'h':
    1240                                         case 'H':
    1241                                         case '?':
    1242                                                 cout << "MoleCuilder suite" << endl << "==================" << endl << endl;
    1243                                                 cout << "Usage: " << argv[0] << "[config file] [-{acefpsthH?vfrp}] [further arguments]" << endl;
    1244                                                 cout << "or simply " << argv[0] << " without arguments for interactive session." << endl;
    1245                                                 cout << "\t-a Z x1 x2 x3\tAdd new atom of element Z at coordinates (x1,x2,x3)." << endl;
    1246                                                 cout << "\t-A <source>\tCreate adjacency list from bonds parsed from 'dbond'-style file." <<endl;
    1247                                                 cout << "\t-b xx xy xz yy yz zz\tCenter atoms in domain with given symmetric matrix of (xx,xy,xz,yy,yz,zz)." << endl;
    1248                                                 cout << "\t-B <basis>\tSetting basis to store to MPQC config files." << endl;
    1249                                                 cout << "\t-c x1 x2 x3\tCenter atoms in domain with a minimum distance to boundary of (x1,x2,x3)." << endl;
    1250                                                 cout << "\t-D <bond distance>\tDepth-First-Search Analysis of the molecule, giving cycles and tree/back edges." << endl;
    1251                                                 cout << "\t-O\tCenter atoms in origin." << endl;
    1252                                                 cout << "\t-d x1 x2 x3\tDuplicate cell along each axis by given factor." << endl;
    1253                                                 cout << "\t-e <file>\tSets the databases path to be parsed (default: ./)." << endl;
    1254                                                 cout << "\t-E <id> <Z>\tChange atom <id>'s element to <Z>, <id> begins at 0." << endl;
    1255                                                 cout << "\t-f/F <dist> <order>\tFragments the molecule in BOSSANOVA manner (with/out rings compressed) and stores config files in same dir as config (return code 0 - fragmented, 2 - no fragmentation necessary)." << endl;
    1256                                                 cout << "\t-h/-H/-?\tGive this help screen." << endl;
    1257                                                 cout << "\t-m <0/1>\tCalculate (0)/ Align in(1) PAS with greatest EV along z axis." << endl;
    1258                                                 cout << "\t-n\tFast parsing (i.e. no trajectories are looked for)." << endl;
    1259                                                 cout << "\t-N\tGet non-convex-envelope." << endl;
    1260                                                 cout << "\t-o <out>\tGet volume of the convex envelope (and store to tecplot file)." << endl;
    1261                                                 cout << "\t-p <file>\tParse given xyz file and create raw config file from it." << endl;
    1262                                                 cout << "\t-P <file>\tParse given forces file and append as an MD step to config file via Verlet." << endl;
    1263                                                 cout << "\t-r\t\tConvert file from an old pcp syntax." << endl;
    1264                                                 cout << "\t-t x1 x2 x3\tTranslate all atoms by this vector (x1,x2,x3)." << endl;
     1227  // simply create a new molecule, wherein the config file is loaded and the manipulation takes place
     1228  molecule *mol = new molecule(periode);
     1229  molecules->insert(mol);
     1230
     1231  if (argc > 1) { // config file specified as option
     1232    // 1. : Parse options that just set variables or print help
     1233    argptr = 1;
     1234    do {
     1235      if (argv[argptr][0] == '-') {
     1236        cout << Verbose(0) << "Recognized command line argument: " << argv[argptr][1] << ".\n";
     1237        argptr++;
     1238        switch(argv[argptr-1][1]) {
     1239          case 'h':
     1240          case 'H':
     1241          case '?':
     1242            cout << "MoleCuilder suite" << endl << "==================" << endl << endl;
     1243            cout << "Usage: " << argv[0] << "[config file] [-{acefpsthH?vfrp}] [further arguments]" << endl;
     1244            cout << "or simply " << argv[0] << " without arguments for interactive session." << endl;
     1245            cout << "\t-a Z x1 x2 x3\tAdd new atom of element Z at coordinates (x1,x2,x3)." << endl;
     1246            cout << "\t-A <source>\tCreate adjacency list from bonds parsed from 'dbond'-style file." <<endl;
     1247            cout << "\t-b xx xy xz yy yz zz\tCenter atoms in domain with given symmetric matrix of (xx,xy,xz,yy,yz,zz)." << endl;
     1248            cout << "\t-B <basis>\tSetting basis to store to MPQC config files." << endl;
     1249            cout << "\t-c x1 x2 x3\tCenter atoms in domain with a minimum distance to boundary of (x1,x2,x3)." << endl;
     1250            cout << "\t-D <bond distance>\tDepth-First-Search Analysis of the molecule, giving cycles and tree/back edges." << endl;
     1251            cout << "\t-O\tCenter atoms in origin." << endl;
     1252            cout << "\t-d x1 x2 x3\tDuplicate cell along each axis by given factor." << endl;
     1253            cout << "\t-e <file>\tSets the databases path to be parsed (default: ./)." << endl;
     1254            cout << "\t-E <id> <Z>\tChange atom <id>'s element to <Z>, <id> begins at 0." << endl;
     1255            cout << "\t-f/F <dist> <order>\tFragments the molecule in BOSSANOVA manner (with/out rings compressed) and stores config files in same dir as config (return code 0 - fragmented, 2 - no fragmentation necessary)." << endl;
     1256            cout << "\t-h/-H/-?\tGive this help screen." << endl;
     1257            cout << "\t-m <0/1>\tCalculate (0)/ Align in(1) PAS with greatest EV along z axis." << endl;
     1258            cout << "\t-n\tFast parsing (i.e. no trajectories are looked for)." << endl;
     1259            cout << "\t-N\tGet non-convex-envelope." << endl;
     1260            cout << "\t-o <out>\tGet volume of the convex envelope (and store to tecplot file)." << endl;
     1261            cout << "\t-p <file>\tParse given xyz file and create raw config file from it." << endl;
     1262            cout << "\t-P <file>\tParse given forces file and append as an MD step to config file via Verlet." << endl;
     1263            cout << "\t-r\t\tConvert file from an old pcp syntax." << endl;
     1264            cout << "\t-t x1 x2 x3\tTranslate all atoms by this vector (x1,x2,x3)." << endl;
    12651265            cout << "\t-T x1 x2 x3\tTranslate periodically all atoms by this vector (x1,x2,x3)." << endl;
    1266                                                 cout << "\t-S <file> Store temperatures from the config file in <file>." << endl;
    1267                                                 cout << "\t-s x1 x2 x3\tScale all atom coordinates by this vector (x1,x2,x3)." << endl;
    1268                                                 cout << "\t-u rho\tsuspend in water solution and output necessary cell lengths, average density rho and repetition." << endl;
    1269                                                 cout << "\t-v/-V\t\tGives version information." << endl;
    1270                                                 cout << "Note: config files must not begin with '-' !" << endl;
    1271                                                 delete(mol);
    1272                                                 delete(periode);
    1273                                                 return (1);
    1274                                                 break;
    1275                                         case 'v':
    1276                                         case 'V':
    1277                                                 cout << argv[0] << " " << VERSIONSTRING << endl;
    1278                                                 cout << "Build your own molecule position set." << endl;
    1279                                                 delete(mol);
    1280                                                 delete(periode);
    1281                                                 return (1);
    1282                                                 break;
    1283                                         case 'e':
    1284                                                 if ((argptr >= argc) || (argv[argptr][0] == '-')) {
    1285                                                         cerr << "Not enough or invalid arguments for specifying element db: -e <db file>" << endl;
    1286                                                 } else {
    1287                                                         cout << "Using " << argv[argptr] << " as elements database." << endl;
    1288                                                         strncpy (configuration.databasepath, argv[argptr], MAXSTRINGSIZE-1);
    1289                                                         argptr+=1;
    1290                                                 }
    1291                                                 break;
    1292                                         case 'n':
    1293                                                 cout << "I won't parse trajectories." << endl;
    1294                                                 configuration.FastParsing = true;
    1295                                                 break;
    1296                                         default:        // no match? Step on
    1297                                                 argptr++;
    1298                                                 break;
    1299                                 }
    1300                         } else
    1301                                 argptr++;
    1302                 } while (argptr < argc);
    1303 
    1304                 // 2. Parse the element database
    1305                 if (periode->LoadPeriodentafel(configuration.databasepath)) {
    1306                         cout << Verbose(0) << "Element list loaded successfully." << endl;
    1307                         //periode->Output((ofstream *)&cout);
    1308                 } else {
    1309                         cout << Verbose(0) << "Element list loading failed." << endl;
    1310                         return 1;
    1311                 }
    1312                 // 3. Find config file name and parse if possible
    1313                 if (argv[1][0] != '-') {
    1314                         cout << Verbose(0) << "Config file given." << endl;
    1315                         test.open(argv[1], ios::in);
    1316                         if (test == NULL) {
    1317                                 //return (1);
    1318                                 output.open(argv[1], ios::out);
    1319                                 if (output == NULL) {
    1320                                         cout << Verbose(1) << "Specified config file " << argv[1] << " not found." << endl;
    1321                                         config_present = absent;
    1322                                 } else {
    1323                                         cout << "Empty configuration file." << endl;
    1324                                         ConfigFileName = argv[1];
    1325                                         config_present = empty;
    1326                                         output.close();
    1327                                 }
    1328                         } else {
    1329                                 test.close();
    1330                                 ConfigFileName = argv[1];
    1331                                 cout << Verbose(1) << "Specified config file found, parsing ... ";
    1332                                 switch (configuration.TestSyntax(ConfigFileName, periode, mol)) {
    1333                                         case 1:
    1334                                                 cout << "new syntax." << endl;
    1335                                                 configuration.Load(ConfigFileName, periode, mol);
    1336                                                 config_present = present;
    1337                                                 break;
    1338                                         case 0:
    1339                                                 cout << "old syntax." << endl;
    1340                                                 configuration.LoadOld(ConfigFileName, periode, mol);
    1341                                                 config_present = present;
    1342                                                 break;
    1343                                         default:
    1344                                                 cout << "Unknown syntax or empty, yet present file." << endl;
    1345                                                 config_present = empty;
    1346                         }
    1347                         }
    1348                 } else
    1349                         config_present = absent;
    1350                 // 4. parse again through options, now for those depending on elements db and config presence
    1351                 argptr = 1;
    1352                 do {
    1353                         cout << "Current Command line argument: " << argv[argptr] << "." << endl;
    1354                         if (argv[argptr][0] == '-') {
    1355                                 argptr++;
    1356                                 if ((config_present == present) || (config_present == empty)) {
    1357                                         switch(argv[argptr-1][1]) {
    1358                                                 case 'p':
    1359                                                         ExitFlag = 1;
    1360                                                         if ((argptr >= argc) || (argv[argptr][0] == '-')) {
    1361                                                                 ExitFlag = 255;
    1362                                                                 cerr << "Not enough arguments for parsing: -p <xyz file>" << endl;
    1363                                                         } else {
    1364                                                                 SaveFlag = true;
    1365                                                                 cout << Verbose(1) << "Parsing xyz file for new atoms." << endl;
    1366                                                                 if (!mol->AddXYZFile(argv[argptr]))
    1367                                                                         cout << Verbose(2) << "File not found." << endl;
    1368                                                                 else {
    1369                                                                         cout << Verbose(2) << "File found and parsed." << endl;
    1370                                                                         config_present = present;
    1371                                                                 }
    1372                                                         }
    1373                                                         break;
    1374                                                 case 'a':
    1375                                                         ExitFlag = 1;
    1376                                                         if ((argptr >= argc) || (argv[argptr][0] == '-') || (!IsValidNumber(argv[argptr+1]))) {
    1377                                                                 ExitFlag = 255;
    1378                                                                 cerr << "Not enough or invalid arguments for adding atom: -a <element> <x> <y> <z>" << endl;
    1379                                                         } else {
    1380                                                                 SaveFlag = true;
    1381                                                                 cout << Verbose(1) << "Adding new atom with element " << argv[argptr] << " at (" << argv[argptr+1] << "," << argv[argptr+2] << "," << argv[argptr+3] << "), ";
    1382                                                                 first = new atom;
    1383                                                                 first->type = periode->FindElement(atoi(argv[argptr]));
    1384                                                                 if (first->type != NULL)
    1385                                                                         cout << Verbose(2) << "found element " << first->type->name << endl;
    1386                                                                 for (int i=NDIM;i--;)
    1387                                                                         first->x.x[i] = atof(argv[argptr+1+i]);
    1388                                                                 if (first->type != NULL) {
    1389                                                                         mol->AddAtom(first);    // add to molecule
    1390                                                                         if ((config_present == empty) && (mol->AtomCount != 0))
    1391                                                                                 config_present = present;
    1392                                                                 } else
    1393                                                                         cerr << Verbose(1) << "Could not find the specified element." << endl;
    1394                                                                 argptr+=4;
    1395                                                         }
    1396                                                         break;
    1397                                                 default:        // no match? Don't step on (this is done in next switch's default)
    1398                                                         break;
    1399                                         }
    1400                                 }
    1401                                 if (config_present == present) {
    1402                                         switch(argv[argptr-1][1]) {
    1403                                                 case 'B':
    1404                                                         if ((argptr >= argc) || (argv[argptr][0] == '-')) {
    1405                                                                 ExitFlag = 255;
    1406                                                                 cerr << "Not enough or invalid arguments given for setting MPQC basis: -B <basis name>" << endl;
    1407                                                         } else {
    1408                                                                 configuration.basis = argv[argptr];
    1409                                                                 cout << Verbose(1) << "Setting MPQC basis to " << configuration.basis << "." << endl;
    1410                                                                 argptr+=1;
    1411                                                         }
    1412                                                         break;
    1413                                                 case 'D':
    1414                                                         ExitFlag = 1;
    1415                                                         {
    1416                                                                 cout << Verbose(1) << "Depth-First-Search Analysis." << endl;
    1417                                                                 MoleculeLeafClass *Subgraphs = NULL;                    // list of subgraphs from DFS analysis
    1418                                                                 int *MinimumRingSize = new int[mol->AtomCount];
    1419                                                                 atom ***ListOfLocalAtoms = NULL;
    1420                                                                 int FragmentCounter = 0;
    1421                                                                 class StackClass<bond *> *BackEdgeStack = NULL;
    1422                                                                 class StackClass<bond *> *LocalBackEdgeStack = NULL;
    1423                                                                 mol->CreateAdjacencyList((ofstream *)&cout, atof(argv[argptr]), configuration.GetIsAngstroem());
    1424                                                                 mol->CreateListOfBondsPerAtom((ofstream *)&cout);
    1425                                                                 Subgraphs = mol->DepthFirstSearchAnalysis((ofstream *)&cout, BackEdgeStack);
    1426                                                                 if (Subgraphs != NULL) {
    1427                                                                         Subgraphs->next->FillBondStructureFromReference((ofstream *)&cout, mol, (FragmentCounter = 0), ListOfLocalAtoms, false);        // we want to keep the created ListOfLocalAtoms
    1428                                                                         while (Subgraphs->next != NULL) {
    1429                                                                                 Subgraphs = Subgraphs->next;
    1430                                                                                 LocalBackEdgeStack = new StackClass<bond *> (Subgraphs->Leaf->BondCount);
    1431                                                                                 Subgraphs->Leaf->PickLocalBackEdges((ofstream *)&cout, ListOfLocalAtoms[FragmentCounter++], BackEdgeStack, LocalBackEdgeStack);
    1432                                                                                 Subgraphs->Leaf->CyclicStructureAnalysis((ofstream *)&cout, BackEdgeStack, MinimumRingSize);
    1433                                                                                 delete(LocalBackEdgeStack);
    1434                                                                                 delete(Subgraphs->previous);
    1435                                                                         }
    1436                                                                         delete(Subgraphs);
    1437                                                                         for (int i=0;i<FragmentCounter;i++)
    1438                                                                                 Free((void **)&ListOfLocalAtoms[FragmentCounter], "ParseCommandLineOptions: **ListOfLocalAtoms[]");
    1439                                                                         Free((void **)&ListOfLocalAtoms, "ParseCommandLineOptions: ***ListOfLocalAtoms");
    1440                                                                 }
    1441                                                                 delete(BackEdgeStack);
    1442                                                                 delete[](MinimumRingSize);
    1443                                                         }
    1444                                                         //argptr+=1;
    1445                                                         break;
    1446                                                 case 'E':
    1447                                                         ExitFlag = 1;
    1448                                                         if ((argptr+1 >= argc) || (!IsValidNumber(argv[argptr])) || (argv[argptr+1][0] == '-')) {
    1449                                                                 ExitFlag = 255;
    1450                                                                 cerr << "Not enough or invalid arguments given for changing element: -E <atom nr.> <element>" << endl;
    1451                                                         } else {
    1452                                                                 SaveFlag = true;
    1453                                                                 cout << Verbose(1) << "Changing atom " << argv[argptr] << " to element " << argv[argptr+1] << "." << endl;
    1454                                                                 first = mol->FindAtom(atoi(argv[argptr]));
    1455                                                                 first->type = periode->FindElement(atoi(argv[argptr+1]));
    1456                                                                 argptr+=2;
    1457                                                         }
    1458                                                         break;
    1459                                                 case 'A':
    1460                                                         ExitFlag = 1;
    1461                                                         if ((argptr >= argc) || (argv[argptr][0] == '-')) {
    1462                                                                 ExitFlag =255;
    1463                                                                 cerr << "Missing source file for bonds in molecule: -A <bond sourcefile>" << endl;
    1464                                                         } else {
    1465                                                                 cout << "Parsing bonds from " << argv[argptr] << "." << endl;
    1466                                                                 ifstream *input = new ifstream(argv[argptr]);
    1467                                                                 mol->CreateAdjacencyList2((ofstream *)&cout, input);
    1468                                                                 input->close();
    1469                                                                 argptr+=1;
    1470                                                         }
    1471                                                         break;
    1472                                                 case 'N':
    1473                                                         ExitFlag = 1;
    1474                                                         if ((argptr+1 >= argc) || (argv[argptr+1][0] == '-')){
    1475                                                                 ExitFlag = 255;
    1476                                                                 cerr << "Not enough or invalid arguments given for non-convex envelope: -o <radius> <tecplot output file>" << endl;
    1477                                                         } else {
    1478                                                                 class Tesselation T;
    1479                                                                 string filename(argv[argptr+1]);
    1480                                                                 filename.append(".csv");
    1481                                                                 cout << Verbose(0) << "Evaluating non-convex envelope.";
    1482                                                                 cout << Verbose(1) << "Using rolling ball of radius " << atof(argv[argptr]) << " and storing tecplot data in " << argv[argptr+1] << "." << endl;
    1483                                                                 LinkedCell LCList(mol, atof(argv[argptr])*2.);
    1484                                                                 Find_non_convex_border((ofstream *)&cout, mol, &T, &LCList, argv[argptr+1], atof(argv[argptr]));
    1485                                                                 //FindDistributionOfEllipsoids((ofstream *)&cout, &T, &LCList, N, number, filename.c_str());
    1486                                                                 argptr+=2;
    1487                                                         }
    1488                                                         break;
    1489                                                 case 'S':
    1490                                                         ExitFlag = 1;
    1491                                                         if ((argptr >= argc) || (argv[argptr][0] == '-')) {
    1492                                                                 ExitFlag = 255;
    1493                                                                 cerr << "Not enough or invalid arguments given for storing tempature: -S <temperature file>" << endl;
    1494                                                         } else {
    1495                                                                 cout << Verbose(1) << "Storing temperatures in " << argv[argptr] << "." << endl;
    1496                                                                 ofstream *output = new ofstream(argv[argptr], ios::trunc);
    1497                                                                 if (!mol->OutputTemperatureFromTrajectories((ofstream *)&cout, 0, mol->MDSteps, output))
    1498                                                                         cout << Verbose(2) << "File could not be written." << endl;
    1499                                                                 else
    1500                                                                         cout << Verbose(2) << "File stored." << endl;
    1501                                                                 output->close();
    1502                                                                 delete(output);
    1503                                                                 argptr+=1;
    1504                                                         }
    1505                                                         break;
    1506                                                 case 'P':
    1507                                                         ExitFlag = 1;
    1508                                                         if ((argptr >= argc) || (argv[argptr][0] == '-')) {
    1509                                                                 ExitFlag = 255;
    1510                                                                 cerr << "Not enough or invalid arguments given for parsing and integrating forces: -P <forces file>" << endl;
    1511                                                         } else {
    1512                                                                 SaveFlag = true;
    1513                                                                 cout << Verbose(1) << "Parsing forces file and Verlet integrating." << endl;
    1514                                                                 if (!mol->VerletForceIntegration(argv[argptr], configuration.Deltat, configuration.GetIsAngstroem()))
    1515                                                                         cout << Verbose(2) << "File not found." << endl;
    1516                                                                 else
    1517                                                                         cout << Verbose(2) << "File found and parsed." << endl;
    1518                                                                 argptr+=1;
    1519                                                         }
    1520                                                         break;
    1521                                                 case 't':
    1522                                                         ExitFlag = 1;
    1523                                                         if ((argptr+2 >= argc) || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1])) || (!IsValidNumber(argv[argptr+2])) ) {
    1524                                                                 ExitFlag = 255;
    1525                                                                 cerr << "Not enough or invalid arguments given for translation: -t <x> <y> <z>" << endl;
    1526                                                         } else {
    1527                                                                 ExitFlag = 1;
    1528                                                                 SaveFlag = true;
    1529                                                                 cout << Verbose(1) << "Translating all ions to new origin." << endl;
    1530                                                                 for (int i=NDIM;i--;)
    1531                                                                         x.x[i] = atof(argv[argptr+i]);
    1532                                                                 mol->Translate((const Vector *)&x);
    1533                                                                 argptr+=3;
    1534                                                         }
    1535                                                         break;
     1266            cout << "\t-S <file> Store temperatures from the config file in <file>." << endl;
     1267            cout << "\t-s x1 x2 x3\tScale all atom coordinates by this vector (x1,x2,x3)." << endl;
     1268            cout << "\t-u rho\tsuspend in water solution and output necessary cell lengths, average density rho and repetition." << endl;
     1269            cout << "\t-v/-V\t\tGives version information." << endl;
     1270            cout << "Note: config files must not begin with '-' !" << endl;
     1271            delete(mol);
     1272            delete(periode);
     1273            return (1);
     1274            break;
     1275          case 'v':
     1276          case 'V':
     1277            cout << argv[0] << " " << VERSIONSTRING << endl;
     1278            cout << "Build your own molecule position set." << endl;
     1279            delete(mol);
     1280            delete(periode);
     1281            return (1);
     1282            break;
     1283          case 'e':
     1284            if ((argptr >= argc) || (argv[argptr][0] == '-')) {
     1285              cerr << "Not enough or invalid arguments for specifying element db: -e <db file>" << endl;
     1286            } else {
     1287              cout << "Using " << argv[argptr] << " as elements database." << endl;
     1288              strncpy (configuration.databasepath, argv[argptr], MAXSTRINGSIZE-1);
     1289              argptr+=1;
     1290            }
     1291            break;
     1292          case 'n':
     1293            cout << "I won't parse trajectories." << endl;
     1294            configuration.FastParsing = true;
     1295            break;
     1296          default:  // no match? Step on
     1297            argptr++;
     1298            break;
     1299        }
     1300      } else
     1301        argptr++;
     1302    } while (argptr < argc);
     1303
     1304    // 2. Parse the element database
     1305    if (periode->LoadPeriodentafel(configuration.databasepath)) {
     1306      cout << Verbose(0) << "Element list loaded successfully." << endl;
     1307      //periode->Output((ofstream *)&cout);
     1308    } else {
     1309      cout << Verbose(0) << "Element list loading failed." << endl;
     1310      return 1;
     1311    }
     1312    // 3. Find config file name and parse if possible
     1313    if (argv[1][0] != '-') {
     1314      cout << Verbose(0) << "Config file given." << endl;
     1315      test.open(argv[1], ios::in);
     1316      if (test == NULL) {
     1317        //return (1);
     1318        output.open(argv[1], ios::out);
     1319        if (output == NULL) {
     1320          cout << Verbose(1) << "Specified config file " << argv[1] << " not found." << endl;
     1321          config_present = absent;
     1322        } else {
     1323          cout << "Empty configuration file." << endl;
     1324          ConfigFileName = argv[1];
     1325          config_present = empty;
     1326          output.close();
     1327        }
     1328      } else {
     1329        test.close();
     1330        ConfigFileName = argv[1];
     1331        cout << Verbose(1) << "Specified config file found, parsing ... ";
     1332        switch (configuration.TestSyntax(ConfigFileName, periode, mol)) {
     1333          case 1:
     1334            cout << "new syntax." << endl;
     1335            configuration.Load(ConfigFileName, periode, mol);
     1336            config_present = present;
     1337            break;
     1338          case 0:
     1339            cout << "old syntax." << endl;
     1340            configuration.LoadOld(ConfigFileName, periode, mol);
     1341            config_present = present;
     1342            break;
     1343          default:
     1344            cout << "Unknown syntax or empty, yet present file." << endl;
     1345            config_present = empty;
     1346      }
     1347      }
     1348    } else
     1349      config_present = absent;
     1350    // 4. parse again through options, now for those depending on elements db and config presence
     1351    argptr = 1;
     1352    do {
     1353      cout << "Current Command line argument: " << argv[argptr] << "." << endl;
     1354      if (argv[argptr][0] == '-') {
     1355        argptr++;
     1356        if ((config_present == present) || (config_present == empty)) {
     1357          switch(argv[argptr-1][1]) {
     1358            case 'p':
     1359              ExitFlag = 1;
     1360              if ((argptr >= argc) || (argv[argptr][0] == '-')) {
     1361                ExitFlag = 255;
     1362                cerr << "Not enough arguments for parsing: -p <xyz file>" << endl;
     1363              } else {
     1364                SaveFlag = true;
     1365                cout << Verbose(1) << "Parsing xyz file for new atoms." << endl;
     1366                if (!mol->AddXYZFile(argv[argptr]))
     1367                  cout << Verbose(2) << "File not found." << endl;
     1368                else {
     1369                  cout << Verbose(2) << "File found and parsed." << endl;
     1370                  config_present = present;
     1371                }
     1372              }
     1373              break;
     1374            case 'a':
     1375              ExitFlag = 1;
     1376              if ((argptr >= argc) || (argv[argptr][0] == '-') || (!IsValidNumber(argv[argptr+1]))) {
     1377                ExitFlag = 255;
     1378                cerr << "Not enough or invalid arguments for adding atom: -a <element> <x> <y> <z>" << endl;
     1379              } else {
     1380                SaveFlag = true;
     1381                cout << Verbose(1) << "Adding new atom with element " << argv[argptr] << " at (" << argv[argptr+1] << "," << argv[argptr+2] << "," << argv[argptr+3] << "), ";
     1382                first = new atom;
     1383                first->type = periode->FindElement(atoi(argv[argptr]));
     1384                if (first->type != NULL)
     1385                  cout << Verbose(2) << "found element " << first->type->name << endl;
     1386                for (int i=NDIM;i--;)
     1387                  first->x.x[i] = atof(argv[argptr+1+i]);
     1388                if (first->type != NULL) {
     1389                  mol->AddAtom(first);  // add to molecule
     1390                  if ((config_present == empty) && (mol->AtomCount != 0))
     1391                    config_present = present;
     1392                } else
     1393                  cerr << Verbose(1) << "Could not find the specified element." << endl;
     1394                argptr+=4;
     1395              }
     1396              break;
     1397            default:  // no match? Don't step on (this is done in next switch's default)
     1398              break;
     1399          }
     1400        }
     1401        if (config_present == present) {
     1402          switch(argv[argptr-1][1]) {
     1403            case 'B':
     1404              if ((argptr >= argc) || (argv[argptr][0] == '-')) {
     1405                ExitFlag = 255;
     1406                cerr << "Not enough or invalid arguments given for setting MPQC basis: -B <basis name>" << endl;
     1407              } else {
     1408                configuration.basis = argv[argptr];
     1409                cout << Verbose(1) << "Setting MPQC basis to " << configuration.basis << "." << endl;
     1410                argptr+=1;
     1411              }
     1412              break;
     1413            case 'D':
     1414              ExitFlag = 1;
     1415              {
     1416                cout << Verbose(1) << "Depth-First-Search Analysis." << endl;
     1417                MoleculeLeafClass *Subgraphs = NULL;      // list of subgraphs from DFS analysis
     1418                int *MinimumRingSize = new int[mol->AtomCount];
     1419                atom ***ListOfLocalAtoms = NULL;
     1420                int FragmentCounter = 0;
     1421                class StackClass<bond *> *BackEdgeStack = NULL;
     1422                class StackClass<bond *> *LocalBackEdgeStack = NULL;
     1423                mol->CreateAdjacencyList((ofstream *)&cout, atof(argv[argptr]), configuration.GetIsAngstroem());
     1424                mol->CreateListOfBondsPerAtom((ofstream *)&cout);
     1425                Subgraphs = mol->DepthFirstSearchAnalysis((ofstream *)&cout, BackEdgeStack);
     1426                if (Subgraphs != NULL) {
     1427                  Subgraphs->next->FillBondStructureFromReference((ofstream *)&cout, mol, (FragmentCounter = 0), ListOfLocalAtoms, false);  // we want to keep the created ListOfLocalAtoms
     1428                  while (Subgraphs->next != NULL) {
     1429                    Subgraphs = Subgraphs->next;
     1430                    LocalBackEdgeStack = new StackClass<bond *> (Subgraphs->Leaf->BondCount);
     1431                    Subgraphs->Leaf->PickLocalBackEdges((ofstream *)&cout, ListOfLocalAtoms[FragmentCounter++], BackEdgeStack, LocalBackEdgeStack);
     1432                    Subgraphs->Leaf->CyclicStructureAnalysis((ofstream *)&cout, BackEdgeStack, MinimumRingSize);
     1433                    delete(LocalBackEdgeStack);
     1434                    delete(Subgraphs->previous);
     1435                  }
     1436                  delete(Subgraphs);
     1437                  for (int i=0;i<FragmentCounter;i++)
     1438                    Free((void **)&ListOfLocalAtoms[FragmentCounter], "ParseCommandLineOptions: **ListOfLocalAtoms[]");
     1439                  Free((void **)&ListOfLocalAtoms, "ParseCommandLineOptions: ***ListOfLocalAtoms");
     1440                }
     1441                delete(BackEdgeStack);
     1442                delete[](MinimumRingSize);
     1443              }
     1444              //argptr+=1;
     1445              break;
     1446            case 'E':
     1447              ExitFlag = 1;
     1448              if ((argptr+1 >= argc) || (!IsValidNumber(argv[argptr])) || (argv[argptr+1][0] == '-')) {
     1449                ExitFlag = 255;
     1450                cerr << "Not enough or invalid arguments given for changing element: -E <atom nr.> <element>" << endl;
     1451              } else {
     1452                SaveFlag = true;
     1453                cout << Verbose(1) << "Changing atom " << argv[argptr] << " to element " << argv[argptr+1] << "." << endl;
     1454                first = mol->FindAtom(atoi(argv[argptr]));
     1455                first->type = periode->FindElement(atoi(argv[argptr+1]));
     1456                argptr+=2;
     1457              }
     1458              break;
     1459            case 'A':
     1460              ExitFlag = 1;
     1461              if ((argptr >= argc) || (argv[argptr][0] == '-')) {
     1462                ExitFlag =255;
     1463                cerr << "Missing source file for bonds in molecule: -A <bond sourcefile>" << endl;
     1464              } else {
     1465                cout << "Parsing bonds from " << argv[argptr] << "." << endl;
     1466                ifstream *input = new ifstream(argv[argptr]);
     1467                mol->CreateAdjacencyList2((ofstream *)&cout, input);
     1468                input->close();
     1469                argptr+=1;
     1470              }
     1471              break;
     1472            case 'N':
     1473              ExitFlag = 1;
     1474              if ((argptr+1 >= argc) || (argv[argptr+1][0] == '-')){
     1475                ExitFlag = 255;
     1476                cerr << "Not enough or invalid arguments given for non-convex envelope: -o <radius> <tecplot output file>" << endl;
     1477              } else {
     1478                class Tesselation T;
     1479                string filename(argv[argptr+1]);
     1480                filename.append(".csv");
     1481                cout << Verbose(0) << "Evaluating non-convex envelope.";
     1482                cout << Verbose(1) << "Using rolling ball of radius " << atof(argv[argptr]) << " and storing tecplot data in " << argv[argptr+1] << "." << endl;
     1483                LinkedCell LCList(mol, atof(argv[argptr])*2.);
     1484                Find_non_convex_border((ofstream *)&cout, mol, &T, &LCList, argv[argptr+1], atof(argv[argptr]));
     1485                //FindDistributionOfEllipsoids((ofstream *)&cout, &T, &LCList, N, number, filename.c_str());
     1486                argptr+=2;
     1487              }
     1488              break;
     1489            case 'S':
     1490              ExitFlag = 1;
     1491              if ((argptr >= argc) || (argv[argptr][0] == '-')) {
     1492                ExitFlag = 255;
     1493                cerr << "Not enough or invalid arguments given for storing tempature: -S <temperature file>" << endl;
     1494              } else {
     1495                cout << Verbose(1) << "Storing temperatures in " << argv[argptr] << "." << endl;
     1496                ofstream *output = new ofstream(argv[argptr], ios::trunc);
     1497                if (!mol->OutputTemperatureFromTrajectories((ofstream *)&cout, 0, mol->MDSteps, output))
     1498                  cout << Verbose(2) << "File could not be written." << endl;
     1499                else
     1500                  cout << Verbose(2) << "File stored." << endl;
     1501                output->close();
     1502                delete(output);
     1503                argptr+=1;
     1504              }
     1505              break;
     1506            case 'P':
     1507              ExitFlag = 1;
     1508              if ((argptr >= argc) || (argv[argptr][0] == '-')) {
     1509                ExitFlag = 255;
     1510                cerr << "Not enough or invalid arguments given for parsing and integrating forces: -P <forces file>" << endl;
     1511              } else {
     1512                SaveFlag = true;
     1513                cout << Verbose(1) << "Parsing forces file and Verlet integrating." << endl;
     1514                if (!mol->VerletForceIntegration(argv[argptr], configuration.Deltat, configuration.GetIsAngstroem()))
     1515                  cout << Verbose(2) << "File not found." << endl;
     1516                else
     1517                  cout << Verbose(2) << "File found and parsed." << endl;
     1518                argptr+=1;
     1519              }
     1520              break;
     1521            case 't':
     1522              ExitFlag = 1;
     1523              if ((argptr+2 >= argc) || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1])) || (!IsValidNumber(argv[argptr+2])) ) {
     1524                ExitFlag = 255;
     1525                cerr << "Not enough or invalid arguments given for translation: -t <x> <y> <z>" << endl;
     1526              } else {
     1527                ExitFlag = 1;
     1528                SaveFlag = true;
     1529                cout << Verbose(1) << "Translating all ions to new origin." << endl;
     1530                for (int i=NDIM;i--;)
     1531                  x.x[i] = atof(argv[argptr+i]);
     1532                mol->Translate((const Vector *)&x);
     1533                argptr+=3;
     1534              }
     1535              break;
    15361536            case 'T':
    15371537              ExitFlag = 1;
    … …  
    15491549              }
    15501550              break;
    1551                                                 case 's':
    1552                                                         ExitFlag = 1;
    1553                                                         if ((argptr >= argc) || (!IsValidNumber(argv[argptr])) ) {
    1554                                                                 ExitFlag = 255;
    1555                                                                 cerr << "Not enough or invalid arguments given for scaling: -s <factor/[factor_x]> [factor_y] [factor_z]" << endl;
    1556                                                         } else {
    1557                                                                 SaveFlag = true;
    1558                                                                 j = -1;
    1559                                                                 cout << Verbose(1) << "Scaling all ion positions by factor." << endl;
    1560                                                                 factor = new double[NDIM];
    1561                                                                 factor[0] = atof(argv[argptr]);
    1562                                                                 if ((argptr < argc) && (IsValidNumber(argv[argptr])))
    1563                                                                         argptr++;
    1564                                                                 factor[1] = atof(argv[argptr]);
    1565                                                                 if ((argptr < argc) && (IsValidNumber(argv[argptr])))
    1566                                                                         argptr++;
    1567                                                                 factor[2] = atof(argv[argptr]);
    1568                                                                 mol->Scale(&factor);
    1569                                                                 for (int i=0;i<NDIM;i++) {
    1570                                                                         j += i+1;
    1571                                                                         x.x[i] = atof(argv[NDIM+i]);
    1572                                                                         mol->cell_size[j]*=factor[i];
    1573                                                                 }
    1574                                                                 delete[](factor);
    1575                                                                 argptr+=1;
    1576                                                         }
    1577                                                         break;
    1578                                                 case 'b':
    1579                                                         ExitFlag = 1;
    1580                                                         if ((argptr+5 >= argc) || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1])) || (!IsValidNumber(argv[argptr+2])) || (!IsValidNumber(argv[argptr+3])) || (!IsValidNumber(argv[argptr+4])) || (!IsValidNumber(argv[argptr+5])) ) {
    1581                                                                 ExitFlag = 255;
    1582                                                                 cerr << "Not enough or invalid arguments given for centering in box: -b <xx> <xy> <xz> <yy> <yz> <zz>" << endl;
    1583                                                         } else {
    1584                                                                 SaveFlag = true;
    1585                                                                 cout << Verbose(1) << "Centering atoms in config file within given simulation box." << endl;
    1586                                                                 for (int i=0;i<6;i++) {
    1587                                                                         mol->cell_size[i] = atof(argv[argptr+i]);
    1588                                                                 }
    1589                                                                 // center
    1590                                                                 mol->CenterInBox((ofstream *)&cout);
     1551            case 's':
     1552              ExitFlag = 1;
     1553              if ((argptr >= argc) || (!IsValidNumber(argv[argptr])) ) {
     1554                ExitFlag = 255;
     1555                cerr << "Not enough or invalid arguments given for scaling: -s <factor/[factor_x]> [factor_y] [factor_z]" << endl;
     1556              } else {
     1557                SaveFlag = true;
     1558                j = -1;
     1559                cout << Verbose(1) << "Scaling all ion positions by factor." << endl;
     1560                factor = new double[NDIM];
     1561                factor[0] = atof(argv[argptr]);
     1562                if ((argptr < argc) && (IsValidNumber(argv[argptr])))
     1563                  argptr++;
     1564                factor[1] = atof(argv[argptr]);
     1565                if ((argptr < argc) && (IsValidNumber(argv[argptr])))
     1566                  argptr++;
     1567                factor[2] = atof(argv[argptr]);
     1568                mol->Scale(&factor);
     1569                for (int i=0;i<NDIM;i++) {
     1570                  j += i+1;
     1571                  x.x[i] = atof(argv[NDIM+i]);
     1572                  mol->cell_size[j]*=factor[i];
     1573                }
     1574                delete[](factor);
     1575                argptr+=1;
     1576              }
     1577              break;
     1578            case 'b':
     1579              ExitFlag = 1;
     1580              if ((argptr+5 >= argc) || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1])) || (!IsValidNumber(argv[argptr+2])) || (!IsValidNumber(argv[argptr+3])) || (!IsValidNumber(argv[argptr+4])) || (!IsValidNumber(argv[argptr+5])) ) {
     1581                ExitFlag = 255;
     1582                cerr << "Not enough or invalid arguments given for centering in box: -b <xx> <xy> <xz> <yy> <yz> <zz>" << endl;
     1583              } else {
     1584                SaveFlag = true;
     1585                cout << Verbose(1) << "Centering atoms in config file within given simulation box." << endl;
     1586                for (int i=0;i<6;i++) {
     1587                  mol->cell_size[i] = atof(argv[argptr+i]);
     1588                }
     1589                // center
     1590                mol->CenterInBox((ofstream *)&cout);
    15911591                argptr+=6;
    1592                                                         }
    1593                                                         break;
    1594                                                 case 'c':
    1595                                                         ExitFlag = 1;
    1596                                                         if ((argptr+2 >= argc) || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1])) || (!IsValidNumber(argv[argptr+2])) ) {
    1597                                                                 ExitFlag = 255;
    1598                                                                 cerr << "Not enough or invalid arguments given for centering with boundary: -c <boundary_x> <boundary_y> <boundary_z>" << endl;
    1599                                                         } else {
    1600                                                                 SaveFlag = true;
    1601                                                                 j = -1;
    1602                                                                 cout << Verbose(1) << "Centering atoms in config file within given additional boundary." << endl;
    1603                                                                 // make every coordinate positive
    1604                                                                 mol->CenterEdge((ofstream *)&cout, &x);
    1605                                                                 // update Box of atoms by boundary
    1606                                                                 mol->SetBoxDimension(&x);
    1607                                                                 // translate each coordinate by boundary
    1608                                                                 j=-1;
    1609                                                                 for (int i=0;i<NDIM;i++) {
    1610                                                                         j += i+1;
    1611                                                                         x.x[i] = atof(argv[argptr++]);
    1612                                                                         mol->cell_size[j] += x.x[i]*2.;
    1613                                                                 }
    1614                                                                 mol->Translate((const Vector *)&x);
     1592              }
     1593              break;
     1594            case 'c':
     1595              ExitFlag = 1;
     1596              if ((argptr+2 >= argc) || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1])) || (!IsValidNumber(argv[argptr+2])) ) {
     1597                ExitFlag = 255;
     1598                cerr << "Not enough or invalid arguments given for centering with boundary: -c <boundary_x> <boundary_y> <boundary_z>" << endl;
     1599              } else {
     1600                SaveFlag = true;
     1601                j = -1;
     1602                cout << Verbose(1) << "Centering atoms in config file within given additional boundary." << endl;
     1603                // make every coordinate positive
     1604                mol->CenterEdge((ofstream *)&cout, &x);
     1605                // update Box of atoms by boundary
     1606                mol->SetBoxDimension(&x);
     1607                // translate each coordinate by boundary
     1608                j=-1;
     1609                for (int i=0;i<NDIM;i++) {
     1610                  j += i+1;
     1611                  x.x[i] = atof(argv[argptr++]);
     1612                  mol->cell_size[j] += x.x[i]*2.;
     1613                }
     1614                mol->Translate((const Vector *)&x);
    16151615                argptr+=3;
    1616                                                         }
    1617                                                         break;
    1618                                                 case 'O':
    1619                                                         ExitFlag = 1;
    1620                                                         SaveFlag = true;
    1621                                                         cout << Verbose(1) << "Centering atoms in origin." << endl;
    1622                                                         mol->CenterOrigin((ofstream *)&cout, &x);
    1623                                                         mol->SetBoxDimension(&x);
     1616              }
     1617              break;
     1618            case 'O':
     1619              ExitFlag = 1;
     1620              SaveFlag = true;
     1621              cout << Verbose(1) << "Centering atoms in origin." << endl;
     1622              mol->CenterOrigin((ofstream *)&cout, &x);
     1623              mol->SetBoxDimension(&x);
    16241624              argptr+=0;
    1625                                                         break;
    1626                                                 case 'r':
    1627                                                         ExitFlag = 1;
    1628                                                         SaveFlag = true;
    1629                                                         cout << Verbose(1) << "Converting config file from supposed old to new syntax." << endl;
    1630                                                         break;
    1631                                                 case 'F':
    1632                                                 case 'f':
    1633                                                         ExitFlag = 1;
    1634                                                         if ((argptr+1 >= argc) || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1]))) {
    1635                                                                 ExitFlag = 255;
    1636                                                                 cerr << "Not enough or invalid arguments for fragmentation: -f <max. bond distance> <bond order>" << endl;
    1637                                                         } else {
    1638                                                                 cout << "Fragmenting molecule with bond distance " << argv[argptr] << " angstroem, order of " << argv[argptr+1] << "." << endl;
    1639                                                                 cout << Verbose(0) << "Creating connection matrix..." << endl;
    1640                                                                 start = clock();
    1641                                                                 mol->CreateAdjacencyList((ofstream *)&cout, atof(argv[argptr++]), configuration.GetIsAngstroem());
    1642                                                                 cout << Verbose(0) << "Fragmenting molecule with current connection matrix ..." << endl;
    1643                                                                 if (mol->first->next != mol->last) {
    1644                                                                         ExitFlag = mol->FragmentMolecule((ofstream *)&cout, atoi(argv[argptr]), &configuration);
    1645                                                                 }
    1646                                                                 end = clock();
    1647                                                                 cout << Verbose(0) << "Clocks for this operation: " << (end-start) << ", time: " << ((double)(end-start)/CLOCKS_PER_SEC) << "s." << endl;
    1648                                                                 argptr+=2;
    1649                                                         }
    1650                                                         break;
    1651                                                 case 'm':
    1652                                                         ExitFlag = 1;
    1653                                                         j = atoi(argv[argptr++]);
    1654                                                         if ((j<0) || (j>1)) {
    1655                                                                 cerr << Verbose(1) << "ERROR: Argument of '-m' should be either 0 for no-rotate or 1 for rotate." << endl;
    1656                                                                 j = 0;
    1657                                                         }
    1658                                                         if (j) {
    1659                                                                 SaveFlag = true;
    1660                                                                 cout << Verbose(0) << "Converting to prinicipal axis system." << endl;
    1661                                                         } else
    1662                                                                 cout << Verbose(0) << "Evaluating prinicipal axis." << endl;
    1663                                                         mol->PrincipalAxisSystem((ofstream *)&cout, (bool)j);
    1664                                                         break;
    1665                                                 case 'o':
    1666                                                         ExitFlag = 1;
    1667                                                         if ((argptr >= argc) || (argv[argptr][0] == '-')){
    1668                                                                 ExitFlag = 255;
    1669                                                                 cerr << "Not enough or invalid arguments given for convex envelope: -o <tecplot output file>" << endl;
    1670                                                         } else {
    1671                                                                 cout << Verbose(0) << "Evaluating volume of the convex envelope.";
    1672                                                                 cout << Verbose(1) << "Storing tecplot data in " << argv[argptr] << "." << endl;
    1673                                                                 VolumeOfConvexEnvelope((ofstream *)&cout, argv[argptr], &configuration, NULL, mol);
    1674                                                                 argptr+=1;
    1675                                                         }
    1676                                                         break;
    1677                                                 case 'U':
    1678                                                         ExitFlag = 1;
    1679                                                         if ((argptr+1 >= argc) || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1])) ) {
    1680                                                                 ExitFlag = 255;
    1681                                                                 cerr << "Not enough or invalid arguments given for suspension with specified volume: -U <volume> <density>" << endl;
    1682                                                                 volume = -1; // for case 'u': don't print error again
    1683                                                         } else {
    1684                                                                 volume = atof(argv[argptr++]);
    1685                                                                 cout << Verbose(0) << "Using " << volume << " angstrom^3 as the volume instead of convex envelope one's." << endl;
    1686                                                         }
    1687                                                 case 'u':
    1688                                                         ExitFlag = 1;
    1689                                                         if ((argptr >= argc) || (!IsValidNumber(argv[argptr])) ) {
    1690                                                                 if (volume != -1)
    1691                                                                         ExitFlag = 255;
    1692                                                                         cerr << "Not enough arguments given for suspension: -u <density>" << endl;
    1693                                                         } else {
    1694                                                                 double density;
    1695                                                                 SaveFlag = true;
    1696                                                                 cout << Verbose(0) << "Evaluating necessary cell volume for a cluster suspended in water.";
    1697                                                                 density = atof(argv[argptr++]);
    1698                                                                 if (density < 1.0) {
    1699                                                                         cerr << Verbose(0) << "Density must be greater than 1.0g/cm^3 !" << endl;
    1700                                                                         density = 1.3;
    1701                                                                 }
    1702 //                                                              for(int i=0;i<NDIM;i++) {
    1703 //                                                                      repetition[i] = atoi(argv[argptr++]);
    1704 //                                                                      if (repetition[i] < 1)
    1705 //                                                                              cerr << Verbose(0) << "ERROR: repetition value must be greater 1!" << endl;
    1706 //                                                                      repetition[i] = 1;
    1707 //                                                              }
    1708                                                                 PrepareClustersinWater((ofstream *)&cout, &configuration, mol, volume, density);        // if volume == 0, will calculate from ConvexEnvelope
    1709                                                         }
    1710                                                         break;
    1711                                                 case 'd':
    1712                                                         ExitFlag = 1;
    1713                                                         if ((argptr+2 >= argc) || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1])) || (!IsValidNumber(argv[argptr+2])) ) {
    1714                                                                 ExitFlag = 255;
    1715                                                                 cerr << "Not enough or invalid arguments given for repeating cells: -d <repeat_x> <repeat_y> <repeat_z>" << endl;
    1716                                                         } else {
    1717                                                                 SaveFlag = true;
    1718                                                                 for (int axis = 1; axis <= NDIM; axis++) {
    1719                                                                         int faktor = atoi(argv[argptr++]);
    1720                                                                         int count;
    1721                                                                         element ** Elements;
    1722                                                                         Vector ** vectors;
    1723                                                                         if (faktor < 1) {
    1724                                                                                 cerr << Verbose(0) << "ERROR: Repetition faktor mus be greater than 1!" << endl;
    1725                                                                                 faktor = 1;
    1726                                                                         }
    1727                                                                         mol->CountAtoms((ofstream *)&cout);     // recount atoms
    1728                                                                         if (mol->AtomCount != 0) {      // if there is more than none
    1729                                                                                 count = mol->AtomCount; // is changed becausing of adding, thus has to be stored away beforehand
    1730                                                                                 Elements = new element *[count];
    1731                                                                                 vectors = new Vector *[count];
    1732                                                                                 j = 0;
    1733                                                                                 first = mol->start;
    1734                                                                                 while (first->next != mol->end) {       // make a list of all atoms with coordinates and element
    1735                                                                                         first = first->next;
    1736                                                                                         Elements[j] = first->type;
    1737                                                                                         vectors[j] = &first->x;
    1738                                                                                         j++;
    1739                                                                                 }
    1740                                                                                 if (count != j)
    1741                                                                                         cout << Verbose(0) << "ERROR: AtomCount " << count << " is not equal to number of atoms in molecule " << j << "!" << endl;
    1742                                                                                 x.Zero();
    1743                                                                                 y.Zero();
    1744                                                                                 y.x[abs(axis)-1] = mol->cell_size[(abs(axis) == 2) ? 2 : ((abs(axis) == 3) ? 5 : 0)] * abs(axis)/axis; // last term is for sign, first is for magnitude
    1745                                                                                 for (int i=1;i<faktor;i++) {    // then add this list with respective translation factor times
    1746                                                                                         x.AddVector(&y); // per factor one cell width further
    1747                                                                                         for (int k=count;k--;) { // go through every atom of the original cell
    1748                                                                                                 first = new atom(); // create a new body
    1749                                                                                                 first->x.CopyVector(vectors[k]);        // use coordinate of original atom
    1750                                                                                                 first->x.AddVector(&x);                 // translate the coordinates
    1751                                                                                                 first->type = Elements[k];      // insert original element
    1752                                                                                                 mol->AddAtom(first);                            // and add to the molecule (which increments ElementsInMolecule, AtomCount, ...)
    1753                                                                                         }
    1754                                                                                 }
    1755                                                                                 // free memory
    1756                                                                                 delete[](Elements);
    1757                                                                                 delete[](vectors);
    1758                                                                                 // correct cell size
    1759                                                                                 if (axis < 0) { // if sign was negative, we have to translate everything
    1760                                                                                         x.Zero();
    1761                                                                                         x.AddVector(&y);
    1762                                                                                         x.Scale(-(faktor-1));
    1763                                                                                         mol->Translate(&x);
    1764                                                                                 }
    1765                                                                                 mol->cell_size[(abs(axis) == 2) ? 2 : ((abs(axis) == 3) ? 5 : 0)] *= faktor;
    1766                                                                         }
    1767                                                                 }
    1768                                                         }
    1769                                                         break;
    1770                                                 default:        // no match? Step on
    1771                                                         if ((argptr < argc) && (argv[argptr][0] != '-')) // if it started with a '-' we've already made a step!
    1772                                                                 argptr++;
    1773                                                         break;
    1774                                         }
    1775                                 }
    1776                         } else argptr++;
    1777                 } while (argptr < argc);
    1778                 if (SaveFlag)
    1779                         SaveConfig(ConfigFileName, &configuration, periode, molecules);
    1780         } else {        // no arguments, hence scan the elements db
    1781                 if (periode->LoadPeriodentafel(configuration.databasepath))
    1782                         cout << Verbose(0) << "Element list loaded successfully." << endl;
    1783                 else
    1784                         cout << Verbose(0) << "Element list loading failed." << endl;
    1785                 configuration.RetrieveConfigPathAndName("main_pcp_linux");
    1786         }
    1787         return(ExitFlag);
     1625              break;
     1626            case 'r':
     1627              ExitFlag = 1;
     1628              SaveFlag = true;
     1629              cout << Verbose(1) << "Converting config file from supposed old to new syntax." << endl;
     1630              break;
     1631            case 'F':
     1632            case 'f':
     1633              ExitFlag = 1;
     1634              if ((argptr+1 >= argc) || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1]))) {
     1635                ExitFlag = 255;
     1636                cerr << "Not enough or invalid arguments for fragmentation: -f <max. bond distance> <bond order>" << endl;
     1637              } else {
     1638                cout << "Fragmenting molecule with bond distance " << argv[argptr] << " angstroem, order of " << argv[argptr+1] << "." << endl;
     1639                cout << Verbose(0) << "Creating connection matrix..." << endl;
     1640                start = clock();
     1641                mol->CreateAdjacencyList((ofstream *)&cout, atof(argv[argptr++]), configuration.GetIsAngstroem());
     1642                cout << Verbose(0) << "Fragmenting molecule with current connection matrix ..." << endl;
     1643                if (mol->first->next != mol->last) {
     1644                  ExitFlag = mol->FragmentMolecule((ofstream *)&cout, atoi(argv[argptr]), &configuration);
     1645                }
     1646                end = clock();
     1647                cout << Verbose(0) << "Clocks for this operation: " << (end-start) << ", time: " << ((double)(end-start)/CLOCKS_PER_SEC) << "s." << endl;
     1648                argptr+=2;
     1649              }
     1650              break;
     1651            case 'm':
     1652              ExitFlag = 1;
     1653              j = atoi(argv[argptr++]);
     1654              if ((j<0) || (j>1)) {
     1655                cerr << Verbose(1) << "ERROR: Argument of '-m' should be either 0 for no-rotate or 1 for rotate." << endl;
     1656                j = 0;
     1657              }
     1658              if (j) {
     1659                SaveFlag = true;
     1660                cout << Verbose(0) << "Converting to prinicipal axis system." << endl;
     1661              } else
     1662                cout << Verbose(0) << "Evaluating prinicipal axis." << endl;
     1663              mol->PrincipalAxisSystem((ofstream *)&cout, (bool)j);
     1664              break;
     1665            case 'o':
     1666              ExitFlag = 1;
     1667              if ((argptr >= argc) || (argv[argptr][0] == '-')){
     1668                ExitFlag = 255;
     1669                cerr << "Not enough or invalid arguments given for convex envelope: -o <tecplot output file>" << endl;
     1670              } else {
     1671                cout << Verbose(0) << "Evaluating volume of the convex envelope.";
     1672                cout << Verbose(1) << "Storing tecplot data in " << argv[argptr] << "." << endl;
     1673                VolumeOfConvexEnvelope((ofstream *)&cout, argv[argptr], &configuration, NULL, mol);
     1674                argptr+=1;
     1675              }
     1676              break;
     1677            case 'U':
     1678              ExitFlag = 1;
     1679              if ((argptr+1 >= argc) || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1])) ) {
     1680                ExitFlag = 255;
     1681                cerr << "Not enough or invalid arguments given for suspension with specified volume: -U <volume> <density>" << endl;
     1682                volume = -1; // for case 'u': don't print error again
     1683              } else {
     1684                volume = atof(argv[argptr++]);
     1685                cout << Verbose(0) << "Using " << volume << " angstrom^3 as the volume instead of convex envelope one's." << endl;
     1686              }
     1687            case 'u':
     1688              ExitFlag = 1;
     1689              if ((argptr >= argc) || (!IsValidNumber(argv[argptr])) ) {
     1690                if (volume != -1)
     1691                  ExitFlag = 255;
     1692                  cerr << "Not enough arguments given for suspension: -u <density>" << endl;
     1693              } else {
     1694                double density;
     1695                SaveFlag = true;
     1696                cout << Verbose(0) << "Evaluating necessary cell volume for a cluster suspended in water.";
     1697                density = atof(argv[argptr++]);
     1698                if (density < 1.0) {
     1699                  cerr << Verbose(0) << "Density must be greater than 1.0g/cm^3 !" << endl;
     1700                  density = 1.3;
     1701                }
     1702//                for(int i=0;i<NDIM;i++) {
     1703//                  repetition[i] = atoi(argv[argptr++]);
     1704//                  if (repetition[i] < 1)
     1705//                    cerr << Verbose(0) << "ERROR: repetition value must be greater 1!" << endl;
     1706//                  repetition[i] = 1;
     1707//                }
     1708                PrepareClustersinWater((ofstream *)&cout, &configuration, mol, volume, density);  // if volume == 0, will calculate from ConvexEnvelope
     1709              }
     1710              break;
     1711            case 'd':
     1712              ExitFlag = 1;
     1713              if ((argptr+2 >= argc) || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1])) || (!IsValidNumber(argv[argptr+2])) ) {
     1714                ExitFlag = 255;
     1715                cerr << "Not enough or invalid arguments given for repeating cells: -d <repeat_x> <repeat_y> <repeat_z>" << endl;
     1716              } else {
     1717                SaveFlag = true;
     1718                for (int axis = 1; axis <= NDIM; axis++) {
     1719                  int faktor = atoi(argv[argptr++]);
     1720                  int count;
     1721                  element ** Elements;
     1722                  Vector ** vectors;
     1723                  if (faktor < 1) {
     1724                    cerr << Verbose(0) << "ERROR: Repetition faktor mus be greater than 1!" << endl;
     1725                    faktor = 1;
     1726                  }
     1727                  mol->CountAtoms((ofstream *)&cout);  // recount atoms
     1728                  if (mol->AtomCount != 0) {  // if there is more than none
     1729                    count = mol->AtomCount;  // is changed becausing of adding, thus has to be stored away beforehand
     1730                    Elements = new element *[count];
     1731                    vectors = new Vector *[count];
     1732                    j = 0;
     1733                    first = mol->start;
     1734                    while (first->next != mol->end) {  // make a list of all atoms with coordinates and element
     1735                      first = first->next;
     1736                      Elements[j] = first->type;
     1737                      vectors[j] = &first->x;
     1738                      j++;
     1739                    }
     1740                    if (count != j)
     1741                      cout << Verbose(0) << "ERROR: AtomCount " << count << " is not equal to number of atoms in molecule " << j << "!" << endl;
     1742                    x.Zero();
     1743                    y.Zero();
     1744                    y.x[abs(axis)-1] = mol->cell_size[(abs(axis) == 2) ? 2 : ((abs(axis) == 3) ? 5 : 0)] * abs(axis)/axis; // last term is for sign, first is for magnitude
     1745                    for (int i=1;i<faktor;i++) {  // then add this list with respective translation factor times
     1746                      x.AddVector(&y); // per factor one cell width further
     1747                      for (int k=count;k--;) { // go through every atom of the original cell
     1748                        first = new atom(); // create a new body
     1749                        first->x.CopyVector(vectors[k]);  // use coordinate of original atom
     1750                        first->x.AddVector(&x);      // translate the coordinates
     1751                        first->type = Elements[k];  // insert original element
     1752                        mol->AddAtom(first);        // and add to the molecule (which increments ElementsInMolecule, AtomCount, ...)
     1753                      }
     1754                    }
     1755                    // free memory
     1756                    delete[](Elements);
     1757                    delete[](vectors);
     1758                    // correct cell size
     1759                    if (axis < 0) { // if sign was negative, we have to translate everything
     1760                      x.Zero();
     1761                      x.AddVector(&y);
     1762                      x.Scale(-(faktor-1));
     1763                      mol->Translate(&x);
     1764                    }
     1765                    mol->cell_size[(abs(axis) == 2) ? 2 : ((abs(axis) == 3) ? 5 : 0)] *= faktor;
     1766                  }
     1767                }
     1768              }
     1769              break;
     1770            default:  // no match? Step on
     1771              if ((argptr < argc) && (argv[argptr][0] != '-')) // if it started with a '-' we've already made a step!
     1772                argptr++;
     1773              break;
     1774          }
     1775        }
     1776      } else argptr++;
     1777    } while (argptr < argc);
     1778    if (SaveFlag)
     1779      SaveConfig(ConfigFileName, &configuration, periode, molecules);
     1780  } else {  // no arguments, hence scan the elements db
     1781    if (periode->LoadPeriodentafel(configuration.databasepath))
     1782      cout << Verbose(0) << "Element list loaded successfully." << endl;
     1783    else
     1784      cout << Verbose(0) << "Element list loading failed." << endl;
     1785    configuration.RetrieveConfigPathAndName("main_pcp_linux");
     1786  }
     1787  return(ExitFlag);
    17881788};
    17891789
    … …  
    17921792int main(int argc, char **argv)
    17931793{
    1794         periodentafel *periode = new periodentafel; // and a period table of all elements
    1795         MoleculeListClass *molecules = new MoleculeListClass;  // list of all molecules
    1796         molecule *mol = NULL;
    1797         config configuration;
    1798         char choice;    // menu choice char
    1799         Vector x,y,z,n; // coordinates for absolute point in cell volume
    1800         ifstream test;
    1801         ofstream output;
    1802         string line;
    1803         char *ConfigFileName = NULL;
    1804         int j, count;
    1805 
    1806         // =========================== PARSE COMMAND LINE OPTIONS ====================================
    1807         j = ParseCommandLineOptions(argc, argv, molecules, periode, configuration, ConfigFileName);
    1808         switch(j) {
    1809           case 0:  // something went wrong
    1810             delete(molecules); // also free's all molecules contained
    1811             delete(periode);
    1812             return j;
    1813             break;
    1814           case 1:  // just for -v and -h options
    1815             delete(molecules); // also free's all molecules contained
    1816             delete(periode);
    1817             return 0;
    1818             break;
    1819           default:
    1820             break;
    1821         }
    1822 
    1823         // General stuff
    1824         if (molecules->ListOfMolecules.size() == 0) {
     1794  periodentafel *periode = new periodentafel; // and a period table of all elements
     1795  MoleculeListClass *molecules = new MoleculeListClass;  // list of all molecules
     1796  molecule *mol = NULL;
     1797  config configuration;
     1798  char choice;  // menu choice char
     1799  Vector x,y,z,n;  // coordinates for absolute point in cell volume
     1800  ifstream test;
     1801  ofstream output;
     1802  string line;
     1803  char *ConfigFileName = NULL;
     1804  int j, count;
     1805
     1806  // =========================== PARSE COMMAND LINE OPTIONS ====================================
     1807  j = ParseCommandLineOptions(argc, argv, molecules, periode, configuration, ConfigFileName);
     1808  switch(j) {
     1809    case 0:  // something went wrong
     1810      delete(molecules); // also free's all molecules contained
     1811      delete(periode);
     1812      return j;
     1813      break;
     1814    case 1:  // just for -v and -h options
     1815      delete(molecules); // also free's all molecules contained
     1816      delete(periode);
     1817      return 0;
     1818      break;
     1819    default:
     1820      break;
     1821  }
     1822
     1823  // General stuff
     1824  if (molecules->ListOfMolecules.size() == 0) {
    18251825    mol = new molecule(periode);
    18261826    if (mol->cell_size[0] == 0.) {
    … …  
    18321832    }
    18331833    molecules->insert(mol);
    1834         }
    1835 
    1836         // =========================== START INTERACTIVE SESSION ====================================
    1837 
    1838         // now the main construction loop
    1839         cout << Verbose(0) << endl << "Now comes the real construction..." << endl;
    1840         do {
    1841                 cout << Verbose(0) << endl << endl;
    1842                 cout << Verbose(0) << "============Molecule list=======================" << endl;
    1843                 molecules->Enumerate((ofstream *)&cout);
    1844                 cout << Verbose(0) << "============Menu===============================" << endl;
     1834  }
     1835
     1836  // =========================== START INTERACTIVE SESSION ====================================
     1837
     1838  // now the main construction loop
     1839  cout << Verbose(0) << endl << "Now comes the real construction..." << endl;
     1840  do {
     1841    cout << Verbose(0) << endl << endl;
     1842    cout << Verbose(0) << "============Molecule list=======================" << endl;
     1843    molecules->Enumerate((ofstream *)&cout);
     1844    cout << Verbose(0) << "============Menu===============================" << endl;
    18451845    cout << Verbose(0) << "a - set molecule (in)active" << endl;
    18461846    cout << Verbose(0) << "e - edit new molecules" << endl;
    … …  
    18581858    cin >> choice;
    18591859
    1860                 switch (choice) {
    1861                         case 'a':  // (in)activate molecule
     1860    switch (choice) {
     1861      case 'a':  // (in)activate molecule
    18621862        {
    18631863          cout << "Enter index of molecule: ";
    … …  
    18691869            (*ListRunner)->ActiveFlag = !(*ListRunner)->ActiveFlag;
    18701870        }
    1871                           break;
    1872 
    1873                         case 'c': // edit each field of the configuration
    1874                         configuration.Edit();
    1875                         break;
     1871        break;
     1872
     1873      case 'c': // edit each field of the configuration
     1874      configuration.Edit();
     1875      break;
    18761876
    18771877      case 'e': // create molecule
    … …  
    18911891        break;
    18921892
    1893                         case 'q': // quit
    1894                                 break;
    1895 
    1896                         case 's': // save to config file
    1897                                 SaveConfig(ConfigFileName, &configuration, periode, molecules);
    1898                                 break;
    1899 
    1900                         case 'T':
    1901                                 testroutine(molecules);
    1902                                 break;
    1903 
    1904                         default:
    1905                           break;
    1906                 };
    1907         } while (choice != 'q');
    1908 
    1909         // save element data base
    1910         if (periode->StorePeriodentafel(configuration.databasepath)) //ElementsFileName
    1911                 cout << Verbose(0) << "Saving of elements.db successful." << endl;
    1912         else
    1913                 cout << Verbose(0) << "Saving of elements.db failed." << endl;
    1914 
    1915         delete(molecules); // also free's all molecules contained
    1916         delete(periode);
    1917         return (0);
     1893      case 'q': // quit
     1894        break;
     1895
     1896      case 's': // save to config file
     1897        SaveConfig(ConfigFileName, &configuration, periode, molecules);
     1898        break;
     1899
     1900      case 'T':
     1901        testroutine(molecules);
     1902        break;
     1903
     1904      default:
     1905        break;
     1906    };
     1907  } while (choice != 'q');
     1908
     1909  // save element data base
     1910  if (periode->StorePeriodentafel(configuration.databasepath)) //ElementsFileName
     1911    cout << Verbose(0) << "Saving of elements.db successful." << endl;
     1912  else
     1913    cout << Verbose(0) << "Saving of elements.db failed." << endl;
     1914
     1915  delete(molecules); // also free's all molecules contained
     1916  delete(periode);
     1917  return (0);
    19181918}
    19191919
Note: See TracChangeset for help on using the changeset viewer.