Changes in src/builder.cpp [375b458:042f82]

File:

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

src/builder.cpp

@@ -308 +308 @@
       }
       // center
-      mol->CenterInBox((ofstream *)&cout, &x);
+      mol->CenterInBox((ofstream *)&cout);
       // update Box of atoms by boundary
       mol->SetBoxDimension(&x);
@@ -436 +436 @@
 static void RemoveAtoms(molecule *mol)
 {
-  atom *first, *second, *third;
+  atom *first, *second;
   int axis;
   double tmp1, tmp2;
@@ -459 +459 @@
       break;
     case 'b':
-      third = mol->AskAtom("Enter number of atom as reference point: ");
+      second = mol->AskAtom("Enter number of atom as reference point: ");
       cout << Verbose(0) << "Enter radius: ";
       cin >> tmp1;
       first = mol->start;
-      second = first->next;
-      while(second != mol->end) {
-        first = second;
-        second = first->next;
-        if (first->x.DistanceSquared((const Vector *)&third->x) > tmp1*tmp1) // distance to first above radius ...
+      while(first->next != mol->end) {
+        first = first->next;
+        if (first->x.DistanceSquared((const Vector *)&second->x) > tmp1*tmp1) // distance to first above radius ...
           mol->RemoveAtom(first);
       }
@@ -474 +472 @@
       cout << Verbose(0) << "Which axis is it: ";
       cin >> axis;
-      cout << Verbose(0) << "Lower boundary: ";
+      cout << Verbose(0) << "Left inward boundary: ";
       cin >> tmp1;
-      cout << Verbose(0) << "Upper boundary: ";
+      cout << Verbose(0) << "Right inward boundary: ";
       cin >> tmp2;
       first = mol->start;
-      second = first->next;
-      while(second != mol->end) {
-        first = second;
-        second = first->next;
-        if ((first->x.x[axis] < tmp1) || (first->x.x[axis] > tmp2)) {// out of boundary ...
-          //cout << "Atom " << *first << " with " << first->x.x[axis] << " on axis " << axis << " is out of bounds [" << tmp1 << "," << tmp2 << "]." << endl;
+      while(first->next != mol->end) {
+        first = first->next;
+        if ((first->x.x[axis] > tmp2) || (first->x.x[axis] < tmp1)) // out of boundary ...
           mol->RemoveAtom(first);
-        }
       }
       break;
@@ -634 +628 @@
 static void ManipulateAtoms(periodentafel *periode, MoleculeListClass *molecules, config *configuration)
 {
-  atom *first, *second, *third, *fourth;
-  Vector **atoms;
+  atom *first, *second;
   molecule *mol = NULL;
   Vector x,y,z,n; // coordinates for absolute point in cell volume
   double *factor; // unit factor if desired
-  double a,b,c;
   double bond, min_bond;
   char choice;  // menu choice char
@@ -749 +741 @@
 static void ManipulateMolecules(periodentafel *periode, MoleculeListClass *molecules, config *configuration)
 {
-  atom *first, *second, *third, *fourth;
-  Vector **atoms;
+  atom *first = NULL;
   Vector x,y,z,n; // coordinates for absolute point in cell volume
-  double a,b,c;
   int j, axis, count, faktor;
   char choice;  // menu choice char
-  bool valid;
   molecule *mol = NULL;
   element **Elements;
@@ -906 +895 @@
 {
   char choice;  // menu choice char
-  bool valid;
   Vector Center;
   int nr, count;
   molecule *mol = NULL;
-  char *molname = NULL;
-  int length;
   char filename[MAXSTRINGSIZE];
 
@@ -1014 +1000 @@
 {
   char choice;  // menu choice char
-  bool valid;
 
   cout << Verbose(0) << "===========MERGE MOLECULES=====================" << endl;
@@ -1223 +1208 @@
  * \return exit code (0 - successful, all else - something's wrong)
  */
-static int ParseCommandLineOptions(int argc, char **argv, MoleculeListClass *&molecules, periodentafel *&periode, config& configuration, char *&ConfigFileName, char *&PathToDatabases)
+static int ParseCommandLineOptions(int argc, char **argv, MoleculeListClass *&molecules, periodentafel *&periode, config& configuration, char *&ConfigFileName)
 {
   Vector x,y,z,n;  // coordinates for absolute point in cell volume
@@ -1238 +1223 @@
   clock_t start,end;
   int argptr;
-  PathToDatabases = LocalPath;
+  strncpy(configuration.databasepath, LocalPath, MAXSTRINGSIZE-1);
 
   // simply create a new molecule, wherein the config file is loaded and the manipulation takes place
@@ -1260 +1245 @@
             cout << "\t-a Z x1 x2 x3\tAdd new atom of element Z at coordinates (x1,x2,x3)." << endl;
             cout << "\t-A <source>\tCreate adjacency 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-b xx xy xz yy yz zz\tCenter atoms in domain with given symmetric matrix of (xx,xy,xz,yy,yz,zz)." << endl;
             cout << "\t-B <basis>\tSetting basis to store to MPQC config files." << endl;
             cout << "\t-c x1 x2 x3\tCenter atoms in domain with a minimum distance to boundary of (x1,x2,x3)." << endl;
@@ -1272 +1257 @@
             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-N <radius> <file>\tGet non-convex-envelope." << endl;
+            cout << "\t-N\tGet non-convex-envelope." << 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-R\t\tRemove all atoms out of sphere around a given one." << 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 x1 x2 x3\tTranslate periodically all atoms by this vector (x1,x2,x3)." << endl;
+            cout << "\t-S <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;
@@ -1301 +1286 @@
             } else {
               cout << "Using " << argv[argptr] << " as elements database." << endl;
-              PathToDatabases = argv[argptr];
+              strncpy (configuration.databasepath, argv[argptr], MAXSTRINGSIZE-1);
               argptr+=1;
             }
@@ -1318 +1303 @@
 
     // 2. Parse the element database
-    if (periode->LoadPeriodentafel(PathToDatabases)) {
+    if (periode->LoadPeriodentafel(configuration.databasepath)) {
       cout << Verbose(0) << "Element list loaded successfully." << endl;
       //periode->Output((ofstream *)&cout);
@@ -1492 +1477 @@
               } else {
                 class Tesselation T;
-                int N = 15;
-                int number = 100;
                 string filename(argv[argptr+1]);
                 filename.append(".csv");
                 cout << Verbose(0) << "Evaluating non-convex envelope.";
                 cout << Verbose(1) << "Using rolling ball of radius " << atof(argv[argptr]) << " and storing tecplot data in " << argv[argptr+1] << "." << endl;
-                LinkedCell LCList(mol, atof(argv[argptr]));    // \NOTE not twice the radius??
+                LinkedCell LCList(mol, atof(argv[argptr])*2.);
                 Find_non_convex_border((ofstream *)&cout, mol, &T, &LCList, argv[argptr+1], atof(argv[argptr]));
                 //FindDistributionOfEllipsoids((ofstream *)&cout, &T, &LCList, N, number, filename.c_str());
@@ -1504 +1487 @@
               }
               break;
-            case 'T':
+            case 'S':
               ExitFlag = 1;
               if ((argptr >= argc) || (argv[argptr][0] == '-')) {
                 ExitFlag = 255;
-                cerr << "Not enough or invalid arguments given for storing tempature: -T <temperature file>" << endl;
+                cerr << "Not enough or invalid arguments given for storing tempature: -S <temperature file>" << endl;
               } else {
                 cout << Verbose(1) << "Storing temperatures in " << argv[argptr] << "." << endl;
@@ -1536 +1519 @@
               }
               break;
-            case 'R':
-              ExitFlag = 1;
-              if ((argptr+1 >= argc) || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1])))  {
-                ExitFlag = 255;
-                cerr << "Not enough or invalid arguments given for removing atoms: -R <id> <distance>" << endl;
-              } else {
-                SaveFlag = true;
-                cout << Verbose(1) << "Removing atoms around " << argv[argptr] << " with radius " << argv[argptr+1] << "." << endl;
-                double tmp1 = atof(argv[argptr+1]);
-                atom *third = mol->FindAtom(atoi(argv[argptr]));
-                atom *first = mol->start;
-                if ((third != NULL) && (first != mol->end)) {
-                  atom *second = first->next;
-                  while(second != mol->end) {
-                    first = second;
-                    second = first->next;
-                    if (first->x.DistanceSquared((const Vector *)&third->x) > tmp1*tmp1) // distance to first above radius ...
-                      mol->RemoveAtom(first);
-                  }
-                } else {
-                  cerr << "Removal failed due to missing atoms on molecule or wrong id." << endl;
-                }
-                argptr+=2;
-              }
-              break;
             case 't':
               ExitFlag = 1;
@@ -1573 +1531 @@
                   x.x[i] = atof(argv[argptr+i]);
                 mol->Translate((const Vector *)&x);
+                argptr+=3;
+              }
+              break;
+            case 'T':
+              ExitFlag = 1;
+              if ((argptr+2 >= argc) || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1])) || (!IsValidNumber(argv[argptr+2])) ) {
+                ExitFlag = 255;
+                cerr << "Not enough or invalid arguments given for periodic translation: -T <x> <y> <z>" << endl;
+              } else {
+                ExitFlag = 1;
+                SaveFlag = true;
+                cout << Verbose(1) << "Translating all ions periodically to new origin." << endl;
+                for (int i=NDIM;i--;)
+                  x.x[i] = atof(argv[argptr+i]);
+                mol->TranslatePeriodically((const Vector *)&x);
                 argptr+=3;
               }
@@ -1605 +1578 @@
             case 'b':
               ExitFlag = 1;
-              if ((argptr+2 >= argc) || (!IsValidNumber(argv[argptr])) || (!IsValidNumber(argv[argptr+1])) || (!IsValidNumber(argv[argptr+2])) ) {
+              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])) ) {
                 ExitFlag = 255;
-                cerr << "Not enough or invalid arguments given for centering in box: -b <length_x> <length_y> <length_z>" << endl;
+                cerr << "Not enough or invalid arguments given for centering in box: -b <xx> <xy> <xz> <yy> <yz> <zz>" << endl;
               } else {
                 SaveFlag = true;
-                j = -1;
                 cout << Verbose(1) << "Centering atoms in config file within given simulation box." << endl;
-                j=-1;
-                for (int i=0;i<NDIM;i++) {
-                  j += i+1;
-                  x.x[i] = atof(argv[argptr++]);
-                  mol->cell_size[j] += x.x[i]*2.;
+                for (int i=0;i<6;i++) {
+                  mol->cell_size[i] = atof(argv[argptr+i]);
                 }
                 // center
-                mol->CenterInBox((ofstream *)&cout, &x);
-                // update Box of atoms by boundary
-                mol->SetBoxDimension(&x);
+                mol->CenterInBox((ofstream *)&cout);
+                argptr+=6;
               }
               break;
@@ -1645 +1613 @@
                 }
                 mol->Translate((const Vector *)&x);
+                argptr+=3;
               }
               break;
@@ -1653 +1622 @@
               mol->CenterOrigin((ofstream *)&cout, &x);
               mol->SetBoxDimension(&x);
+              argptr+=0;
               break;
             case 'r':
@@ -1808 +1778 @@
     if (SaveFlag)
       SaveConfig(ConfigFileName, &configuration, periode, molecules);
-    if ((ExitFlag >= 1)) {
-      delete(mol);
-      delete(periode);
-      return (ExitFlag);
-    }
   } else {  // no arguments, hence scan the elements db
-    if (periode->LoadPeriodentafel(PathToDatabases))
+    if (periode->LoadPeriodentafel(configuration.databasepath))
       cout << Verbose(0) << "Element list loaded successfully." << endl;
     else
@@ -1820 +1785 @@
     configuration.RetrieveConfigPathAndName("main_pcp_linux");
   }
-  return(0);
+  return(ExitFlag);
 };
 
@@ -1831 +1796 @@
   molecule *mol = NULL;
   config configuration;
-  double tmp1;
-  atom *first, *second;
   char choice;  // menu choice char
   Vector x,y,z,n;  // coordinates for absolute point in cell volume
-  bool valid; // flag if input was valid or not
   ifstream test;
   ofstream output;
   string line;
   char *ConfigFileName = NULL;
-  char *ElementsFileName = NULL;
-  int Z;
-  int j, axis, count, faktor;
+  int j, count;
 
   // =========================== PARSE COMMAND LINE OPTIONS ====================================
-  j = ParseCommandLineOptions(argc, argv, molecules, periode, configuration, ConfigFileName, ElementsFileName);
-  if (j == 1) return 0; // just for -v and -h options
-  if (j) return j;  // something went wrong
+  j = ParseCommandLineOptions(argc, argv, molecules, periode, configuration, ConfigFileName);
+  switch(j) {
+    case 0:  // something went wrong
+      delete(molecules); // also free's all molecules contained
+      delete(periode);
+      return j;
+      break;
+    case 1:  // just for -v and -h options
+      delete(molecules); // also free's all molecules contained
+      delete(periode);
+      return 0;
+      break;
+    default:
+      break;
+  }
 
   // General stuff
@@ -1936 +1908 @@
 
   // save element data base
-  if (periode->StorePeriodentafel(ElementsFileName)) //ElementsFileName
+  if (periode->StorePeriodentafel(configuration.databasepath)) //ElementsFileName
     cout << Verbose(0) << "Saving of elements.db successful." << endl;
   else
     cout << Verbose(0) << "Saving of elements.db failed." << endl;
 
-  delete(molecules);
+  delete(molecules); // also free's all molecules contained
   delete(periode);
   return (0);