Changes in src/builder.cpp [fedd5f:14d4d4]

File:

• src/builder.cpp (modified) (12 diffs)

src/builder.cpp

@@ -52 +52 @@
 #include "helpers.hpp"
 #include "molecules.hpp"
+#include "boundary.hpp"
 
 /********************************************** Submenu routine **************************************/
@@ -484 +485 @@
  * \param *mol the molecule with all the atoms
  */
-static void MeasureAtoms(periodentafel *periode, molecule *mol)
+static void MeasureAtoms(periodentafel *periode, molecule *mol, config *configuration)
 {
   atom *first, *second, *third;
@@ -496 +497 @@
   cout << Verbose(0) << " b - calculate bond length between two atoms" << endl;
   cout << Verbose(0) << " c - calculate bond angle" << endl;
+  cout << Verbose(0) << " d - calculate principal axis of the system" << endl;
+  cout << Verbose(0) << " e - calculate volume of the convex envelope" << endl;
   cout << Verbose(0) << "all else - go back" << endl;
   cout << Verbose(0) << "===============================================" << endl;
@@ -553 +556 @@
       cout << Verbose(0) << (acos(x.ScalarProduct((const vector *)&y)/(y.Norm()*x.Norm()))/M_PI*180.) << " degrees" << endl;         
       break;
+    case 'd':
+        cout << Verbose(0) << "Evaluating prinicipal axis." << endl;
+        cout << Verbose(0) << "Shall we rotate? [0/1]: ";
+        cin >> Z;
+        if ((Z >=0) && (Z <=1))
+          mol->PrincipalAxisSystem((ofstream *)&cout, (bool)Z);
+        else
+          mol->PrincipalAxisSystem((ofstream *)&cout, false);
+        break;
+    case 'e':
+        cout << Verbose(0) << "Evaluating volume of the convex envelope.";
+        VolumeOfConvexEnvelope((ofstream *)&cout, configuration, NULL, mol);
+        break;
   }
 };
@@ -720 +736 @@
   int ExitFlag = 0;
   int j;
+  double volume = 0.;
   enum ConfigStatus config_present = absent;
   clock_t start,end;
@@ -742 +759 @@
             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-d x1 x2 x3\tDuplicate cell along each axis by given factor." << endl;
             cout << "\t-e <file>\tSets the databases path to be parsed (default: ./)." << endl;
             cout << "\t-f <dist> <order>\tFragments the molecule in BOSSANOVA manner and stores config files in same dir as config." << endl;
             cout << "\t-h/-H/-?\tGive this help screen." << endl;
+            cout << "\t-m\tAlign in 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-p <file>\tParse given xyz file and create raw config file from it." << 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-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-v/-V\t\tGives version information." << endl;
             cout << "Note: config files must not begin with '-' !" << endl;
@@ -768 +790 @@
             argptr+=1;
             break;
+          case 'n':
+            cout << "I won't parse trajectories." << endl;
+            configuration.FastParsing = true;
           default:   // no match? Step on
             argptr++;
@@ -780 +805 @@
       cout << Verbose(0) << "Element list loaded successfully." << endl;
       periode->Output((ofstream *)&cout);
-    } else
+    } else {
       cout << Verbose(0) << "Element list loading failed." << endl;
+      return 1;
+    }
     
     // 3. Find config file name and parse if possible
@@ -945 +972 @@
               }
               break;
+            case 'm':
+              ExitFlag = 1;
+              j = atoi(argv[argptr++]);
+              if ((j<0) || (j>1)) {
+                cerr << Verbose(1) << "ERROR: Argument of '-m' should be either 0 for no-rotate or 1 for rotate." << endl;
+                j = 0;
+              }
+              if (j)
+                cout << Verbose(0) << "Converting to prinicipal axis system." << endl;
+              else
+                cout << Verbose(0) << "Evaluating prinicipal axis." << endl;
+              mol->PrincipalAxisSystem((ofstream *)&cout, (bool)j);
+              break;
+            case 'o':
+              ExitFlag = 1;
+              cout << Verbose(0) << "Evaluating volume of the convex envelope.";
+              VolumeOfConvexEnvelope((ofstream *)&cout, &configuration, NULL, mol);
+              break;
+            case 'U':
+              volume = atof(argv[argptr++]);
+              cout << Verbose(0) << "Using " << volume << " angstrom^3 as the volume instead of convex envelope one's." << endl;  
+            case 'u':
+              {
+                double density;
+                ExitFlag = 1;
+                cout << Verbose(0) << "Evaluating necessary cell volume for a cluster suspended in water.";
+                density = atof(argv[argptr++]);
+                if (density < 1.0) {
+                  cerr << Verbose(0) << "Density must be greater than 1.0g/cm^3 !" << endl;
+                  density = 1.3;
+                }
+//                for(int i=0;i<NDIM;i++) {
+//                  repetition[i] = atoi(argv[argptr++]);
+//                  if (repetition[i] < 1)
+//                    cerr << Verbose(0) << "ERROR: repetition value must be greater 1!" << endl;
+//                  repetition[i] = 1;
+//                }
+                PrepareClustersinWater((ofstream *)&cout, &configuration, mol, volume, density);
+              }
+              break;
+            case 'd':
+              for (int axis = 1; axis <= NDIM; axis++) {
+                int faktor = atoi(argv[argptr++]);
+                int count;
+                element ** Elements;
+                vector ** Vectors;
+                if (faktor < 1) {
+                  cerr << Verbose(0) << "ERROR: Repetition faktor mus be greater than 1!" << endl;
+                  faktor = 1;
+                }
+                mol->CountAtoms((ofstream *)&cout);  // recount atoms
+                if (mol->AtomCount != 0) {  // if there is more than none
+                  count = mol->AtomCount;   // is changed becausing of adding, thus has to be stored away beforehand
+                  Elements = new element *[count];
+                  Vectors = new vector *[count];
+                  j = 0;
+                  first = mol->start;
+                  while (first->next != mol->end) {  // make a list of all atoms with coordinates and element
+                    first = first->next;
+                    Elements[j] = first->type;
+                    Vectors[j] = &first->x;
+                    j++;
+                  }
+                  if (count != j)
+                    cout << Verbose(0) << "ERROR: AtomCount " << count << " is not equal to number of atoms in molecule " << j << "!" << endl;
+                  x.Zero();
+                  y.Zero();
+                  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
+                  for (int i=1;i<faktor;i++) {  // then add this list with respective translation factor times
+                    x.AddVector(&y); // per factor one cell width further
+                    for (int k=count;k--;) { // go through every atom of the original cell
+                      first = new atom(); // create a new body
+                      first->x.CopyVector(Vectors[k]);  // use coordinate of original atom
+                      first->x.AddVector(&x);      // translate the coordinates
+                      first->type = Elements[k];  // insert original element
+                      mol->AddAtom(first);        // and add to the molecule (which increments ElementsInMolecule, AtomCount, ...)
+                    }
+                  }
+                  // free memory
+                  delete[](Elements);
+                  delete[](Vectors);
+                  // correct cell size
+                  if (axis < 0) { // if sign was negative, we have to translate everything
+                    x.Zero();
+                    x.AddVector(&y);
+                    x.Scale(-(faktor-1));
+                    mol->Translate(&x);
+                  }
+                  mol->cell_size[(abs(axis) == 2) ? 2 : ((abs(axis) == 3) ? 5 : 0)] *= faktor; 
+                }
+              }
+              break;
             default:   // no match? Step on
-              argptr++;
+              if (argv[argptr][0] != '-') // if it started with a '-' we've already made a step!
+                argptr++;
               break;
           }
@@ -1151 +1271 @@
 
       case 'l': // measure distances or angles
-        MeasureAtoms(periode, mol);
+        MeasureAtoms(periode, mol, &configuration);
         break;
 
@@ -1164 +1284 @@
         mol->CreateAdjacencyList((ofstream *)&cout, tmp1, configuration.GetIsAngstroem());
         //mol->CreateListOfBondsPerAtom((ofstream *)&cout);
-        Subgraphs = mol->DepthFirstSearchAnalysis((ofstream *)&cout, false, MinimumRingSize);
+        Subgraphs = mol->DepthFirstSearchAnalysis((ofstream *)&cout, MinimumRingSize);
         while (Subgraphs->next != NULL) {
           Subgraphs = Subgraphs->next;
@@ -1220 +1340 @@
   
   // save element data base
-  if (periode->StorePeriodentafel()) //ElementsFileName
+  if (periode->StorePeriodentafel(ElementsFileName)) //ElementsFileName
     cout << Verbose(0) << "Saving of elements.db successful." << endl;
   else