source: molecuilder/src/atom.cpp@ 0f508c

Last change on this file since 0f508c was 644ba1, checked in by Frederik Heber <heber@…>, 17 years ago

Adaptivity fixes, MD by VerletForceIntegration introduced, MD molecule::Trajectories, atom Max::Order, no more recursive going down the fragmentation level

MD
==
molecule::Trajectories is now a map to a struct trajectory list of all the MD steps.
struct Trajectory contains STL vectors of coordinates, velocities and forces. Both are needed for the new VerletForceIntegration.
Parsing of coordinates, velocities and forces from the config file was completely rewritten:

  • in the FastParsing case, we just scan for IonType1_1 to find the last step, set the file pointer there and scan all remaining ones
  • in the other case, we create the atoms in the first step and put them in a hash for lookup on later steps and read in sequentially (with file pointer moving on).
  • This is a lot faster than the old variant.

VerletForceIntegration() implemented in a working manner with force smoothing (mean of actual and last one).
OutputTrajectoriesXYZ() now concatenates the single MD steps into one xyz file, so that the animation can be viewed with e.g. jmol or vmd
config:Deltat is now public (lazy me) and set to 1 instead of 0 initially, also Deltat is parsed accordingly (if not present, defaults to 1)
MatrixContainer::ParseMatrix from parser.cpp is used during VerletForceIntegration() to parse the forces file. Consequently, we have included parser.* in the Makefile.am.
Fix: MoleculeListClass::OutputConfigForListOfFragments() stores config file under config::configpath, before it backup'd the path twice to PathBackup.

Adaptivity
==========
Adaptivity (CheckOrderAtSite()) had compared not absolute, but real values which caused sign problems and faulty behaviour.
Adapatvity (CheckOrderAtSite()) had included atoms by Order (desired one) not by FragOrder (current one) into the list of candidates, which caused faulty behaviour.
CheckOrderAtSite() did single stepping wrong as the mask bit (last in AtomMask) was checked for true not false! Also bit was not set to false initially in FragmentMolecule().
Adaptivity: FragmentMolecule now returns 1 - continue, 2 - don't ... to tell whether we still have to continue with the adaptive cycle (is given as return value from molecuilder)
introduced atom::MaxOrder
StoreOrderAtSiteFile() now also stores the MaxOrder and ParseOrderAtSiteFromFile() parses it back into atom class

Removed Fragmentation Recursion
===============================
As we switched in the prelude of the adaptivity to a monotonous increase from order 1 to the desired one, we don't need to recursively go down each level from a given current bond order, as all these fragments have already been created on the lower orders. Consequently, FragmentBOSSANOVA does not contain NumLevels or alike anymore. This simplified the code a bit, but probably is not yet completely done. (but is working, checked on heptan).
SPFragmentGenerator() does not need Labels anymore, global ones are used for checks. Consequently, PowerSetGenerator() and FragmentSearch structure don't initialise/contain them anymore. We always compare against ...->GetTrueFather()->nr.
  • Property mode set to 100644
File size: 2.7 KB
Line 
1/** \file atom.cpp
2 *
3 * Function implementations for the class atom.
4 *
5 */
6
7#include "molecules.hpp"
8
9/************************************* Functions for class atom *************************************/
10
11/** Constructor of class atom.
12 */
13atom::atom()
14{
15 Name = NULL;
16 previous = NULL;
17 next = NULL;
18 father = this; // generally, father is itself
19 Ancestor = NULL;
20 type = NULL;
21 sort = NULL;
22 FixedIon = 0;
23 nr = -1;
24 GraphNr = -1;
25 ComponentNr = NULL;
26 IsCyclic = false;
27 SeparationVertex = false;
28 LowpointNr = -1;
29 AdaptiveOrder = 0;
30 MaxOrder = false;
31};
32
33/** Destructor of class atom.
34 */
35atom::~atom()
36{
37 Free((void **)&Name, "atom::~atom: *Name");
38 Free((void **)&ComponentNr, "atom::~atom: *ComponentNr");
39};
40
41
42/** Climbs up the father list until NULL, last is returned.
43 * \return true father, i.e. whose father points to itself, NULL if it could not be found or has none (added hydrogen)
44 */
45atom *atom::GetTrueFather()
46{
47 atom *walker = this;
48 do {
49 if (walker == walker->father) // top most father is the one that points on itself
50 break;
51 walker = walker->father;
52 } while (walker != NULL);
53 return walker;
54};
55
56/** Output of a single atom.
57 * \param ElementNo cardinal number of the element
58 * \param AtomNo cardinal number among these atoms of the same element
59 * \param *out stream to output to
60 */
61bool atom::Output(int ElementNo, int AtomNo, ofstream *out) const
62{
63 if (out != NULL) {
64 *out << "Ion_Type" << ElementNo << "_" << AtomNo << "\t" << fixed << setprecision(9) << showpoint;
65 *out << x.x[0] << "\t" << x.x[1] << "\t" << x.x[2];
66 *out << "\t" << FixedIon;
67 if (v.Norm() > MYEPSILON)
68 *out << "\t" << scientific << setprecision(6) << v.x[0] << "\t" << v.x[1] << "\t" << v.x[2] << "\t";
69 *out << " # Number in molecule " << nr << endl;
70 return true;
71 } else
72 return false;
73};
74
75/** Output of a single atom as one lin in xyz file.
76 * \param *out stream to output to
77 */
78bool atom::OutputXYZLine(ofstream *out) const
79{
80 if (out != NULL) {
81 *out << type->symbol << "\t" << x.x[0] << "\t" << x.x[1] << "\t" << x.x[2] << "\t" << endl;
82 return true;
83 } else
84 return false;
85};
86
87ostream & operator << (ostream &ost, atom &a)
88{
89 ost << "[" << a.Name << "|" << &a << "]";
90 return ost;
91};
92
93/** Compares the indices of \a this atom with a given \a ptr.
94 * \param ptr atom to compare index against
95 * \return true - this one's is smaller, false - not
96 */
97bool atom::Compare(atom &ptr)
98{
99 if (nr < ptr.nr)
100 return true;
101 else
102 return false;
103};
104
105bool operator < (atom &a, atom &b)
106{
107 return a.Compare(b);
108};
109
Note: See TracBrowser for help on using the repository browser.