Changeset 3f0c46 for molecuilder/src/builder.cpp

Timestamp:

Jun 13, 2008, 9:18:57 AM (17 years ago)

Author:

Frederik Heber <heber@…>

Children:

0779a9

Parents:

41c2f1

Message:

Final touches for Create Clusters in water

repetition in CreateClustersinWater(), some stuff thrown out in VolumeOfConvexEnvelope (loop to throw out some more points, which does not happen if we take any appearing in at least one BoundaryPoints), some variable declarations shifted to the start of a function, BoundaryPoints are given to functions and if NULL is created by calling GetBoundaryPoints() and free'd subsequently.

File:

• molecuilder/src/builder.cpp (modified) (3 diffs)

molecuilder/src/builder.cpp

@@ -758 +758 @@
             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;
@@ -768 +769 @@
             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;
@@ -983 +985 @@
               ExitFlag = 1;
               cout << Verbose(0) << "Evaluating volume of the convex envelope.";
-//              tmp = atof(argv[argptr++]);
-//              if (tmp < 1.0) {
-//                cerr << Verbose(0) << "Density must be greater than 1.0g/cm^3 !" << endl;
-//                tmp = 1.3;
-//              }
               VolumeOfConvexEnvelope((ofstream *)&cout, &configuration, NULL, mol);
+              break;
+            case 'u':
+              {
+                double tmp;
+                ExitFlag = 1;
+                cout << Verbose(0) << "Evaluating necessary cell volume for a cluster suspended in water.";
+                tmp = atof(argv[argptr++]);
+                if (tmp < 1.0) {
+                  cerr << Verbose(0) << "Density must be greater than 1.0g/cm^3 !" << endl;
+                  tmp = 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, tmp);
+              }
+              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