Changeset a98603 for src/builder.cpp


Ignore:
Timestamp:
Feb 9, 2009, 2:18:13 PM (16 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:
5bc4d0
Parents:
674220 (diff), cc2ee5 (diff)
Note: this is a merge changeset, the changes displayed below correspond to the merge itself.
Use the (diff) links above to see all the changes relative to each parent.
Message:

Merge ../espack3

File:
1 edited

Legend:

Unmodified
Added
Removed

  • src/builder.cpp

    • Property mode changed from 100644 to 100755
    r674220 ra98603  
    11/** \file builder.cpp
    2  * 
     2 *
    33 * By stating absolute positions or binding angles and distances atomic positions of a molecule can be constructed.
    44 * The output is the complete configuration file for PCP for direct use.
    … …  
    66 * -# Atomic data is retrieved from a file, if not found requested and stored there for later re-use
    77 * -# step-by-step construction of the molecule beginning either at a centre of with a certain atom
    8  *   
     8 *
    99 */
    1010
    1111/*! \mainpage Molecuilder - a molecular set builder
    12  * 
     12 *
    1313 * This introductory shall briefly make aquainted with the program, helping in installing and a first run.
    14  * 
     14 *
    1515 * \section about About the Program
    16  * 
     16 *
    1717 *  Molecuilder is a short program, written in C++, that enables the construction of a coordinate set for the
    1818 *  atoms making up an molecule by the successive statement of binding angles and distances and referencing to
    1919 *  already constructed atoms.
    20  * 
     20 *
    2121 *  A configuration file may be written that is compatible to the format used by PCP - a parallel Car-Parrinello
    2222 *  molecular dynamics implementation.
    23  * 
     23 *
    2424 * \section install Installation
    25  * 
     25 *
    2626 *  Installation should without problems succeed as follows:
    2727 *  -# ./configure (or: mkdir build;mkdir run;cd build; ../configure --bindir=../run)
    2828 *  -# make
    2929 *  -# make install
    30  * 
     30 *
    3131 *  Further useful commands are
    3232 *  -# 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 
    34  * 
     33 *  -# make doxygen-doc: Creates these html pages out of the documented source
     34 *
    3535 * \section run Running
    36  * 
     36 *
    3737 *  The program can be executed by running: ./molecuilder
    38  * 
     38 *
    3939 *  Note, that it uses a database, called "elements.db", in the executable's directory. If the file is not found,
    4040 *  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. 
    42  * 
     41 *  later re-execution.
     42 *
    4343 * \section ref References
    44  * 
     44 *
    4545 *  For the special configuration file format, see the documentation of pcp.
    46  * 
     46 *
    4747 */
    4848
    … …  
    8080  cout << Verbose(0) << "INPUT: ";
    8181  cin >> choice;
    82  
     82
    8383  switch (choice) {
    8484      case 'a': // absolute coordinates of atom
    … …  
    8989        mol->AddAtom(first);  // add to molecule
    9090        break;
    91        
     91
    9292      case 'b': // relative coordinates of atom wrt to reference point
    9393        first = new atom;
    … …  
    105105        mol->AddAtom(first);  // add to molecule
    106106        break;
    107        
     107
    108108      case 'c': // relative coordinates of atom wrt to already placed atom
    109109        first = new atom;
    … …  
    111111        do {
    112112          if (!valid) cout << Verbose(0) << "Resulting position out of cell." << endl;
    113           second = mol->AskAtom("Enter atom number: ");               
     113          second = mol->AskAtom("Enter atom number: ");
    114114          cout << Verbose(0) << "Enter relative coordinates." << endl;
    115115          first->x.AskPosition(mol->cell_size, false);
    … …  
    121121        mol->AddAtom(first);  // add to molecule
    122122        break;
    123      
     123
    124124      case 'd': // two atoms, two angles and a distance
    125125        first = new atom;
    … …  
    152152          x.Copyvector(&fourth->x);
    153153          x.SubtractVector(&third->x);
    154          
     154
    155155          if (!z.SolveSystem(&x,&y,&n, b, c, a)) {
    156156            cout << Verbose(0) << "Failure solving self-dependent linear system!" << endl;
    … …  
    167167          cout << "x: ",
    168168          x.Output((ofstream *)&cout);
    169           cout << endl;         
     169          cout << endl;
    170170          z.MakeNormalVector(&second->x,&third->x,&fourth->x);
    171171          cout << "z: ",
    172172          z.Output((ofstream *)&cout);
    173           cout << endl;         
     173          cout << endl;
    174174          y.MakeNormalVector(&x,&z);
    175175          cout << "y: ",
    176176          y.Output((ofstream *)&cout);
    177           cout << endl;         
    178          
     177          cout << endl;
     178
    179179          // rotate vector around first angle
    180180          first->x.CopyVector(&x);
    … …  
    182182          cout << "Rotated vector: ",
    183183          first->x.Output((ofstream *)&cout);
    184           cout << endl;         
     184          cout << endl;
    185185          // remove the projection onto the rotation plane of the second angle
    186186          n.CopyVector(&y);
    … …  
    188188          cout << "N1: ",
    189189          n.Output((ofstream *)&cout);
    190           cout << endl;         
     190          cout << endl;
    191191          first->x.SubtractVector(&n);
    192192          cout << "Subtracted vector: ",
    193193          first->x.Output((ofstream *)&cout);
    194           cout << endl;         
     194          cout << endl;
    195195          n.CopyVector(&z);
    196196          n.Scale(first->x.Projection(&z));
    197197          cout << "N2: ",
    198198          n.Output((ofstream *)&cout);
    199           cout << endl;         
     199          cout << endl;
    200200          first->x.SubtractVector(&n);
    201201          cout << "2nd subtracted vector: ",
    202202          first->x.Output((ofstream *)&cout);
    203           cout << endl;         
    204          
     203          cout << endl;
     204
    205205          // rotate another vector around second angle
    206206          n.CopyVector(&y);
    … …  
    208208          cout << "2nd Rotated vector: ",
    209209          n.Output((ofstream *)&cout);
    210           cout << endl;         
    211          
     210          cout << endl;
     211
    212212          // add the two linear independent vectors
    213213          first->x.AddVector(&n);
    214           first->x.Normalize();         
     214          first->x.Normalize();
    215215          first->x.Scale(a);
    216216          first->x.AddVector(&second->x);
    217          
     217
    218218          cout << Verbose(0) << "resulting coordinates: ";
    219219          first->x.Output((ofstream *)&cout);
    … …  
    241241        } while ((j != -1) && (i<128));
    242242        if (i >= 2) {
    243           first->x.LSQdistance(atoms, i);             
     243          first->x.LSQdistance(atoms, i);
    244244
    245245          first->x.Output((ofstream *)&cout);
    … …  
    259259static void CenterAtoms(molecule *mol)
    260260{
    261   Vector x, y;
     261  Vector x, y, helper;
    262262  char choice;  // menu choice char
    263  
     263
    264264  cout << Verbose(0) << "===========CENTER ATOMS=========================" << endl;
    265265  cout << Verbose(0) << " a - on origin" << endl;
    … …  
    271271  cout << Verbose(0) << "INPUT: ";
    272272  cin >> choice;
    273  
     273
    274274  switch (choice) {
    275275    default:
    … …  
    292292      mol->CenterEdge((ofstream *)&cout, &x);  // make every coordinate positive
    293293      mol->Translate(&y); // translate by boundary
    294       mol->SetBoxDimension(&(x+y*2));  // update Box of atoms by boundary
     294      helper.CopyVector(&y);
     295      helper.Scale(2.);
     296      helper.AddVector(&x);
     297      mol->SetBoxDimension(&helper);  // update Box of atoms by boundary
    295298      break;
    296299    case 'd':
    … …  
    327330  cout << Verbose(0) << "INPUT: ";
    328331  cin >> choice;
    329  
     332
    330333  switch (choice) {
    331334    default:
    … …  
    346349      second = mol->AskAtom("Enter second atom: ");
    347350
    348       n.CopyVector((const Vector *)&first->x); 
    349       n.SubtractVector((const Vector *)&second->x); 
     351      n.CopyVector((const Vector *)&first->x);
     352      n.SubtractVector((const Vector *)&second->x);
    350353      n.Normalize();
    351       break;       
     354      break;
    352355    case 'd':
    353356      char shorthand[4];
    … …  
    363366        x.x[i] = gsl_vector_get(param.x,i);
    364367        n.x[i] = gsl_vector_get(param.x,i+NDIM);
    365       } 
     368      }
    366369      gsl_vector_free(param.x);
    367370      cout << Verbose(0) << "Offset vector: ";
    … …  
    369372      cout << Verbose(0) << endl;
    370373      n.Normalize();
    371       break;       
     374      break;
    372375  };
    373376  cout << Verbose(0) << "Alignment vector: ";
    … …  
    385388  Vector n;
    386389  char choice;  // menu choice char
    387  
     390
    388391  cout << Verbose(0) << "===========MIRROR ATOMS=========================" << endl;
    389392  cout << Verbose(0) << " a - state three atoms defining mirror plane" << endl;
    … …  
    394397  cout << Verbose(0) << "INPUT: ";
    395398  cin >> choice;
    396  
     399
    397400  switch (choice) {
    398401    default:
    … …  
    413416      second = mol->AskAtom("Enter second atom: ");
    414417
    415       n.CopyVector((const Vector *)&first->x); 
    416       n.SubtractVector((const Vector *)&second->x); 
     418      n.CopyVector((const Vector *)&first->x);
     419      n.SubtractVector((const Vector *)&second->x);
    417420      n.Normalize();
    418       break;         
     421      break;
    419422  };
    420423  cout << Verbose(0) << "Normal vector: ";
    … …  
    433436  double tmp1, tmp2;
    434437  char choice;  // menu choice char
    435  
     438
    436439  cout << Verbose(0) << "===========REMOVE ATOMS=========================" << endl;
    437440  cout << Verbose(0) << " a - state atom for removal by number" << endl;
    … …  
    442445  cout << Verbose(0) << "INPUT: ";
    443446  cin >> choice;
    444  
     447
    445448  switch (choice) {
    446449    default:
    … …  
    475478          mol->RemoveAtom(first);
    476479      }
    477       break;         
     480      break;
    478481  };
    479482  //mol->Output((ofstream *)&cout);
    … …  
    492495  int Z;
    493496  char choice;  // menu choice char
    494  
     497
    495498  cout << Verbose(0) << "===========MEASURE ATOMS=========================" << endl;
    496499  cout << Verbose(0) << " a - calculate bond length between one atom and all others" << endl;
    … …  
    514517      for (int i=MAX_ELEMENTS;i--;)
    515518        min[i] = 0.;
    516        
    517       second = mol->start;   
     519
     520      second = mol->start;
    518521      while ((second->next != mol->end)) {
    519522        second = second->next; // advance
    … …  
    526529        }
    527530        if ((tmp1 != 0.) && ((min[Z] == 0.) || (tmp1 < min[Z]))) min[Z] = tmp1;
    528         //cout << Verbose(0) << "Bond length between Atom " << first->nr << " and " << second->nr << ": " << tmp1 << " a.u." << endl;         
     531        //cout << Verbose(0) << "Bond length between Atom " << first->nr << " and " << second->nr << ": " << tmp1 << " a.u." << endl;
    529532      }
    530533      for (int i=MAX_ELEMENTS;i--;)
    531534        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;
    532535      break;
    533      
     536
    534537    case 'b':
    535538      first = mol->AskAtom("Enter first atom: ");
    … …  
    556559      y.SubtractVector((const Vector *)&second->x);
    557560      cout << Verbose(0) << "Bond angle between first atom Nr." << first->nr << ", central atom Nr." << second->nr << " and last atom Nr." << third->nr << ": ";
    558       cout << Verbose(0) << (acos(x.ScalarProduct((const Vector *)&y)/(y.Norm()*x.Norm()))/M_PI*180.) << " degrees" << endl;         
     561      cout << Verbose(0) << (acos(x.ScalarProduct((const Vector *)&y)/(y.Norm()*x.Norm()))/M_PI*180.) << " degrees" << endl;
    559562      break;
    560563    case 'd':
    … …  
    600603  int Order1;
    601604  clock_t start, end;
    602  
     605
    603606  cout << Verbose(0) << "Fragmenting molecule with current connection matrix ..." << endl;
    604607  cout << Verbose(0) << "What's the desired bond order: ";
    … …  
    609612    end = clock();
    610613    cout << Verbose(0) << "Clocks for this operation: " << (end-start) << ", time: " << ((double)(end-start)/CLOCKS_PER_SEC) << "s." << endl;
    611   } else 
     614  } else
    612615    cout << Verbose(0) << "Connection matrix has not yet been generated!" << endl;
    613616};
    … …  
    623626  atom *Walker = mol->start;
    624627  int i, comp, counter=0;
    625  
     628
    626629  // generate some KeySets
    627630  cout << "Generating KeySets." << endl;
    … …  
    637640  cout << "Testing insertion of already present item in KeySets." << endl;
    638641  KeySetTestPair test;
    639   test = TestSets[mol->AtomCount-1].insert(Walker->nr); 
     642  test = TestSets[mol->AtomCount-1].insert(Walker->nr);
    640643  if (test.second) {
    641644    cout << Verbose(1) << "Insertion worked?!" << endl;
    … …  
    646649  TestSets[mol->AtomCount].insert(mol->end->previous->previous->previous->nr);
    647650
    648   // constructing Graph structure 
     651  // constructing Graph structure
    649652  cout << "Generating Subgraph class." << endl;
    650653  Graph Subgraphs;
    … …  
    657660  cout << "Testing insertion of already present item in Subgraph." << endl;
    658661  GraphTestPair test2;
    659   test2 = Subgraphs.insert(GraphPair (TestSets[mol->AtomCount],pair<int, double>(counter++, 1.))); 
     662  test2 = Subgraphs.insert(GraphPair (TestSets[mol->AtomCount],pair<int, double>(counter++, 1.)));
    660663  if (test2.second) {
    661664    cout << Verbose(1) << "Insertion worked?!" << endl;
    … …  
    663666    cout << Verbose(1) << "Insertion rejected: Present object is " << (*(test2.first)).second.first << "." << endl;
    664667  }
    665  
     668
    666669  // show graphs
    667670  cout << "Showing Subgraph's contents, checking that it's sorted." << endl;
    … …  
    674677      if ((*key) > comp)
    675678        cout << (*key) << " ";
    676       else 
     679      else
    677680        cout << (*key) << "! ";
    678681      comp = (*key);
    … …  
    716719  else
    717720    cout << "failed." << endl;
    718  
     721
    719722  // and save to xyz file
    720723  if (ConfigFileName != NULL) {
    … …  
    727730    strcat(filename, ".xyz");
    728731    output.open(filename, ios::trunc);
    729   } 
     732  }
    730733  cout << Verbose(0) << "Saving of XYZ file ";
    731734  if (mol->MDSteps <= 1) {
    … …  
    742745  output.close();
    743746  output.clear();
    744  
     747
    745748  // and save as MPQC configuration
    746749  if (ConfigFileName != NULL)
    … …  
    753756  else
    754757    cout << "failed." << endl;
    755  
     758
    756759  if (!strcmp(configuration->configpath, configuration->GetDefaultPath())) {
    757760    cerr << "WARNING: config is found under different path then stated in config file::defaultpath!" << endl;
    … …  
    785788  int argptr;
    786789  PathToDatabases = LocalPath;
    787  
     790
    788791  if (argc > 1) { // config file specified as option
    789792    // 1. : Parse options that just set variables or print help
    … …  
    798801          case '?':
    799802            cout << "MoleCuilder suite" << endl << "==================" << endl << endl;
    800             cout << "Usage: " << argv[0] << "[config file] [-{acefpsthH?vfrp}] [further arguments]" << endl; 
     803            cout << "Usage: " << argv[0] << "[config file] [-{acefpsthH?vfrp}] [further arguments]" << endl;
    801804            cout << "or simply " << argv[0] << " without arguments for interactive session." << endl;
    802805            cout << "\t-a Z x1 x2 x3\tAdd new atom of element Z at coordinates (x1,x2,x3)." << endl;
     806            cout << "\t-A <source>\tCreate adjacency list from bonds parsed from 'dbond'-style file." <<endl;
    803807            cout << "\t-b x1 x2 x3\tCenter atoms in domain with given edge lengths of (x1,x2,x3)." << endl;
    804808            cout << "\t-c x1 x2 x3\tCenter atoms in domain with a minimum distance to boundary of (x1,x2,x3)." << endl;
    … …  
    812816            cout << "\t-m <0/1>\tCalculate (0)/ Align in(1) PAS with greatest EV along z axis." << endl;
    813817            cout << "\t-n\tFast parsing (i.e. no trajectories are looked for)." << endl;
    814             cout << "\t-o\tGet volume of the convex envelope (and store to tecplot file)." << endl;
     818            cout << "\t-N\tGet non-convex-envelope." << endl;
     819            cout << "\t-o <out>\tGet volume of the convex envelope (and store to tecplot file)." << endl;
    815820            cout << "\t-p <file>\tParse given xyz file and create raw config file from it." << endl;
    816821            cout << "\t-P <file>\tParse given forces file and append as an MD step to config file via Verlet." << endl;
    817822            cout << "\t-r\t\tConvert file from an old pcp syntax." << endl;
    818823            cout << "\t-t x1 x2 x3\tTranslate all atoms by this vector (x1,x2,x3)." << endl;
    819             cout << "\t-T <file> Store temperatures from the config file in <file>." << endl; 
     824            cout << "\t-T <file> Store temperatures from the config file in <file>." << endl;
    820825            cout << "\t-s x1 x2 x3\tScale all atom coordinates by this vector (x1,x2,x3)." << endl;
    821             cout << "\t-u rho\tsuspend in water solution and output necessary cell lengths, average density rho and repetition." << endl; 
     826            cout << "\t-u rho\tsuspend in water solution and output necessary cell lengths, average density rho and repetition." << endl;
    822827            cout << "\t-v/-V\t\tGives version information." << endl;
    823828            cout << "Note: config files must not begin with '-' !" << endl;
    … …  
    854859        argptr++;
    855860    } while (argptr < argc);
    856    
     861
    857862    // 2. Parse the element database
    858863    if (periode->LoadPeriodentafel(PathToDatabases)) {
    … …  
    863868      return 1;
    864869    }
    865    
     870
    866871    // 3. Find config file name and parse if possible
    867872    if (argv[1][0] != '-') {
    … …  
    902907    } else
    903908      config_present = absent;
    904    
     909
    905910    // 4. parse again through options, now for those depending on elements db and config presence
    906911    argptr = 1;
    … …  
    946951                    config_present = present;
    947952                } else
    948                   cerr << Verbose(1) << "Could not find the specified element." << endl; 
     953                  cerr << Verbose(1) << "Could not find the specified element." << endl;
    949954                argptr+=4;
    950955              }
    … …  
    956961        if (config_present == present) {
    957962          switch(argv[argptr-1][1]) {
    958             case 'D': 
     963            case 'D':
    959964              ExitFlag = 1;
    960965              {
    … …  
    10021007              }
    10031008              break;
     1009                                                case 'A':
     1010                                                        ExitFlag = 1;
     1011                                                        if ((argptr >= argc) || (!IsValidNumber(argv[argptr])) || (argv[argptr][0] == '-')) {
     1012                                                                ExitFlag =255;
     1013                                                                cerr << "Missing source file for bonds in molecule: -A <bond sourcefile>" << endl;
     1014                                                        } else {
     1015                                                                cout << "Parsing bonds from " << argv[argptr] << "." << endl;
     1016                                                                ifstream *input = new ifstream(argv[argptr]);
     1017                                                                mol->CreateAdjacencyList2((ofstream *)&cout, input);
     1018                                                                input->close();
     1019                                                        }
     1020                                                        break;
     1021                                                case 'N':
     1022                                                        ExitFlag = 1;
     1023                                                        if ((argptr+1 >= argc) || (argv[argptr+1][0] == '-')){
     1024                                                                ExitFlag = 255;
     1025                                                                cerr << "Not enough or invalid arguments given for non-convex envelope: -o <radius> <tecplot output file>" << endl;
     1026                                                        } else {
     1027                                                                cout << Verbose(0) << "Evaluating npn-convex envelope.";
     1028                                                                string TempName(argv[argptr+1]);
     1029                                                                TempName.append(".r3d");
     1030                                                                ofstream *output = new ofstream(TempName.c_str(), ios::trunc);
     1031                                                                cout << Verbose(1) << "Using rolling ball of radius " << atof(argv[argptr]) << " and storing tecplot data in " << argv[argptr+1] << "." << endl;
     1032                                                                Find_non_convex_border((ofstream *)&cout, output, mol, argv[argptr+1], atof(argv[argptr]));
     1033                                                                output->close();
     1034                                                                delete(output);
     1035                                                                argptr+=2;
     1036                                                        }
     1037                                                        break;
    10041038            case 'T':
    10051039              ExitFlag = 1;
    … …  
    11681202            case 'o':
    11691203              ExitFlag = 1;
    1170               if ((argptr >= argc) || (argv[argptr][0] == '-')) {
     1204              if ((argptr >= argc) || (argv[argptr][0] == '-')){
    11711205                ExitFlag = 255;
    11721206                cerr << "Not enough or invalid arguments given for convex envelope: -o <tecplot output file>" << endl;
    11731207              } else {
    11741208                cout << Verbose(0) << "Evaluating volume of the convex envelope.";
     1209                ofstream *output = new ofstream(argv[argptr], ios::trunc);
     1210                //$$$
    11751211                cout << Verbose(1) << "Storing tecplot data in " << argv[argptr] << "." << endl;
    1176                 ofstream *output = new ofstream(argv[argptr], ios::trunc);
    11771212                VolumeOfConvexEnvelope((ofstream *)&cout, output, &configuration, NULL, mol);
    11781213                output->close();
    … …  
    12691304                      mol->Translate(&x);
    12701305                    }
    1271                     mol->cell_size[(abs(axis) == 2) ? 2 : ((abs(axis) == 3) ? 5 : 0)] *= faktor; 
     1306                    mol->cell_size[(abs(axis) == 2) ? 2 : ((abs(axis) == 3) ? 5 : 0)] *= faktor;
    12721307                  }
    12731308                }
    … …  
    13321367  if (j == 1) return 0; // just for -v and -h options
    13331368  if (j) return j;  // something went wrong
    1334  
     1369
    13351370  // General stuff
    13361371  if (mol->cell_size[0] == 0.) {
    … …  
    13461381  // now the main construction loop
    13471382  cout << Verbose(0) << endl << "Now comes the real construction..." << endl;
    1348   do {   
     1383  do {
    13491384    cout << Verbose(0) << endl << endl;
    13501385    cout << Verbose(0) << "============Element list=======================" << endl;
    … …  
    13651400    cout << Verbose(0) << "-----------------------------------------------" << endl;
    13661401    cout << Verbose(0) << "d - duplicate molecule/periodic cell" << endl;
    1367     cout << Verbose(0) << "i - realign molecule" << endl; 
    1368     cout << Verbose(0) << "m - mirror all molecules" << endl; 
     1402    cout << Verbose(0) << "i - realign molecule" << endl;
     1403    cout << Verbose(0) << "m - mirror all molecules" << endl;
    13691404    cout << Verbose(0) << "t - translate molecule by vector" << endl;
    13701405    cout << Verbose(0) << "c - scale by unit transformation" << endl;
    … …  
    13771412    cout << Verbose(0) << "Input: ";
    13781413    cin >> choice;
    1379    
     1414
    13801415    switch (choice) {
    13811416      default:
    13821417      case 'a': // add atom
    13831418        AddAtoms(periode, mol);
    1384         choice = 'a'; 
    1385         break;
    1386      
     1419        choice = 'a';
     1420        break;
     1421
    13871422      case 'b': // scale a bond
    13881423        cout << Verbose(0) << "Scaling bond length between two atoms." << endl;
    … …  
    14001435        }
    14011436        //cout << Verbose(0) << "New coordinates of Atom " << second->nr << " are: ";
    1402         //second->Output(second->type->No, 1, (ofstream *)&cout);       
    1403         break;
    1404 
    1405       case 'c': // unit scaling of the metric 
     1437        //second->Output(second->type->No, 1, (ofstream *)&cout);
     1438        break;
     1439
     1440      case 'c': // unit scaling of the metric
    14061441       cout << Verbose(0) << "Angstroem -> Bohrradius: 1.8897261\t\tBohrradius -> Angstroem: 0.52917721" << endl;
    14071442       cout << Verbose(0) << "Enter three factors: ";
    … …  
    14141449       delete[](factor);
    14151450       break;
    1416        
     1451
    14171452      case 'd': // duplicate the periodic cell along a given axis, given times
    14181453        cout << Verbose(0) << "State the axis [(+-)123]: ";
    … …  
    14201455        cout << Verbose(0) << "State the factor: ";
    14211456        cin >> faktor;
    1422        
     1457
    14231458        mol->CountAtoms((ofstream *)&cout);  // recount atoms
    14241459        if (mol->AtomCount != 0) {  // if there is more than none
    … …  
    14611496            mol->Translate(&x);
    14621497          }
    1463           mol->cell_size[(abs(axis) == 2) ? 2 : ((abs(axis) == 3) ? 5 : 0)] *= faktor; 
     1498          mol->cell_size[(abs(axis) == 2) ? 2 : ((abs(axis) == 3) ? 5 : 0)] *= faktor;
    14641499        }
    14651500        break;
    1466      
     1501
    14671502      case 'e': // edit each field of the configuration
    14681503       configuration.Edit(mol);
    14691504       break;
    1470  
     1505
    14711506      case 'f':
    14721507        FragmentAtoms(mol, &configuration);
    14731508        break;
    1474        
     1509
    14751510      case 'g': // center the atoms
    14761511        CenterAtoms(mol);
    14771512        break;
    1478        
    1479       case 'i': // align all atoms 
     1513
     1514      case 'i': // align all atoms
    14801515        AlignAtoms(periode, mol);
    14811516        break;
    … …  
    14881523        MirrorAtoms(mol);
    14891524        break;
    1490        
     1525
    14911526      case 'o': // create the connection matrix
    14921527        {
    … …  
    15091544        }
    15101545        break;
    1511        
     1546
    15121547      case 'p': // parse and XYZ file
    15131548        cout << Verbose(0) << "Format should be XYZ with: ShorthandOfElement\tX\tY\tZ" << endl;
    … …  
    15181553        break;
    15191554
    1520       case 'q': // quit 
    1521         break;
    1522        
     1555      case 'q': // quit
     1556        break;
     1557
    15231558      case 'r': // remove atom
    1524         RemoveAtoms(mol);       
    1525         break;
    1526        
     1559        RemoveAtoms(mol);
     1560        break;
     1561
    15271562      case 's': // save to config file
    15281563        SaveConfig(ConfigFileName, &configuration, periode, mol);
    … …  
    15301565
    15311566      case 't': // translate all atoms
    1532        cout << Verbose(0) << "Enter translation vector." << endl;       
     1567       cout << Verbose(0) << "Enter translation vector." << endl;
    15331568       x.AskPosition(mol->cell_size,0);
    15341569       mol->Translate((const Vector *)&x);
    15351570       break;
    1536  
     1571
    15371572      case 'T':
    15381573        testroutine(mol);
    15391574        break;
    1540      
     1575
    15411576      case 'u': // change an atom's element
    15421577        first = NULL;
    … …  
    15451580          cin >> Z;
    15461581        } while ((first = mol->FindAtom(Z)) == NULL);
    1547         cout << Verbose(0) << "New element by atomic number Z: ";       
     1582        cout << Verbose(0) << "New element by atomic number Z: ";
    15481583        cin >> Z;
    15491584        first->type = periode->FindElement(Z);
    1550         cout << Verbose(0) << "Atom " << first->nr << "'s element is " << first->type->name << "." << endl;   
     1585        cout << Verbose(0) << "Atom " << first->nr << "'s element is " << first->type->name << "." << endl;
    15511586        break;
    15521587    };
    15531588  } while (choice != 'q');
    1554  
     1589
    15551590  // save element data base
    15561591  if (periode->StorePeriodentafel(ElementsFileName)) //ElementsFileName
Note: See TracChangeset for help on using the changeset viewer.