Context Navigation

-              rf683fe
+              ra8bcea6
 /** \file builder.cpp
+ *
+ *
  * By stating absolute positions or binding angles and distances atomic positions of a molecule can be constructed.
  * The output is the complete configuration file for PCP for direct use.
 …
  * -# Atomic data is retrieved from a file, if not found requested and stored there for later re-use
  * -# step-by-step construction of the molecule beginning either at a centre of with a certain atom
+ *
+ *
  */
 /*! \mainpage Molecuilder - a molecular set builder
+ *
+ *
  * This introductory shall briefly make aquainted with the program, helping in installing and a first run.
+ *
+ *
  * \section about About the Program
+ *
+ *
  *  Molecuilder is a short program, written in C++, that enables the construction of a coordinate set for the
  *  atoms making up an molecule by the successive statement of binding angles and distances and referencing to
  *  already constructed atoms.
+ *
+ *
  *  A configuration file may be written that is compatible to the format used by PCP - a parallel Car-Parrinello
  *  molecular dynamics implementation.
+ *
+ *
  * \section install Installation
+ *
+ *
  *  Installation should without problems succeed as follows:
  *  -# ./configure (or: mkdir build;mkdir run;cd build; ../configure --bindir=../run)
  *  -# make
  *  -# make install
+ *
+ *
  *  Further useful commands are
  *  -# make clean uninstall: deletes .o-files and removes executable from the given binary directory\n
  *  -# make doxygen-doc: Creates these html pages out of the documented source
+ *
+ *  -# make doxygen-doc: Creates these html pages out of the documented source
+ *
  * \section run Running
+ *
+ *
  *  The program can be executed by running: ./molecuilder
+ *
+ *
  *  Note, that it uses a database, called "elements.db", in the executable's directory. If the file is not found,
  *  it is created and any given data on elements of the periodic table will be stored therein and re-used on
  *  later re-execution.
+ *
+ *  later re-execution.
+ *
  * \section ref References
+ *
+ *
  *  For the special configuration file format, see the documentation of pcp.
+ *
+ *
  */
 …
   cout << Verbose(0) << "INPUT: ";
   cin >> choice;
   switch (choice) {
       case 'a': // absolute coordinates of atom
 …
         mol->AddAtom(first);  // add to molecule
         break;
       case 'b': // relative coordinates of atom wrt to reference point
         first = new atom;
 …
         mol->AddAtom(first);  // add to molecule
         break;
       case 'c': // relative coordinates of atom wrt to already placed atom
         first = new atom;
 …
         do {
           if (!valid) cout << Verbose(0) << "Resulting position out of cell." << endl;
           second = mol->AskAtom("Enter atom number: ");
+          second = mol->AskAtom("Enter atom number: ");
           cout << Verbose(0) << "Enter relative coordinates." << endl;
           first->x.AskPosition(mol->cell_size, false);
 …
         mol->AddAtom(first);  // add to molecule
         break;
       case 'd': // two atoms, two angles and a distance
         first = new atom;
 …
           x.Copyvector(&fourth->x);
           x.SubtractVector(&third->x);
           if (!z.SolveSystem(&x,&y,&n, b, c, a)) {
             cout << Verbose(0) << "Failure solving self-dependent linear system!" << endl;
 …
           cout << "x: ",
           x.Output((ofstream *)&cout);
           cout << endl;
+          cout << endl;
           z.MakeNormalVector(&second->x,&third->x,&fourth->x);
           cout << "z: ",
           z.Output((ofstream *)&cout);
           cout << endl;
+          cout << endl;
           y.MakeNormalVector(&x,&z);
           cout << "y: ",
           y.Output((ofstream *)&cout);
           cout << endl;
+          cout << endl;
           // rotate vector around first angle
           first->x.CopyVector(&x);
 …
           cout << "Rotated vector: ",
           first->x.Output((ofstream *)&cout);
           cout << endl;
+          cout << endl;
           // remove the projection onto the rotation plane of the second angle
           n.CopyVector(&y);
 …
           cout << "N1: ",
           n.Output((ofstream *)&cout);
           cout << endl;
+          cout << endl;
           first->x.SubtractVector(&n);
           cout << "Subtracted vector: ",
           first->x.Output((ofstream *)&cout);
           cout << endl;
+          cout << endl;
           n.CopyVector(&z);
           n.Scale(first->x.Projection(&z));
           cout << "N2: ",
           n.Output((ofstream *)&cout);
           cout << endl;
+          cout << endl;
           first->x.SubtractVector(&n);
           cout << "2nd subtracted vector: ",
           first->x.Output((ofstream *)&cout);
           cout << endl;
+          cout << endl;
           // rotate another vector around second angle
           n.CopyVector(&y);
 …
           cout << "2nd Rotated vector: ",
           n.Output((ofstream *)&cout);
           cout << endl;
+          cout << endl;
           // add the two linear independent vectors
           first->x.AddVector(&n);
           first->x.Normalize();
+          first->x.Normalize();
           first->x.Scale(a);
           first->x.AddVector(&second->x);
           cout << Verbose(0) << "resulting coordinates: ";
           first->x.Output((ofstream *)&cout);
 …
         } while ((j != -1) && (i<128));
         if (i >= 2) {
           first->x.LSQdistance(atoms, i);
+          first->x.LSQdistance(atoms, i);
           first->x.Output((ofstream *)&cout);
 …
   Vector x, y;
   char choice;  // menu choice char
   cout << Verbose(0) << "===========CENTER ATOMS=========================" << endl;
   cout << Verbose(0) << " a - on origin" << endl;
 …
   cout << Verbose(0) << "INPUT: ";
   cin >> choice;
   switch (choice) {
     default:
 …
   cout << Verbose(0) << "INPUT: ";
   cin >> choice;
   switch (choice) {
     default:
 …
       second = mol->AskAtom("Enter second atom: ");
       n.CopyVector((const Vector *)&first->x);
       n.SubtractVector((const Vector *)&second->x);
+      n.CopyVector((const Vector *)&first->x);
+      n.SubtractVector((const Vector *)&second->x);
       n.Normalize();
       break;
+      break;
     case 'd':
       char shorthand[4];
 …
         x.x[i] = gsl_vector_get(param.x,i);
         n.x[i] = gsl_vector_get(param.x,i+NDIM);
+      }
+      }
       gsl_vector_free(param.x);
       cout << Verbose(0) << "Offset vector: ";
 …
       cout << Verbose(0) << endl;
       n.Normalize();
       break;
+      break;
   };
   cout << Verbose(0) << "Alignment vector: ";
 …
   Vector n;
   char choice;  // menu choice char
   cout << Verbose(0) << "===========MIRROR ATOMS=========================" << endl;
   cout << Verbose(0) << " a - state three atoms defining mirror plane" << endl;
 …
   cout << Verbose(0) << "INPUT: ";
   cin >> choice;
   switch (choice) {
     default:
 …
       second = mol->AskAtom("Enter second atom: ");
       n.CopyVector((const Vector *)&first->x);
       n.SubtractVector((const Vector *)&second->x);
+      n.CopyVector((const Vector *)&first->x);
+      n.SubtractVector((const Vector *)&second->x);
       n.Normalize();
       break;
+      break;
   };
   cout << Verbose(0) << "Normal vector: ";
 …
   double tmp1, tmp2;
   char choice;  // menu choice char
   cout << Verbose(0) << "===========REMOVE ATOMS=========================" << endl;
   cout << Verbose(0) << " a - state atom for removal by number" << endl;
 …
   cout << Verbose(0) << "INPUT: ";
   cin >> choice;
   switch (choice) {
     default:
 …
           mol->RemoveAtom(first);
+      }
       break;
+      break;
   };
   //mol->Output((ofstream *)&cout);
 …
   int Z;
   char choice;  // menu choice char
   cout << Verbose(0) << "===========MEASURE ATOMS=========================" << endl;
   cout << Verbose(0) << " a - calculate bond length between one atom and all others" << endl;
 …
       for (int i=MAX_ELEMENTS;i--;)
         min[i] = 0.;
       second = mol->start;
+      second = mol->start;
       while ((second->next != mol->end)) {
         second = second->next; // advance
 …
+        }
         if ((tmp1 != 0.) && ((min[Z] == 0.) || (tmp1 < min[Z]))) min[Z] = tmp1;
         //cout << Verbose(0) << "Bond length between Atom " << first->nr << " and " << second->nr << ": " << tmp1 << " a.u." << endl;
+        //cout << Verbose(0) << "Bond length between Atom " << first->nr << " and " << second->nr << ": " << tmp1 << " a.u." << endl;
+      }
       for (int i=MAX_ELEMENTS;i--;)
         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;
       break;
     case 'b':
       first = mol->AskAtom("Enter first atom: ");
 …
       y.SubtractVector((const Vector *)&second->x);
       cout << Verbose(0) << "Bond angle between first atom Nr." << first->nr << ", central atom Nr." << second->nr << " and last atom Nr." << third->nr << ": ";
       cout << Verbose(0) << (acos(x.ScalarProduct((const Vector *)&y)/(y.Norm()*x.Norm()))/M_PI*180.) << " degrees" << endl;
+      cout << Verbose(0) << (acos(x.ScalarProduct((const Vector *)&y)/(y.Norm()*x.Norm()))/M_PI*180.) << " degrees" << endl;
       break;
     case 'd':
 …
   int Order1;
   clock_t start, end;
   cout << Verbose(0) << "Fragmenting molecule with current connection matrix ..." << endl;
   cout << Verbose(0) << "What's the desired bond order: ";
 …
     end = clock();
     cout << Verbose(0) << "Clocks for this operation: " << (end-start) << ", time: " << ((double)(end-start)/CLOCKS_PER_SEC) << "s." << endl;
   } else
+  } else
     cout << Verbose(0) << "Connection matrix has not yet been generated!" << endl;
 };
 …
   atom *Walker = mol->start;
   int i, comp, counter=0;
   // generate some KeySets
   cout << "Generating KeySets." << endl;
 …
   cout << "Testing insertion of already present item in KeySets." << endl;
   KeySetTestPair test;
   test = TestSets[mol->AtomCount-1].insert(Walker->nr);
+  test = TestSets[mol->AtomCount-1].insert(Walker->nr);
   if (test.second) {
     cout << Verbose(1) << "Insertion worked?!" << endl;
 …
   TestSets[mol->AtomCount].insert(mol->end->previous->previous->previous->nr);
   // constructing Graph structure
+  // constructing Graph structure
   cout << "Generating Subgraph class." << endl;
   Graph Subgraphs;
 …
   cout << "Testing insertion of already present item in Subgraph." << endl;
   GraphTestPair test2;
   test2 = Subgraphs.insert(GraphPair (TestSets[mol->AtomCount],pair<int, double>(counter++, 1.)));
+  test2 = Subgraphs.insert(GraphPair (TestSets[mol->AtomCount],pair<int, double>(counter++, 1.)));
   if (test2.second) {
     cout << Verbose(1) << "Insertion worked?!" << endl;
 …
     cout << Verbose(1) << "Insertion rejected: Present object is " << (*(test2.first)).second.first << "." << endl;
+  }
   // show graphs
   cout << "Showing Subgraph's contents, checking that it's sorted." << endl;
 …
       if ((*key) > comp)
         cout << (*key) << " ";
       else
+      else
         cout << (*key) << "! ";
       comp = (*key);
 …
   else
     cout << "failed." << endl;
   // and save to xyz file
   if (ConfigFileName != NULL) {
 …
     strcat(filename, ".xyz");
     output.open(filename, ios::trunc);
+  }
+  }
   cout << Verbose(0) << "Saving of XYZ file ";
   if (mol->MDSteps <= 1) {
 …
   output.close();
   output.clear();
   // and save as MPQC configuration
   if (ConfigFileName != NULL)
 …
   else
     cout << "failed." << endl;
   if (!strcmp(configuration->configpath, configuration->GetDefaultPath())) {
     cerr << "WARNING: config is found under different path then stated in config file::defaultpath!" << endl;
 …
   int argptr;
   PathToDatabases = LocalPath;
   if (argc > 1) { // config file specified as option
     // 1. : Parse options that just set variables or print help
 …
           case '?':
             cout << "MoleCuilder suite" << endl << "==================" << endl << endl;
             cout << "Usage: " << argv[0] << "[config file] [-{acefpsthH?vfrp}] [further arguments]" << endl;
+            cout << "Usage: " << argv[0] << "[config file] [-{acefpsthH?vfrp}] [further arguments]" << endl;
             cout << "or simply " << argv[0] << " without arguments for interactive session." << endl;
             cout << "\t-a Z x1 x2 x3\tAdd new atom of element Z at coordinates (x1,x2,x3)." << endl;
+            cout << "\t-A <source>\tCreate adjacenzy list from bonds parsed from 'dbond'-style file." <<endl;
             cout << "\t-b x1 x2 x3\tCenter atoms in domain with given edge lengths of (x1,x2,x3)." << endl;
             cout << "\t-c x1 x2 x3\tCenter atoms in domain with a minimum distance to boundary of (x1,x2,x3)." << endl;
 …
             cout << "\t-m <0/1>\tCalculate (0)/ Align in(1) PAS with greatest EV along z axis." << endl;
             cout << "\t-n\tFast parsing (i.e. no trajectories are looked for)." << endl;
             cout << "\t-o\tGet volume of the convex envelope (and store to tecplot file)." << endl;
+            cout << "\t-o <out>\tGet volume of the convex envelope (and store to tecplot file)." << endl;
             cout << "\t-p <file>\tParse given xyz file and create raw config file from it." << endl;
             cout << "\t-P <file>\tParse given forces file and append as an MD step to config file via Verlet." << endl;
             cout << "\t-r\t\tConvert file from an old pcp syntax." << endl;
             cout << "\t-t x1 x2 x3\tTranslate all atoms by this vector (x1,x2,x3)." << endl;
             cout << "\t-T <file> Store temperatures from the config file in <file>." << endl;
+            cout << "\t-T <file> Store temperatures from the config file in <file>." << endl;
             cout << "\t-s x1 x2 x3\tScale all atom coordinates by this vector (x1,x2,x3)." << endl;
             cout << "\t-u rho\tsuspend in water solution and output necessary cell lengths, average density rho and repetition." << endl;
+            cout << "\t-u rho\tsuspend in water solution and output necessary cell lengths, average density rho and repetition." << endl;
             cout << "\t-v/-V\t\tGives version information." << endl;
             cout << "Note: config files must not begin with '-' !" << endl;
 …
         argptr++;
     } while (argptr < argc);
     // 2. Parse the element database
     if (periode->LoadPeriodentafel(PathToDatabases)) {
 …
       return 1;
+    }
     // 3. Find config file name and parse if possible
     if (argv[1][0] != '-') {
 …
     } else
       config_present = absent;
     // 4. parse again through options, now for those depending on elements db and config presence
     argptr = 1;
 …
                     config_present = present;
                 } else
                   cerr << Verbose(1) << "Could not find the specified element." << endl;
+                  cerr << Verbose(1) << "Could not find the specified element." << endl;
                 argptr+=4;
+              }
 …
         if (config_present == present) {
           switch(argv[argptr-1][1]) {
             case 'D':
+            case 'D':
               ExitFlag = 1;
+              {
 …
+              }
               break;
+            case 'A':
+                ExitFlag = 1;
+                if ((argptr >= argc) || (argv[argptr][0] == '-')){
+                        ExitFlag =255;
+                        cerr << "Missing source file for bonds in molecule: -A <bond sourcefile>" << endl;
+                }
+                else{
+                        cout << "Parsing bonds from " << argv[argptr] << "." << endl;
+                    ifstream *input = new ifstream(argv[argptr]);
+                        mol->CreateAdjacencyList2((ofstream *)&cout, input);
+                        input->close();
+                }
+                break;
             case 'T':
               ExitFlag = 1;
 …
             case 'o':
               ExitFlag = 1;
               if ((argptr >= argc) || (argv[argptr][0] == '-')) {
+              if ((argptr >= argc) || (argv[argptr][0] == '-')){
                 ExitFlag = 255;
                 cerr << "Not enough or invalid arguments given for convex envelope: -o <tecplot output file>" << endl;
               } else {
                 cout << Verbose(0) << "Evaluating volume of the convex envelope.";
+                ofstream *output = new ofstream(argv[argptr], ios::trunc);
+                //$$$
                 cout << Verbose(1) << "Storing tecplot data in " << argv[argptr] << "." << endl;
-                ofstream *output = new ofstream(argv[argptr], ios::trunc);
                 VolumeOfConvexEnvelope((ofstream *)&cout, output, &configuration, NULL, mol);
                 output->close();
 …
                       mol->Translate(&x);
+                    }
                     mol->cell_size[(abs(axis) == 2) ? 2 : ((abs(axis) == 3) ? 5 : 0)] *= faktor;
+                    mol->cell_size[(abs(axis) == 2) ? 2 : ((abs(axis) == 3) ? 5 : 0)] *= faktor;
+                  }
+                }
 …
   if (j == 1) return 0; // just for -v and -h options
   if (j) return j;  // something went wrong
   // General stuff
   if (mol->cell_size[0] == 0.) {
 …
   // now the main construction loop
   cout << Verbose(0) << endl << "Now comes the real construction..." << endl;
   do {
+  do {
     cout << Verbose(0) << endl << endl;
     cout << Verbose(0) << "============Element list=======================" << endl;
 …
     cout << Verbose(0) << "-----------------------------------------------" << endl;
     cout << Verbose(0) << "d - duplicate molecule/periodic cell" << endl;
     cout << Verbose(0) << "i - realign molecule" << endl;
     cout << Verbose(0) << "m - mirror all molecules" << endl;
+    cout << Verbose(0) << "i - realign molecule" << endl;
+    cout << Verbose(0) << "m - mirror all molecules" << endl;
     cout << Verbose(0) << "t - translate molecule by vector" << endl;
     cout << Verbose(0) << "c - scale by unit transformation" << endl;
 …
     cout << Verbose(0) << "Input: ";
     cin >> choice;
     switch (choice) {
       default:
       case 'a': // add atom
         AddAtoms(periode, mol);
         choice = 'a';
         break;
+        choice = 'a';
+        break;
       case 'b': // scale a bond
         cout << Verbose(0) << "Scaling bond length between two atoms." << endl;
 …
+        }
         //cout << Verbose(0) << "New coordinates of Atom " << second->nr << " are: ";
         //second->Output(second->type->No, 1, (ofstream *)&cout);
         break;
       case 'c': // unit scaling of the metric
+        //second->Output(second->type->No, 1, (ofstream *)&cout);
+        break;
+      case 'c': // unit scaling of the metric
        cout << Verbose(0) << "Angstroem -> Bohrradius: 1.8897261\t\tBohrradius -> Angstroem: 0.52917721" << endl;
        cout << Verbose(0) << "Enter three factors: ";
 …
        delete[](factor);
        break;
       case 'd': // duplicate the periodic cell along a given axis, given times
         cout << Verbose(0) << "State the axis [(+-)123]: ";
 …
         cout << Verbose(0) << "State the factor: ";
         cin >> faktor;
         mol->CountAtoms((ofstream *)&cout);  // recount atoms
         if (mol->AtomCount != 0) {  // if there is more than none
 …
             mol->Translate(&x);
+          }
           mol->cell_size[(abs(axis) == 2) ? 2 : ((abs(axis) == 3) ? 5 : 0)] *= faktor;
+          mol->cell_size[(abs(axis) == 2) ? 2 : ((abs(axis) == 3) ? 5 : 0)] *= faktor;
+        }
         break;
       case 'e': // edit each field of the configuration
        configuration.Edit(mol);
        break;
       case 'f':
         FragmentAtoms(mol, &configuration);
         break;
       case 'g': // center the atoms
         CenterAtoms(mol);
         break;
       case 'i': // align all atoms
+      case 'i': // align all atoms
         AlignAtoms(periode, mol);
         break;
 …
         MirrorAtoms(mol);
         break;
       case 'o': // create the connection matrix
+        {
 …
+        }
         break;
       case 'p': // parse and XYZ file
         cout << Verbose(0) << "Format should be XYZ with: ShorthandOfElement\tX\tY\tZ" << endl;
 …
         break;
       case 'q': // quit
         break;
+      case 'q': // quit
+        break;
       case 'r': // remove atom
         RemoveAtoms(mol);
         break;
+        RemoveAtoms(mol);
+        break;
       case 's': // save to config file
         SaveConfig(ConfigFileName, &configuration, periode, mol);
 …
       case 't': // translate all atoms
        cout << Verbose(0) << "Enter translation vector." << endl;
+       cout << Verbose(0) << "Enter translation vector." << endl;
        x.AskPosition(mol->cell_size,0);
        mol->Translate((const Vector *)&x);
        break;
       case 'T':
         testroutine(mol);
         break;
       case 'u': // change an atom's element
         first = NULL;
 …
           cin >> Z;
         } while ((first = mol->FindAtom(Z)) == NULL);
         cout << Verbose(0) << "New element by atomic number Z: ";
+        cout << Verbose(0) << "New element by atomic number Z: ";
         cin >> Z;
         first->type = periode->FindElement(Z);
         cout << Verbose(0) << "Atom " << first->nr << "'s element is " << first->type->name << "." << endl;
+        cout << Verbose(0) << "Atom " << first->nr << "'s element is " << first->type->name << "." << endl;
         break;
     };
   } while (choice != 'q');
   // save element data base
   if (periode->StorePeriodentafel(ElementsFileName)) //ElementsFileName

Note: See TracChangeset for help on using the changeset viewer.

Context Navigation

Changeset a8bcea6 for src/builder.cpp

Legend:

src/builder.cpp

Download in other formats: