| 1 | /** \file config.cpp
 | 
|---|
| 2 |  *
 | 
|---|
| 3 |  * Function implementations for the class config.
 | 
|---|
| 4 |  *
 | 
|---|
| 5 |  */
 | 
|---|
| 6 | 
 | 
|---|
| 7 | #include "molecules.hpp"
 | 
|---|
| 8 | 
 | 
|---|
| 9 | /************************************* Functions for class config ***************************/
 | 
|---|
| 10 | 
 | 
|---|
| 11 | /** Constructor for config file class.
 | 
|---|
| 12 |  */
 | 
|---|
| 13 | config::config()
 | 
|---|
| 14 | {
 | 
|---|
| 15 |         mainname = (char *) MallocString(sizeof(char)*MAXSTRINGSIZE,"config constructor: mainname");
 | 
|---|
| 16 |         defaultpath = (char *) MallocString(sizeof(char)*MAXSTRINGSIZE,"config constructor: mainname");
 | 
|---|
| 17 |         pseudopotpath = (char *) MallocString(sizeof(char)*MAXSTRINGSIZE,"config constructor: mainname");
 | 
|---|
| 18 |         configpath = (char *) MallocString(sizeof(char)*MAXSTRINGSIZE,"config constructor: mainname");
 | 
|---|
| 19 |         configname = (char *) MallocString(sizeof(char)*MAXSTRINGSIZE,"config constructor: mainname");
 | 
|---|
| 20 |         strcpy(mainname,"pcp");
 | 
|---|
| 21 |         strcpy(defaultpath,"not specified");
 | 
|---|
| 22 |         strcpy(pseudopotpath,"not specified");
 | 
|---|
| 23 |         configpath[0]='\0';
 | 
|---|
| 24 |         configname[0]='\0';
 | 
|---|
| 25 |         basis = "3-21G";
 | 
|---|
| 26 | 
 | 
|---|
| 27 |         FastParsing = false;
 | 
|---|
| 28 |         ProcPEGamma=8;
 | 
|---|
| 29 |         ProcPEPsi=1;
 | 
|---|
| 30 |         DoOutVis=0;
 | 
|---|
| 31 |         DoOutMes=1;
 | 
|---|
| 32 |         DoOutNICS=0;
 | 
|---|
| 33 |         DoOutOrbitals=0;
 | 
|---|
| 34 |         DoOutCurrent=0;
 | 
|---|
| 35 |         DoPerturbation=0;
 | 
|---|
| 36 |         DoFullCurrent=0;
 | 
|---|
| 37 |         DoWannier=0;
 | 
|---|
| 38 |         CommonWannier=0;
 | 
|---|
| 39 |         SawtoothStart=0.01;
 | 
|---|
| 40 |         VectorPlane=0;
 | 
|---|
| 41 |         VectorCut=0;
 | 
|---|
| 42 |         UseAddGramSch=1;
 | 
|---|
| 43 |         Seed=1;
 | 
|---|
| 44 | 
 | 
|---|
| 45 |         MaxOuterStep=0;
 | 
|---|
| 46 |         Deltat=1;
 | 
|---|
| 47 |         OutVisStep=10;
 | 
|---|
| 48 |         OutSrcStep=5;
 | 
|---|
| 49 |         TargetTemp=0.00095004455;
 | 
|---|
| 50 |         ScaleTempStep=25;
 | 
|---|
| 51 |         MaxPsiStep=0;
 | 
|---|
| 52 |         EpsWannier=1e-7;
 | 
|---|
| 53 | 
 | 
|---|
| 54 |         MaxMinStep=100;
 | 
|---|
| 55 |         RelEpsTotalEnergy=1e-7;
 | 
|---|
| 56 |         RelEpsKineticEnergy=1e-5;
 | 
|---|
| 57 |         MaxMinStopStep=1;
 | 
|---|
| 58 |         MaxMinGapStopStep=0;
 | 
|---|
| 59 |         MaxInitMinStep=100;
 | 
|---|
| 60 |         InitRelEpsTotalEnergy=1e-5;
 | 
|---|
| 61 |         InitRelEpsKineticEnergy=1e-4;
 | 
|---|
| 62 |         InitMaxMinStopStep=1;
 | 
|---|
| 63 |         InitMaxMinGapStopStep=0;
 | 
|---|
| 64 | 
 | 
|---|
| 65 |         //BoxLength[NDIM*NDIM];
 | 
|---|
| 66 | 
 | 
|---|
| 67 |         ECut=128.;
 | 
|---|
| 68 |         MaxLevel=5;
 | 
|---|
| 69 |         RiemannTensor=0;
 | 
|---|
| 70 |         LevRFactor=0;
 | 
|---|
| 71 |         RiemannLevel=0;
 | 
|---|
| 72 |         Lev0Factor=2;
 | 
|---|
| 73 |         RTActualUse=0;
 | 
|---|
| 74 |         PsiType=0;
 | 
|---|
| 75 |         MaxPsiDouble=0;
 | 
|---|
| 76 |         PsiMaxNoUp=0;
 | 
|---|
| 77 |         PsiMaxNoDown=0;
 | 
|---|
| 78 |         AddPsis=0;
 | 
|---|
| 79 | 
 | 
|---|
| 80 |         RCut=20.;
 | 
|---|
| 81 |         StructOpt=0;
 | 
|---|
| 82 |         IsAngstroem=1;
 | 
|---|
| 83 |         RelativeCoord=0;
 | 
|---|
| 84 |         MaxTypes=0;
 | 
|---|
| 85 | };
 | 
|---|
| 86 | 
 | 
|---|
| 87 | 
 | 
|---|
| 88 | /** Destructor for config file class.
 | 
|---|
| 89 |  */
 | 
|---|
| 90 | config::~config()
 | 
|---|
| 91 | {
 | 
|---|
| 92 |         Free((void **)&mainname, "config::~config: *mainname");
 | 
|---|
| 93 |         Free((void **)&defaultpath, "config::~config: *defaultpath");
 | 
|---|
| 94 |         Free((void **)&pseudopotpath, "config::~config: *pseudopotpath");
 | 
|---|
| 95 |         Free((void **)&configpath, "config::~config: *configpath");
 | 
|---|
| 96 |         Free((void **)&configname, "config::~config: *configname");
 | 
|---|
| 97 | };
 | 
|---|
| 98 | 
 | 
|---|
| 99 | /** Displays menu for editing each entry of the config file.
 | 
|---|
| 100 |  * Nothing fancy here, just lots of cout << Verbose(0)s for the menu and a switch/case
 | 
|---|
| 101 |  * for each entry of the config file structure.
 | 
|---|
| 102 |  */
 | 
|---|
| 103 | void config::Edit()
 | 
|---|
| 104 | {
 | 
|---|
| 105 |         char choice;
 | 
|---|
| 106 | 
 | 
|---|
| 107 |         do {
 | 
|---|
| 108 |                 cout << Verbose(0) << "===========EDIT CONFIGURATION============================" << endl;
 | 
|---|
| 109 |                 cout << Verbose(0) << " A - mainname (prefix for all runtime files)" << endl;
 | 
|---|
| 110 |                 cout << Verbose(0) << " B - Default path (for runtime files)" << endl;
 | 
|---|
| 111 |                 cout << Verbose(0) << " C - Path of pseudopotential files" << endl;
 | 
|---|
| 112 |                 cout << Verbose(0) << " D - Number of coefficient sharing processes" << endl;
 | 
|---|
| 113 |                 cout << Verbose(0) << " E - Number of wave function sharing processes" << endl;
 | 
|---|
| 114 |                 cout << Verbose(0) << " F - 0: Don't output density for OpenDX, 1: do" << endl;
 | 
|---|
| 115 |                 cout << Verbose(0) << " G - 0: Don't output physical data, 1: do" << endl;
 | 
|---|
| 116 |                 cout << Verbose(0) << " H - 0: Don't output densities of each unperturbed orbital for OpenDX, 1: do" << endl;
 | 
|---|
| 117 |                 cout << Verbose(0) << " I - 0: Don't output current density for OpenDX, 1: do" << endl;
 | 
|---|
| 118 |                 cout << Verbose(0) << " J - 0: Don't do the full current calculation, 1: do" << endl;
 | 
|---|
| 119 |                 cout << Verbose(0) << " K - 0: Don't do perturbation calculation to obtain susceptibility and shielding, 1: do" << endl;
 | 
|---|
| 120 |                 cout << Verbose(0) << " L - 0: Wannier centres as calculated, 1: common centre for all, 2: unite centres according to spread, 3: cell centre, 4: shifted to nearest grid point" << endl;
 | 
|---|
| 121 |                 cout << Verbose(0) << " M - Absolute begin of unphysical sawtooth transfer for position operator within cell" << endl;
 | 
|---|
| 122 |                 cout << Verbose(0) << " N - (0,1,2) x,y,z-plane to do two-dimensional current vector cut" << endl;
 | 
|---|
| 123 |                 cout << Verbose(0) << " O - Absolute position along vector cut axis for cut plane" << endl;
 | 
|---|
| 124 |                 cout << Verbose(0) << " P - Additional Gram-Schmidt-Orthonormalization to stabilize numerics" << endl;
 | 
|---|
| 125 |                 cout << Verbose(0) << " Q - Initial integer value of random number generator" << endl;
 | 
|---|
| 126 |                 cout << Verbose(0) << " R - for perturbation 0, for structure optimization defines upper limit of iterations" << endl;
 | 
|---|
| 127 |                 cout << Verbose(0) << " T - Output visual after ...th step" << endl;
 | 
|---|
| 128 |                 cout << Verbose(0) << " U - Output source densities of wave functions after ...th step" << endl;
 | 
|---|
| 129 |                 cout << Verbose(0) << " X - minimization iterations per wave function, if unsure leave at default value 0" << endl;
 | 
|---|
| 130 |                 cout << Verbose(0) << " Y - tolerance value for total spread in iterative Jacobi diagonalization" << endl;
 | 
|---|
| 131 |                 cout << Verbose(0) << " Z - Maximum number of minimization iterations" << endl;
 | 
|---|
| 132 |                 cout << Verbose(0) << " a - Relative change in total energy to stop min. iteration" << endl;
 | 
|---|
| 133 |                 cout << Verbose(0) << " b - Relative change in kinetic energy to stop min. iteration" << endl;
 | 
|---|
| 134 |                 cout << Verbose(0) << " c - Check stop conditions every ..th step during min. iteration" << endl;
 | 
|---|
| 135 |                 cout << Verbose(0) << " e - Maximum number of minimization iterations during initial level" << endl;
 | 
|---|
| 136 |                 cout << Verbose(0) << " f - Relative change in total energy to stop min. iteration during initial level" << endl;
 | 
|---|
| 137 |                 cout << Verbose(0) << " g - Relative change in kinetic energy to stop min. iteration during initial level" << endl;
 | 
|---|
| 138 |                 cout << Verbose(0) << " h - Check stop conditions every ..th step during min. iteration during initial level" << endl;
 | 
|---|
| 139 | //              cout << Verbose(0) << " j - six lower diagonal entries of matrix, defining the unit cell" << endl;
 | 
|---|
| 140 |                 cout << Verbose(0) << " k - Energy cutoff of plane wave basis in Hartree" << endl;
 | 
|---|
| 141 |                 cout << Verbose(0) << " l - Maximum number of levels in multi-level-ansatz" << endl;
 | 
|---|
| 142 |                 cout << Verbose(0) << " m - Factor by which grid nodes increase between standard and upper level" << endl;
 | 
|---|
| 143 |                 cout << Verbose(0) << " n - 0: Don't use RiemannTensor, 1: Do" << endl;
 | 
|---|
| 144 |                 cout << Verbose(0) << " o - Factor by which grid nodes increase between Riemann and standard(?) level" << endl;
 | 
|---|
| 145 |                 cout << Verbose(0) << " p - Number of Riemann levels" << endl;
 | 
|---|
| 146 |                 cout << Verbose(0) << " r - 0: Don't Use RiemannTensor, 1: Do" << endl;
 | 
|---|
| 147 |                 cout << Verbose(0) << " s - 0: Doubly occupied orbitals, 1: Up-/Down-Orbitals" << endl;
 | 
|---|
| 148 |                 cout << Verbose(0) << " t - Number of orbitals (depends pn SpinType)" << endl;
 | 
|---|
| 149 |                 cout << Verbose(0) << " u - Number of SpinUp orbitals (depends on SpinType)" << endl;
 | 
|---|
| 150 |                 cout << Verbose(0) << " v - Number of SpinDown orbitals (depends on SpinType)" << endl;
 | 
|---|
| 151 |                 cout << Verbose(0) << " w - Number of additional, unoccupied orbitals" << endl;
 | 
|---|
| 152 |                 cout << Verbose(0) << " x - radial cutoff for ewald summation in Bohrradii" << endl;
 | 
|---|
| 153 |                 cout << Verbose(0) << " y - 0: Don't do structure optimization beforehand, 1: Do" << endl;
 | 
|---|
| 154 |                 cout << Verbose(0) << " z - 0: Units are in Bohr radii, 1: units are in Aengstrom" << endl;
 | 
|---|
| 155 |                 cout << Verbose(0) << " i - 0: Coordinates given in file are absolute, 1: ... are relative to unit cell" << endl;
 | 
|---|
| 156 |                 cout << Verbose(0) << "=========================================================" << endl;
 | 
|---|
| 157 |                 cout << Verbose(0) << "INPUT: ";
 | 
|---|
| 158 |                 cin >> choice;
 | 
|---|
| 159 | 
 | 
|---|
| 160 |                 switch (choice) {
 | 
|---|
| 161 |                                 case 'A': // mainname
 | 
|---|
| 162 |                                         cout << Verbose(0) << "Old: " << config::mainname << "\t new: ";
 | 
|---|
| 163 |                                         cin >> config::mainname;
 | 
|---|
| 164 |                                         break;
 | 
|---|
| 165 |                                 case 'B': // defaultpath
 | 
|---|
| 166 |                                         cout << Verbose(0) << "Old: " << config::defaultpath << "\t new: ";
 | 
|---|
| 167 |                                         cin >> config::defaultpath;
 | 
|---|
| 168 |                                         break;
 | 
|---|
| 169 |                                 case 'C': // pseudopotpath
 | 
|---|
| 170 |                                         cout << Verbose(0) << "Old: " << config::pseudopotpath << "\t new: ";
 | 
|---|
| 171 |                                         cin >> config::pseudopotpath;
 | 
|---|
| 172 |                                         break;
 | 
|---|
| 173 | 
 | 
|---|
| 174 |                                 case 'D': // ProcPEGamma
 | 
|---|
| 175 |                                         cout << Verbose(0) << "Old: " << config::ProcPEGamma << "\t new: ";
 | 
|---|
| 176 |                                         cin >> config::ProcPEGamma;
 | 
|---|
| 177 |                                         break;
 | 
|---|
| 178 |                                 case 'E': // ProcPEPsi
 | 
|---|
| 179 |                                         cout << Verbose(0) << "Old: " << config::ProcPEPsi << "\t new: ";
 | 
|---|
| 180 |                                         cin >> config::ProcPEPsi;
 | 
|---|
| 181 |                                         break;
 | 
|---|
| 182 |                                 case 'F': // DoOutVis
 | 
|---|
| 183 |                                         cout << Verbose(0) << "Old: " << config::DoOutVis << "\t new: ";
 | 
|---|
| 184 |                                         cin >> config::DoOutVis;
 | 
|---|
| 185 |                                         break;
 | 
|---|
| 186 |                                 case 'G': // DoOutMes
 | 
|---|
| 187 |                                         cout << Verbose(0) << "Old: " << config::DoOutMes << "\t new: ";
 | 
|---|
| 188 |                                         cin >> config::DoOutMes;
 | 
|---|
| 189 |                                         break;
 | 
|---|
| 190 |                                 case 'H': // DoOutOrbitals
 | 
|---|
| 191 |                                         cout << Verbose(0) << "Old: " << config::DoOutOrbitals << "\t new: ";
 | 
|---|
| 192 |                                         cin >> config::DoOutOrbitals;
 | 
|---|
| 193 |                                         break;
 | 
|---|
| 194 |                                 case 'I': // DoOutCurrent
 | 
|---|
| 195 |                                         cout << Verbose(0) << "Old: " << config::DoOutCurrent << "\t new: ";
 | 
|---|
| 196 |                                         cin >> config::DoOutCurrent;
 | 
|---|
| 197 |                                         break;
 | 
|---|
| 198 |                                 case 'J': // DoFullCurrent
 | 
|---|
| 199 |                                         cout << Verbose(0) << "Old: " << config::DoFullCurrent << "\t new: ";
 | 
|---|
| 200 |                                         cin >> config::DoFullCurrent;
 | 
|---|
| 201 |                                         break;
 | 
|---|
| 202 |                                 case 'K': // DoPerturbation
 | 
|---|
| 203 |                                         cout << Verbose(0) << "Old: " << config::DoPerturbation << "\t new: ";
 | 
|---|
| 204 |                                         cin >> config::DoPerturbation;
 | 
|---|
| 205 |                                         break;
 | 
|---|
| 206 |                                 case 'L': // CommonWannier
 | 
|---|
| 207 |                                         cout << Verbose(0) << "Old: " << config::CommonWannier << "\t new: ";
 | 
|---|
| 208 |                                         cin >> config::CommonWannier;
 | 
|---|
| 209 |                                         break;
 | 
|---|
| 210 |                                 case 'M': // SawtoothStart
 | 
|---|
| 211 |                                         cout << Verbose(0) << "Old: " << config::SawtoothStart << "\t new: ";
 | 
|---|
| 212 |                                         cin >> config::SawtoothStart;
 | 
|---|
| 213 |                                         break;
 | 
|---|
| 214 |                                 case 'N': // VectorPlane
 | 
|---|
| 215 |                                         cout << Verbose(0) << "Old: " << config::VectorPlane << "\t new: ";
 | 
|---|
| 216 |                                         cin >> config::VectorPlane;
 | 
|---|
| 217 |                                         break;
 | 
|---|
| 218 |                                 case 'O': // VectorCut
 | 
|---|
| 219 |                                         cout << Verbose(0) << "Old: " << config::VectorCut << "\t new: ";
 | 
|---|
| 220 |                                         cin >> config::VectorCut;
 | 
|---|
| 221 |                                         break;
 | 
|---|
| 222 |                                 case 'P': // UseAddGramSch
 | 
|---|
| 223 |                                         cout << Verbose(0) << "Old: " << config::UseAddGramSch << "\t new: ";
 | 
|---|
| 224 |                                         cin >> config::UseAddGramSch;
 | 
|---|
| 225 |                                         break;
 | 
|---|
| 226 |                                 case 'Q': // Seed
 | 
|---|
| 227 |                                         cout << Verbose(0) << "Old: " << config::Seed << "\t new: ";
 | 
|---|
| 228 |                                         cin >> config::Seed;
 | 
|---|
| 229 |                                         break;
 | 
|---|
| 230 | 
 | 
|---|
| 231 |                                 case 'R': // MaxOuterStep
 | 
|---|
| 232 |                                         cout << Verbose(0) << "Old: " << config::MaxOuterStep << "\t new: ";
 | 
|---|
| 233 |                                         cin >> config::MaxOuterStep;
 | 
|---|
| 234 |                                         break;
 | 
|---|
| 235 |                                 case 'T': // OutVisStep
 | 
|---|
| 236 |                                         cout << Verbose(0) << "Old: " << config::OutVisStep << "\t new: ";
 | 
|---|
| 237 |                                         cin >> config::OutVisStep;
 | 
|---|
| 238 |                                         break;
 | 
|---|
| 239 |                                 case 'U': // OutSrcStep
 | 
|---|
| 240 |                                         cout << Verbose(0) << "Old: " << config::OutSrcStep << "\t new: ";
 | 
|---|
| 241 |                                         cin >> config::OutSrcStep;
 | 
|---|
| 242 |                                         break;
 | 
|---|
| 243 |                                 case 'X': // MaxPsiStep
 | 
|---|
| 244 |                                         cout << Verbose(0) << "Old: " << config::MaxPsiStep << "\t new: ";
 | 
|---|
| 245 |                                         cin >> config::MaxPsiStep;
 | 
|---|
| 246 |                                         break;
 | 
|---|
| 247 |                                 case 'Y': // EpsWannier
 | 
|---|
| 248 |                                         cout << Verbose(0) << "Old: " << config::EpsWannier << "\t new: ";
 | 
|---|
| 249 |                                         cin >> config::EpsWannier;
 | 
|---|
| 250 |                                         break;
 | 
|---|
| 251 | 
 | 
|---|
| 252 |                                 case 'Z': // MaxMinStep
 | 
|---|
| 253 |                                         cout << Verbose(0) << "Old: " << config::MaxMinStep << "\t new: ";
 | 
|---|
| 254 |                                         cin >> config::MaxMinStep;
 | 
|---|
| 255 |                                         break;
 | 
|---|
| 256 |                                 case 'a': // RelEpsTotalEnergy
 | 
|---|
| 257 |                                         cout << Verbose(0) << "Old: " << config::RelEpsTotalEnergy << "\t new: ";
 | 
|---|
| 258 |                                         cin >> config::RelEpsTotalEnergy;
 | 
|---|
| 259 |                                         break;
 | 
|---|
| 260 |                                 case 'b': // RelEpsKineticEnergy
 | 
|---|
| 261 |                                         cout << Verbose(0) << "Old: " << config::RelEpsKineticEnergy << "\t new: ";
 | 
|---|
| 262 |                                         cin >> config::RelEpsKineticEnergy;
 | 
|---|
| 263 |                                         break;
 | 
|---|
| 264 |                                 case 'c': // MaxMinStopStep
 | 
|---|
| 265 |                                         cout << Verbose(0) << "Old: " << config::MaxMinStopStep << "\t new: ";
 | 
|---|
| 266 |                                         cin >> config::MaxMinStopStep;
 | 
|---|
| 267 |                                         break;
 | 
|---|
| 268 |                                 case 'e': // MaxInitMinStep
 | 
|---|
| 269 |                                         cout << Verbose(0) << "Old: " << config::MaxInitMinStep << "\t new: ";
 | 
|---|
| 270 |                                         cin >> config::MaxInitMinStep;
 | 
|---|
| 271 |                                         break;
 | 
|---|
| 272 |                                 case 'f': // InitRelEpsTotalEnergy
 | 
|---|
| 273 |                                         cout << Verbose(0) << "Old: " << config::InitRelEpsTotalEnergy << "\t new: ";
 | 
|---|
| 274 |                                         cin >> config::InitRelEpsTotalEnergy;
 | 
|---|
| 275 |                                         break;
 | 
|---|
| 276 |                                 case 'g': // InitRelEpsKineticEnergy
 | 
|---|
| 277 |                                         cout << Verbose(0) << "Old: " << config::InitRelEpsKineticEnergy << "\t new: ";
 | 
|---|
| 278 |                                         cin >> config::InitRelEpsKineticEnergy;
 | 
|---|
| 279 |                                         break;
 | 
|---|
| 280 |                                 case 'h': // InitMaxMinStopStep
 | 
|---|
| 281 |                                         cout << Verbose(0) << "Old: " << config::InitMaxMinStopStep << "\t new: ";
 | 
|---|
| 282 |                                         cin >> config::InitMaxMinStopStep;
 | 
|---|
| 283 |                                         break;
 | 
|---|
| 284 | 
 | 
|---|
| 285 | //                              case 'j': // BoxLength
 | 
|---|
| 286 | //                                      cout << Verbose(0) << "enter lower triadiagonalo form of basis matrix" << endl << endl;
 | 
|---|
| 287 | //                                      for (int i=0;i<6;i++) {
 | 
|---|
| 288 | //                                              cout << Verbose(0) << "Cell size" << i << ": ";
 | 
|---|
| 289 | //                                              cin >> mol->cell_size[i];
 | 
|---|
| 290 | //                                      }
 | 
|---|
| 291 | //                                      break;
 | 
|---|
| 292 | 
 | 
|---|
| 293 |                                 case 'k': // ECut
 | 
|---|
| 294 |                                         cout << Verbose(0) << "Old: " << config::ECut << "\t new: ";
 | 
|---|
| 295 |                                         cin >> config::ECut;
 | 
|---|
| 296 |                                         break;
 | 
|---|
| 297 |                                 case 'l': // MaxLevel
 | 
|---|
| 298 |                                         cout << Verbose(0) << "Old: " << config::MaxLevel << "\t new: ";
 | 
|---|
| 299 |                                         cin >> config::MaxLevel;
 | 
|---|
| 300 |                                         break;
 | 
|---|
| 301 |                                 case 'm': // RiemannTensor
 | 
|---|
| 302 |                                         cout << Verbose(0) << "Old: " << config::RiemannTensor << "\t new: ";
 | 
|---|
| 303 |                                         cin >> config::RiemannTensor;
 | 
|---|
| 304 |                                         break;
 | 
|---|
| 305 |                                 case 'n': // LevRFactor
 | 
|---|
| 306 |                                         cout << Verbose(0) << "Old: " << config::LevRFactor << "\t new: ";
 | 
|---|
| 307 |                                         cin >> config::LevRFactor;
 | 
|---|
| 308 |                                         break;
 | 
|---|
| 309 |                                 case 'o': // RiemannLevel
 | 
|---|
| 310 |                                         cout << Verbose(0) << "Old: " << config::RiemannLevel << "\t new: ";
 | 
|---|
| 311 |                                         cin >> config::RiemannLevel;
 | 
|---|
| 312 |                                         break;
 | 
|---|
| 313 |                                 case 'p': // Lev0Factor
 | 
|---|
| 314 |                                         cout << Verbose(0) << "Old: " << config::Lev0Factor << "\t new: ";
 | 
|---|
| 315 |                                         cin >> config::Lev0Factor;
 | 
|---|
| 316 |                                         break;
 | 
|---|
| 317 |                                 case 'r': // RTActualUse
 | 
|---|
| 318 |                                         cout << Verbose(0) << "Old: " << config::RTActualUse << "\t new: ";
 | 
|---|
| 319 |                                         cin >> config::RTActualUse;
 | 
|---|
| 320 |                                         break;
 | 
|---|
| 321 |                                 case 's': // PsiType
 | 
|---|
| 322 |                                         cout << Verbose(0) << "Old: " << config::PsiType << "\t new: ";
 | 
|---|
| 323 |                                         cin >> config::PsiType;
 | 
|---|
| 324 |                                         break;
 | 
|---|
| 325 |                                 case 't': // MaxPsiDouble
 | 
|---|
| 326 |                                         cout << Verbose(0) << "Old: " << config::MaxPsiDouble << "\t new: ";
 | 
|---|
| 327 |                                         cin >> config::MaxPsiDouble;
 | 
|---|
| 328 |                                         break;
 | 
|---|
| 329 |                                 case 'u': // PsiMaxNoUp
 | 
|---|
| 330 |                                         cout << Verbose(0) << "Old: " << config::PsiMaxNoUp << "\t new: ";
 | 
|---|
| 331 |                                         cin >> config::PsiMaxNoUp;
 | 
|---|
| 332 |                                         break;
 | 
|---|
| 333 |                                 case 'v': // PsiMaxNoDown
 | 
|---|
| 334 |                                         cout << Verbose(0) << "Old: " << config::PsiMaxNoDown << "\t new: ";
 | 
|---|
| 335 |                                         cin >> config::PsiMaxNoDown;
 | 
|---|
| 336 |                                         break;
 | 
|---|
| 337 |                                 case 'w': // AddPsis
 | 
|---|
| 338 |                                         cout << Verbose(0) << "Old: " << config::AddPsis << "\t new: ";
 | 
|---|
| 339 |                                         cin >> config::AddPsis;
 | 
|---|
| 340 |                                         break;
 | 
|---|
| 341 | 
 | 
|---|
| 342 |                                 case 'x': // RCut
 | 
|---|
| 343 |                                         cout << Verbose(0) << "Old: " << config::RCut << "\t new: ";
 | 
|---|
| 344 |                                         cin >> config::RCut;
 | 
|---|
| 345 |                                         break;
 | 
|---|
| 346 |                                 case 'y': // StructOpt
 | 
|---|
| 347 |                                         cout << Verbose(0) << "Old: " << config::StructOpt << "\t new: ";
 | 
|---|
| 348 |                                         cin >> config::StructOpt;
 | 
|---|
| 349 |                                         break;
 | 
|---|
| 350 |                                 case 'z': // IsAngstroem
 | 
|---|
| 351 |                                         cout << Verbose(0) << "Old: " << config::IsAngstroem << "\t new: ";
 | 
|---|
| 352 |                                         cin >> config::IsAngstroem;
 | 
|---|
| 353 |                                         break;
 | 
|---|
| 354 |                                 case 'i': // RelativeCoord
 | 
|---|
| 355 |                                         cout << Verbose(0) << "Old: " << config::RelativeCoord << "\t new: ";
 | 
|---|
| 356 |                                         cin >> config::RelativeCoord;
 | 
|---|
| 357 |                                         break;
 | 
|---|
| 358 |                 };
 | 
|---|
| 359 |         } while (choice != 'q');
 | 
|---|
| 360 | };
 | 
|---|
| 361 | 
 | 
|---|
| 362 | /** Tests whether a given configuration file adhears to old or new syntax.
 | 
|---|
| 363 |  * \param *filename filename of config file to be tested
 | 
|---|
| 364 |  * \param *periode pointer to a periodentafel class with all elements
 | 
|---|
| 365 |  * \param *mol pointer to molecule containing all atoms of the molecule
 | 
|---|
| 366 |  * \return 0 - old syntax, 1 - new syntax, -1 - unknown syntax
 | 
|---|
| 367 |  */
 | 
|---|
| 368 | int config::TestSyntax(char *filename, periodentafel *periode, molecule *mol)
 | 
|---|
| 369 | {
 | 
|---|
| 370 |         int test;
 | 
|---|
| 371 |         ifstream file(filename);
 | 
|---|
| 372 | 
 | 
|---|
| 373 |         // search file for keyword: ProcPEGamma (new syntax)
 | 
|---|
| 374 |         if (ParseForParameter(1,&file,"ProcPEGamma", 0, 1, 1, int_type, &test, 1, optional)) {
 | 
|---|
| 375 |                 file.close();
 | 
|---|
| 376 |                 return 1;
 | 
|---|
| 377 |         }
 | 
|---|
| 378 |         // search file for keyword: ProcsGammaPsi (old syntax)
 | 
|---|
| 379 |         if (ParseForParameter(1,&file,"ProcsGammaPsi", 0, 1, 1, int_type, &test, 1, optional)) {
 | 
|---|
| 380 |                 file.close();
 | 
|---|
| 381 |                 return 0;
 | 
|---|
| 382 |         }
 | 
|---|
| 383 |         file.close();
 | 
|---|
| 384 |         return -1;
 | 
|---|
| 385 | }
 | 
|---|
| 386 | 
 | 
|---|
| 387 | /** Returns private config::IsAngstroem.
 | 
|---|
| 388 |  * \return IsAngstroem
 | 
|---|
| 389 |  */
 | 
|---|
| 390 | bool config::GetIsAngstroem() const
 | 
|---|
| 391 | {
 | 
|---|
| 392 |         return (IsAngstroem == 1);
 | 
|---|
| 393 | };
 | 
|---|
| 394 | 
 | 
|---|
| 395 | /** Returns private config::*defaultpath.
 | 
|---|
| 396 |  * \return *defaultpath
 | 
|---|
| 397 |  */
 | 
|---|
| 398 | char * config::GetDefaultPath() const
 | 
|---|
| 399 | {
 | 
|---|
| 400 |         return defaultpath;
 | 
|---|
| 401 | };
 | 
|---|
| 402 | 
 | 
|---|
| 403 | 
 | 
|---|
| 404 | /** Returns private config::*defaultpath.
 | 
|---|
| 405 |  * \return *defaultpath
 | 
|---|
| 406 |  */
 | 
|---|
| 407 | void config::SetDefaultPath(const char *path)
 | 
|---|
| 408 | {
 | 
|---|
| 409 |         strcpy(defaultpath, path);
 | 
|---|
| 410 | };
 | 
|---|
| 411 | 
 | 
|---|
| 412 | /** Retrieves the path in the given config file name.
 | 
|---|
| 413 |  * \param filename config file string
 | 
|---|
| 414 |  */
 | 
|---|
| 415 | void config::RetrieveConfigPathAndName(string filename)
 | 
|---|
| 416 | {
 | 
|---|
| 417 |         char *ptr = NULL;
 | 
|---|
| 418 |         char *buffer = new char[MAXSTRINGSIZE];
 | 
|---|
| 419 |         strncpy(buffer, filename.c_str(), MAXSTRINGSIZE);
 | 
|---|
| 420 |         int last = -1;
 | 
|---|
| 421 |         for(last=MAXSTRINGSIZE;last--;) {
 | 
|---|
| 422 |                 if (buffer[last] == '/')
 | 
|---|
| 423 |                         break;
 | 
|---|
| 424 |         }
 | 
|---|
| 425 |         if (last == -1) { // no path in front, set to local directory.
 | 
|---|
| 426 |                 strcpy(configpath, "./");
 | 
|---|
| 427 |                 ptr = buffer;
 | 
|---|
| 428 |         } else {
 | 
|---|
| 429 |                 strncpy(configpath, buffer, last+1);
 | 
|---|
| 430 |                 ptr = &buffer[last+1];
 | 
|---|
| 431 |                 if (last < 254)
 | 
|---|
| 432 |                         configpath[last+1]='\0';
 | 
|---|
| 433 |         }
 | 
|---|
| 434 |         strcpy(configname, ptr);
 | 
|---|
| 435 |         cout << "Found configpath: " << configpath << ", dir slash was found at " << last << ", config name is " << configname << "." << endl;
 | 
|---|
| 436 |         delete[](buffer);
 | 
|---|
| 437 | };
 | 
|---|
| 438 | 
 | 
|---|
| 439 | 
 | 
|---|
| 440 | /** Initializes config file structure by loading elements from a give file.
 | 
|---|
| 441 |  * \param *file input file stream being the opened config file
 | 
|---|
| 442 |  * \param *periode pointer to a periodentafel class with all elements
 | 
|---|
| 443 |  * \param *mol pointer to molecule containing all atoms of the molecule
 | 
|---|
| 444 |  */
 | 
|---|
| 445 | void config::Load(char *filename, periodentafel *periode, molecule *mol)
 | 
|---|
| 446 | {
 | 
|---|
| 447 |         ifstream *file = new ifstream(filename);
 | 
|---|
| 448 |         if (file == NULL) {
 | 
|---|
| 449 |                 cerr << "ERROR: config file " << filename << " missing!" << endl;
 | 
|---|
| 450 |                 return;
 | 
|---|
| 451 |         }
 | 
|---|
| 452 |         RetrieveConfigPathAndName(filename);
 | 
|---|
| 453 |         // ParseParameters
 | 
|---|
| 454 | 
 | 
|---|
| 455 |         /* Oeffne Hauptparameterdatei */
 | 
|---|
| 456 |         int di;
 | 
|---|
| 457 |         double BoxLength[9];
 | 
|---|
| 458 |         string zeile;
 | 
|---|
| 459 |         string dummy;
 | 
|---|
| 460 |         element *elementhash[MAX_ELEMENTS];
 | 
|---|
| 461 |         char name[MAX_ELEMENTS];
 | 
|---|
| 462 |         char keyword[MAX_ELEMENTS];
 | 
|---|
| 463 |         int Z, No[MAX_ELEMENTS];
 | 
|---|
| 464 |         int verbose = 0;
 | 
|---|
| 465 |         double value[3];
 | 
|---|
| 466 | 
 | 
|---|
| 467 |         /* Namen einlesen */
 | 
|---|
| 468 | 
 | 
|---|
| 469 |         ParseForParameter(verbose,file, "mainname", 0, 1, 1, string_type, (config::mainname), 1, critical);
 | 
|---|
| 470 |         ParseForParameter(verbose,file, "defaultpath", 0, 1, 1, string_type, (config::defaultpath), 1, critical);
 | 
|---|
| 471 |         ParseForParameter(verbose,file, "pseudopotpath", 0, 1, 1, string_type, (config::pseudopotpath), 1, critical);
 | 
|---|
| 472 |         ParseForParameter(verbose,file,"ProcPEGamma", 0, 1, 1, int_type, &(config::ProcPEGamma), 1, critical);
 | 
|---|
| 473 |         ParseForParameter(verbose,file,"ProcPEPsi", 0, 1, 1, int_type, &(config::ProcPEPsi), 1, critical);
 | 
|---|
| 474 | 
 | 
|---|
| 475 |         if (!ParseForParameter(verbose,file,"Seed", 0, 1, 1, int_type, &(config::Seed), 1, optional))
 | 
|---|
| 476 |                 config::Seed = 1;
 | 
|---|
| 477 | 
 | 
|---|
| 478 |         if(!ParseForParameter(verbose,file,"DoOutOrbitals", 0, 1, 1, int_type, &(config::DoOutOrbitals), 1, optional)) {
 | 
|---|
| 479 |                 config::DoOutOrbitals = 0;
 | 
|---|
| 480 |         } else {
 | 
|---|
| 481 |                 if (config::DoOutOrbitals < 0) config::DoOutOrbitals = 0;
 | 
|---|
| 482 |                 if (config::DoOutOrbitals > 1) config::DoOutOrbitals = 1;
 | 
|---|
| 483 |         }
 | 
|---|
| 484 |         ParseForParameter(verbose,file,"DoOutVis", 0, 1, 1, int_type, &(config::DoOutVis), 1, critical);
 | 
|---|
| 485 |         if (config::DoOutVis < 0) config::DoOutVis = 0;
 | 
|---|
| 486 |         if (config::DoOutVis > 1) config::DoOutVis = 1;
 | 
|---|
| 487 |         if (!ParseForParameter(verbose,file,"VectorPlane", 0, 1, 1, int_type, &(config::VectorPlane), 1, optional))
 | 
|---|
| 488 |                 config::VectorPlane = -1;
 | 
|---|
| 489 |         if (!ParseForParameter(verbose,file,"VectorCut", 0, 1, 1, double_type, &(config::VectorCut), 1, optional))
 | 
|---|
| 490 |                 config::VectorCut = 0.;
 | 
|---|
| 491 |         ParseForParameter(verbose,file,"DoOutMes", 0, 1, 1, int_type, &(config::DoOutMes), 1, critical);
 | 
|---|
| 492 |         if (config::DoOutMes < 0) config::DoOutMes = 0;
 | 
|---|
| 493 |         if (config::DoOutMes > 1) config::DoOutMes = 1;
 | 
|---|
| 494 |         if (!ParseForParameter(verbose,file,"DoOutCurr", 0, 1, 1, int_type, &(config::DoOutCurrent), 1, optional))
 | 
|---|
| 495 |                 config::DoOutCurrent = 0;
 | 
|---|
| 496 |         if (config::DoOutCurrent < 0) config::DoOutCurrent = 0;
 | 
|---|
| 497 |         if (config::DoOutCurrent > 1) config::DoOutCurrent = 1;
 | 
|---|
| 498 |         ParseForParameter(verbose,file,"AddGramSch", 0, 1, 1, int_type, &(config::UseAddGramSch), 1, critical);
 | 
|---|
| 499 |         if (config::UseAddGramSch < 0) config::UseAddGramSch = 0;
 | 
|---|
| 500 |         if (config::UseAddGramSch > 2) config::UseAddGramSch = 2;
 | 
|---|
| 501 |         if(!ParseForParameter(verbose,file,"DoWannier", 0, 1, 1, int_type, &(config::DoWannier), 1, optional)) {
 | 
|---|
| 502 |                 config::DoWannier = 0;
 | 
|---|
| 503 |         } else {
 | 
|---|
| 504 |                 if (config::DoWannier < 0) config::DoWannier = 0;
 | 
|---|
| 505 |                 if (config::DoWannier > 1) config::DoWannier = 1;
 | 
|---|
| 506 |         }
 | 
|---|
| 507 |         if(!ParseForParameter(verbose,file,"CommonWannier", 0, 1, 1, int_type, &(config::CommonWannier), 1, optional)) {
 | 
|---|
| 508 |                 config::CommonWannier = 0;
 | 
|---|
| 509 |         } else {
 | 
|---|
| 510 |                 if (config::CommonWannier < 0) config::CommonWannier = 0;
 | 
|---|
| 511 |                 if (config::CommonWannier > 4) config::CommonWannier = 4;
 | 
|---|
| 512 |         }
 | 
|---|
| 513 |         if(!ParseForParameter(verbose,file,"SawtoothStart", 0, 1, 1, double_type, &(config::SawtoothStart), 1, optional)) {
 | 
|---|
| 514 |                 config::SawtoothStart = 0.01;
 | 
|---|
| 515 |         } else {
 | 
|---|
| 516 |                 if (config::SawtoothStart < 0.) config::SawtoothStart = 0.;
 | 
|---|
| 517 |                 if (config::SawtoothStart > 1.) config::SawtoothStart = 1.;
 | 
|---|
| 518 |         }
 | 
|---|
| 519 | 
 | 
|---|
| 520 |         ParseForParameter(verbose,file,"MaxOuterStep", 0, 1, 1, int_type, &(config::MaxOuterStep), 1, critical);
 | 
|---|
| 521 |         if (!ParseForParameter(verbose,file,"Deltat", 0, 1, 1, double_type, &(config::Deltat), 1, optional))
 | 
|---|
| 522 |                 config::Deltat = 1;
 | 
|---|
| 523 |         ParseForParameter(verbose,file,"OutVisStep", 0, 1, 1, int_type, &(config::OutVisStep), 1, optional);
 | 
|---|
| 524 |         ParseForParameter(verbose,file,"OutSrcStep", 0, 1, 1, int_type, &(config::OutSrcStep), 1, optional);
 | 
|---|
| 525 |         ParseForParameter(verbose,file,"TargetTemp", 0, 1, 1, double_type, &(config::TargetTemp), 1, optional);
 | 
|---|
| 526 |         //ParseForParameter(verbose,file,"Thermostat", 0, 1, 1, int_type, &(config::ScaleTempStep), 1, optional);
 | 
|---|
| 527 |         if (!ParseForParameter(verbose,file,"EpsWannier", 0, 1, 1, double_type, &(config::EpsWannier), 1, optional))
 | 
|---|
| 528 |                 config::EpsWannier = 1e-8;
 | 
|---|
| 529 | 
 | 
|---|
| 530 |         // stop conditions
 | 
|---|
| 531 |         //if (config::MaxOuterStep <= 0) config::MaxOuterStep = 1;
 | 
|---|
| 532 |         ParseForParameter(verbose,file,"MaxPsiStep", 0, 1, 1, int_type, &(config::MaxPsiStep), 1, critical);
 | 
|---|
| 533 |         if (config::MaxPsiStep <= 0) config::MaxPsiStep = 3;
 | 
|---|
| 534 | 
 | 
|---|
| 535 |         ParseForParameter(verbose,file,"MaxMinStep", 0, 1, 1, int_type, &(config::MaxMinStep), 1, critical);
 | 
|---|
| 536 |         ParseForParameter(verbose,file,"RelEpsTotalE", 0, 1, 1, double_type, &(config::RelEpsTotalEnergy), 1, critical);
 | 
|---|
| 537 |         ParseForParameter(verbose,file,"RelEpsKineticE", 0, 1, 1, double_type, &(config::RelEpsKineticEnergy), 1, critical);
 | 
|---|
| 538 |         ParseForParameter(verbose,file,"MaxMinStopStep", 0, 1, 1, int_type, &(config::MaxMinStopStep), 1, critical);
 | 
|---|
| 539 |         ParseForParameter(verbose,file,"MaxMinGapStopStep", 0, 1, 1, int_type, &(config::MaxMinGapStopStep), 1, critical);
 | 
|---|
| 540 |         if (config::MaxMinStep <= 0) config::MaxMinStep = config::MaxPsiStep;
 | 
|---|
| 541 |         if (config::MaxMinStopStep < 1) config::MaxMinStopStep = 1;
 | 
|---|
| 542 |         if (config::MaxMinGapStopStep < 1) config::MaxMinGapStopStep = 1;
 | 
|---|
| 543 | 
 | 
|---|
| 544 |         ParseForParameter(verbose,file,"MaxInitMinStep", 0, 1, 1, int_type, &(config::MaxInitMinStep), 1, critical);
 | 
|---|
| 545 |         ParseForParameter(verbose,file,"InitRelEpsTotalE", 0, 1, 1, double_type, &(config::InitRelEpsTotalEnergy), 1, critical);
 | 
|---|
| 546 |         ParseForParameter(verbose,file,"InitRelEpsKineticE", 0, 1, 1, double_type, &(config::InitRelEpsKineticEnergy), 1, critical);
 | 
|---|
| 547 |         ParseForParameter(verbose,file,"InitMaxMinStopStep", 0, 1, 1, int_type, &(config::InitMaxMinStopStep), 1, critical);
 | 
|---|
| 548 |         ParseForParameter(verbose,file,"InitMaxMinGapStopStep", 0, 1, 1, int_type, &(config::InitMaxMinGapStopStep), 1, critical);
 | 
|---|
| 549 |         if (config::MaxInitMinStep <= 0) config::MaxInitMinStep = config::MaxPsiStep;
 | 
|---|
| 550 |         if (config::InitMaxMinStopStep < 1) config::InitMaxMinStopStep = 1;
 | 
|---|
| 551 |         if (config::InitMaxMinGapStopStep < 1) config::InitMaxMinGapStopStep = 1;
 | 
|---|
| 552 | 
 | 
|---|
| 553 |         // Unit cell and magnetic field
 | 
|---|
| 554 |         ParseForParameter(verbose,file, "BoxLength", 0, 3, 3, lower_trigrid, BoxLength, 1, critical); /* Lattice->RealBasis */
 | 
|---|
| 555 |         mol->cell_size[0] = BoxLength[0];
 | 
|---|
| 556 |         mol->cell_size[1] = BoxLength[3];
 | 
|---|
| 557 |         mol->cell_size[2] = BoxLength[4];
 | 
|---|
| 558 |         mol->cell_size[3] = BoxLength[6];
 | 
|---|
| 559 |         mol->cell_size[4] = BoxLength[7];
 | 
|---|
| 560 |         mol->cell_size[5] = BoxLength[8];
 | 
|---|
| 561 |         if (1) fprintf(stderr,"\n");
 | 
|---|
| 562 | 
 | 
|---|
| 563 |         ParseForParameter(verbose,file,"DoPerturbation", 0, 1, 1, int_type, &(config::DoPerturbation), 1, optional);
 | 
|---|
| 564 |         ParseForParameter(verbose,file,"DoOutNICS", 0, 1, 1, int_type, &(config::DoOutNICS), 1, optional);
 | 
|---|
| 565 |         if (!ParseForParameter(verbose,file,"DoFullCurrent", 0, 1, 1, int_type, &(config::DoFullCurrent), 1, optional))
 | 
|---|
| 566 |                 config::DoFullCurrent = 0;
 | 
|---|
| 567 |         if (config::DoFullCurrent < 0) config::DoFullCurrent = 0;
 | 
|---|
| 568 |         if (config::DoFullCurrent > 2) config::DoFullCurrent = 2;
 | 
|---|
| 569 |         if (config::DoOutNICS < 0) config::DoOutNICS = 0;
 | 
|---|
| 570 |         if (config::DoOutNICS > 2) config::DoOutNICS = 2;
 | 
|---|
| 571 |         if (config::DoPerturbation == 0) {
 | 
|---|
| 572 |                 config::DoFullCurrent = 0;
 | 
|---|
| 573 |                 config::DoOutNICS = 0;
 | 
|---|
| 574 |         }
 | 
|---|
| 575 | 
 | 
|---|
| 576 |         ParseForParameter(verbose,file,"ECut", 0, 1, 1, double_type, &(config::ECut), 1, critical);
 | 
|---|
| 577 |         ParseForParameter(verbose,file,"MaxLevel", 0, 1, 1, int_type, &(config::MaxLevel), 1, critical);
 | 
|---|
| 578 |         ParseForParameter(verbose,file,"Level0Factor", 0, 1, 1, int_type, &(config::Lev0Factor), 1, critical);
 | 
|---|
| 579 |         if (config::Lev0Factor < 2) {
 | 
|---|
| 580 |                 config::Lev0Factor = 2;
 | 
|---|
| 581 |         }
 | 
|---|
| 582 |         ParseForParameter(verbose,file,"RiemannTensor", 0, 1, 1, int_type, &di, 1, critical);
 | 
|---|
| 583 |         if (di >= 0 && di < 2) {
 | 
|---|
| 584 |                 config::RiemannTensor = di;
 | 
|---|
| 585 |         } else {
 | 
|---|
| 586 |                 fprintf(stderr, "0 <= RiemanTensor < 2: 0 UseNotRT, 1 UseRT");
 | 
|---|
| 587 |                 exit(1);
 | 
|---|
| 588 |         }
 | 
|---|
| 589 |         switch (config::RiemannTensor) {
 | 
|---|
| 590 |                 case 0: //UseNoRT
 | 
|---|
| 591 |                         if (config::MaxLevel < 2) {
 | 
|---|
| 592 |                                 config::MaxLevel = 2;
 | 
|---|
| 593 |                         }
 | 
|---|
| 594 |                         config::LevRFactor = 2;
 | 
|---|
| 595 |                         config::RTActualUse = 0;
 | 
|---|
| 596 |                         break;
 | 
|---|
| 597 |                 case 1: // UseRT
 | 
|---|
| 598 |                         if (config::MaxLevel < 3) {
 | 
|---|
| 599 |                                 config::MaxLevel = 3;
 | 
|---|
| 600 |                         }
 | 
|---|
| 601 |                         ParseForParameter(verbose,file,"RiemannLevel", 0, 1, 1, int_type, &(config::RiemannLevel), 1, critical);
 | 
|---|
| 602 |                         if (config::RiemannLevel < 2) {
 | 
|---|
| 603 |                                 config::RiemannLevel = 2;
 | 
|---|
| 604 |                         }
 | 
|---|
| 605 |                         if (config::RiemannLevel > config::MaxLevel-1) {
 | 
|---|
| 606 |                                 config::RiemannLevel = config::MaxLevel-1;
 | 
|---|
| 607 |                         }
 | 
|---|
| 608 |                         ParseForParameter(verbose,file,"LevRFactor", 0, 1, 1, int_type, &(config::LevRFactor), 1, critical);
 | 
|---|
| 609 |                         if (config::LevRFactor < 2) {
 | 
|---|
| 610 |                                 config::LevRFactor = 2;
 | 
|---|
| 611 |                         }
 | 
|---|
| 612 |                         config::Lev0Factor = 2;
 | 
|---|
| 613 |                         config::RTActualUse = 2;
 | 
|---|
| 614 |                         break;
 | 
|---|
| 615 |         }
 | 
|---|
| 616 |         ParseForParameter(verbose,file,"PsiType", 0, 1, 1, int_type, &di, 1, critical);
 | 
|---|
| 617 |         if (di >= 0 && di < 2) {
 | 
|---|
| 618 |                 config::PsiType = di;
 | 
|---|
| 619 |         } else {
 | 
|---|
| 620 |                 fprintf(stderr, "0 <= PsiType < 2: 0 UseSpinDouble, 1 UseSpinUpDown");
 | 
|---|
| 621 |                 exit(1);
 | 
|---|
| 622 |         }
 | 
|---|
| 623 |         switch (config::PsiType) {
 | 
|---|
| 624 |         case 0: // SpinDouble
 | 
|---|
| 625 |                 ParseForParameter(verbose,file,"MaxPsiDouble", 0, 1, 1, int_type, &(config::MaxPsiDouble), 1, critical);
 | 
|---|
| 626 |                 ParseForParameter(verbose,file,"AddPsis", 0, 1, 1, int_type, &(config::AddPsis), 1, optional);
 | 
|---|
| 627 |                 break;
 | 
|---|
| 628 |         case 1: // SpinUpDown
 | 
|---|
| 629 |                 if (config::ProcPEGamma % 2) config::ProcPEGamma*=2;
 | 
|---|
| 630 |                 ParseForParameter(verbose,file,"PsiMaxNoUp", 0, 1, 1, int_type, &(config::PsiMaxNoUp), 1, critical);
 | 
|---|
| 631 |                 ParseForParameter(verbose,file,"PsiMaxNoDown", 0, 1, 1, int_type, &(config::PsiMaxNoDown), 1, critical);
 | 
|---|
| 632 |                 ParseForParameter(verbose,file,"AddPsis", 0, 1, 1, int_type, &(config::AddPsis), 1, optional);
 | 
|---|
| 633 |                 break;
 | 
|---|
| 634 |         }
 | 
|---|
| 635 | 
 | 
|---|
| 636 |         // IonsInitRead
 | 
|---|
| 637 | 
 | 
|---|
| 638 |         ParseForParameter(verbose,file,"RCut", 0, 1, 1, double_type, &(config::RCut), 1, critical);
 | 
|---|
| 639 |         ParseForParameter(verbose,file,"IsAngstroem", 0, 1, 1, int_type, &(config::IsAngstroem), 1, critical);
 | 
|---|
| 640 |         ParseForParameter(verbose,file,"MaxTypes", 0, 1, 1, int_type, &(config::MaxTypes), 1, critical);
 | 
|---|
| 641 |         if (!ParseForParameter(verbose,file,"RelativeCoord", 0, 1, 1, int_type, &(config::RelativeCoord) , 1, optional))
 | 
|---|
| 642 |                 config::RelativeCoord = 0;
 | 
|---|
| 643 |         if (!ParseForParameter(verbose,file,"StructOpt", 0, 1, 1, int_type, &(config::StructOpt), 1, optional))
 | 
|---|
| 644 |                 config::StructOpt = 0;
 | 
|---|
| 645 |         if (MaxTypes == 0) {
 | 
|---|
| 646 |                 cerr << "There are no atoms according to MaxTypes in this config file." << endl;
 | 
|---|
| 647 |         } else {
 | 
|---|
| 648 |                 // prescan number of ions per type
 | 
|---|
| 649 |                 cout << Verbose(0) << "Prescanning ions per type: " << endl;
 | 
|---|
| 650 |                 for (int i=0; i < config::MaxTypes; i++) {
 | 
|---|
| 651 |                         sprintf(name,"Ion_Type%i",i+1);
 | 
|---|
| 652 |                         ParseForParameter(verbose,file, (const char*)name, 0, 1, 1, int_type, &No[i], 1, critical);
 | 
|---|
| 653 |                         ParseForParameter(verbose,file, name, 0, 2, 1, int_type, &Z, 1, critical);
 | 
|---|
| 654 |                         elementhash[i] = periode->FindElement(Z);
 | 
|---|
| 655 |                         cout << Verbose(1) << i << ". Z = " << elementhash[i]->Z << " with " << No[i] << " ions." << endl;
 | 
|---|
| 656 |                 }
 | 
|---|
| 657 |                 int repetition = 0; // which repeated keyword shall be read
 | 
|---|
| 658 | 
 | 
|---|
| 659 |                 map<int, atom *> AtomList[config::MaxTypes];
 | 
|---|
| 660 |                 if (!FastParsing) {
 | 
|---|
| 661 |                         // parse in trajectories
 | 
|---|
| 662 |                         bool status = true;
 | 
|---|
| 663 |                         atom *neues = NULL;
 | 
|---|
| 664 |                         while (status) {
 | 
|---|
| 665 |                                 cout << "Currently parsing MD step " << repetition << "." << endl;
 | 
|---|
| 666 |                                 for (int i=0; i < config::MaxTypes; i++) {
 | 
|---|
| 667 |                                         sprintf(name,"Ion_Type%i",i+1);
 | 
|---|
| 668 |                                         for(int j=0;j<No[i];j++) {
 | 
|---|
| 669 |                                                 sprintf(keyword,"%s_%i",name, j+1);
 | 
|---|
| 670 |                                                 if (repetition == 0) {
 | 
|---|
| 671 |                                                         neues = new atom();
 | 
|---|
| 672 |                                                         AtomList[i][j] = neues;
 | 
|---|
| 673 |                                                         neues->type = elementhash[i]; // find element type
 | 
|---|
| 674 |                                                         mol->AddAtom(neues);
 | 
|---|
| 675 |                                                 } else
 | 
|---|
| 676 |                                                         neues = AtomList[i][j];
 | 
|---|
| 677 |                                                 status = (status &&
 | 
|---|
| 678 |                                                                                 ParseForParameter(verbose,file, keyword, 0, 1, 1, double_type, &neues->x.x[0], 1, (repetition == 0) ? critical : optional) &&
 | 
|---|
| 679 |                                                                                 ParseForParameter(verbose,file, keyword, 0, 2, 1, double_type, &neues->x.x[1], 1, (repetition == 0) ? critical : optional) &&
 | 
|---|
| 680 |                                                                                 ParseForParameter(verbose,file, keyword, 0, 3, 1, double_type, &neues->x.x[2], 1, (repetition == 0) ? critical : optional) &&
 | 
|---|
| 681 |                                                                                 ParseForParameter(verbose,file, keyword, 0, 4, 1, int_type, &neues->FixedIon, 1, (repetition == 0) ? critical : optional));
 | 
|---|
| 682 |                                                 if (!status) break;
 | 
|---|
| 683 | 
 | 
|---|
| 684 |                                                 // check size of vectors
 | 
|---|
| 685 |                                                 if (mol->Trajectories[neues].R.size() <= (unsigned int)(repetition)) {
 | 
|---|
| 686 |                                                         //cout << "Increasing size for trajectory array of " << keyword << " to " << (repetition+10) << "." << endl;
 | 
|---|
| 687 |                                                         mol->Trajectories[neues].R.resize(repetition+10);
 | 
|---|
| 688 |                                                         mol->Trajectories[neues].U.resize(repetition+10);
 | 
|---|
| 689 |                                                         mol->Trajectories[neues].F.resize(repetition+10);
 | 
|---|
| 690 |                                                 }
 | 
|---|
| 691 | 
 | 
|---|
| 692 |                                                 // put into trajectories list
 | 
|---|
| 693 |                                                 for (int d=0;d<NDIM;d++)
 | 
|---|
| 694 |                                                         mol->Trajectories[neues].R.at(repetition).x[d] = neues->x.x[d];
 | 
|---|
| 695 | 
 | 
|---|
| 696 |                                                 // parse velocities if present
 | 
|---|
| 697 |                                                 if(!ParseForParameter(verbose,file, keyword, 0, 5, 1, double_type, &neues->v.x[0], 1,optional))
 | 
|---|
| 698 |                                                         neues->v.x[0] = 0.;
 | 
|---|
| 699 |                                                 if(!ParseForParameter(verbose,file, keyword, 0, 6, 1, double_type, &neues->v.x[1], 1,optional))
 | 
|---|
| 700 |                                                         neues->v.x[1] = 0.;
 | 
|---|
| 701 |                                                 if(!ParseForParameter(verbose,file, keyword, 0, 7, 1, double_type, &neues->v.x[2], 1,optional))
 | 
|---|
| 702 |                                                         neues->v.x[2] = 0.;
 | 
|---|
| 703 |                                                 for (int d=0;d<NDIM;d++)
 | 
|---|
| 704 |                                                         mol->Trajectories[neues].U.at(repetition).x[d] = neues->v.x[d];
 | 
|---|
| 705 | 
 | 
|---|
| 706 |                                                 // parse forces if present
 | 
|---|
| 707 |                                                 if(!ParseForParameter(verbose,file, keyword, 0, 8, 1, double_type, &value[0], 1,optional))
 | 
|---|
| 708 |                                                         value[0] = 0.;
 | 
|---|
| 709 |                                                 if(!ParseForParameter(verbose,file, keyword, 0, 9, 1, double_type, &value[1], 1,optional))
 | 
|---|
| 710 |                                                         value[1] = 0.;
 | 
|---|
| 711 |                                                 if(!ParseForParameter(verbose,file, keyword, 1, 10, 1, double_type, &value[2], 1,optional))
 | 
|---|
| 712 |                                                         value[2] = 0.;
 | 
|---|
| 713 |                                                 for (int d=0;d<NDIM;d++)
 | 
|---|
| 714 |                                                         mol->Trajectories[neues].F.at(repetition).x[d] = value[d];
 | 
|---|
| 715 | 
 | 
|---|
| 716 |         //                                              cout << "Parsed position of step " << (repetition) << ": (";
 | 
|---|
| 717 |         //                                              for (int d=0;d<NDIM;d++)
 | 
|---|
| 718 |         //                                                      cout << mol->Trajectories[neues].R.at(repetition).x[d] << " ";                                  // next step
 | 
|---|
| 719 |         //                                              cout << ")\t(";
 | 
|---|
| 720 |         //                                              for (int d=0;d<NDIM;d++)
 | 
|---|
| 721 |         //                                                      cout << mol->Trajectories[neues].U.at(repetition).x[d] << " ";                                  // next step
 | 
|---|
| 722 |         //                                              cout << ")\t(";
 | 
|---|
| 723 |         //                                              for (int d=0;d<NDIM;d++)
 | 
|---|
| 724 |         //                                                      cout << mol->Trajectories[neues].F.at(repetition).x[d] << " ";                                  // next step
 | 
|---|
| 725 |         //                                              cout << ")" << endl;
 | 
|---|
| 726 |                                         }
 | 
|---|
| 727 |                                 }
 | 
|---|
| 728 |                                 repetition++;
 | 
|---|
| 729 |                         }
 | 
|---|
| 730 |                         repetition--;
 | 
|---|
| 731 |                         cout << "Found " << repetition << " trajectory steps." << endl;
 | 
|---|
| 732 |                         mol->MDSteps = repetition;
 | 
|---|
| 733 |                 } else {
 | 
|---|
| 734 |                         // find the maximum number of MD steps so that we may parse last one (Ion_Type1_1 must always be present, because is the first atom)
 | 
|---|
| 735 |                         repetition = 0;
 | 
|---|
| 736 |                         while ( ParseForParameter(verbose,file, "Ion_Type1_1", 0, 1, 1, double_type, &value[0], repetition, (repetition == 0) ? critical : optional) &&
 | 
|---|
| 737 |                                                         ParseForParameter(verbose,file, "Ion_Type1_1", 0, 2, 1, double_type, &value[1], repetition, (repetition == 0) ? critical : optional) &&
 | 
|---|
| 738 |                                                         ParseForParameter(verbose,file, "Ion_Type1_1", 0, 3, 1, double_type, &value[2], repetition, (repetition == 0) ? critical : optional))
 | 
|---|
| 739 |                                 repetition++;
 | 
|---|
| 740 |                         cout << "I found " << repetition << " times the keyword Ion_Type1_1." << endl;
 | 
|---|
| 741 |                         // parse in molecule coordinates
 | 
|---|
| 742 |                         for (int i=0; i < config::MaxTypes; i++) {
 | 
|---|
| 743 |                                 sprintf(name,"Ion_Type%i",i+1);
 | 
|---|
| 744 |                                 for(int j=0;j<No[i];j++) {
 | 
|---|
| 745 |                                         sprintf(keyword,"%s_%i",name, j+1);
 | 
|---|
| 746 |                                         atom *neues = new atom();
 | 
|---|
| 747 |                                         // then parse for each atom the coordinates as often as present
 | 
|---|
| 748 |                                         ParseForParameter(verbose,file, keyword, 0, 1, 1, double_type, &neues->x.x[0], repetition,critical);
 | 
|---|
| 749 |                                         ParseForParameter(verbose,file, keyword, 0, 2, 1, double_type, &neues->x.x[1], repetition,critical);
 | 
|---|
| 750 |                                         ParseForParameter(verbose,file, keyword, 0, 3, 1, double_type, &neues->x.x[2], repetition,critical);
 | 
|---|
| 751 |                                         ParseForParameter(verbose,file, keyword, 0, 4, 1, int_type, &neues->FixedIon, repetition,critical);
 | 
|---|
| 752 |                                         if(!ParseForParameter(verbose,file, keyword, 0, 5, 1, double_type, &neues->v.x[0], repetition,optional))
 | 
|---|
| 753 |                                                 neues->v.x[0] = 0.;
 | 
|---|
| 754 |                                         if(!ParseForParameter(verbose,file, keyword, 0, 6, 1, double_type, &neues->v.x[1], repetition,optional))
 | 
|---|
| 755 |                                                 neues->v.x[1] = 0.;
 | 
|---|
| 756 |                                         if(!ParseForParameter(verbose,file, keyword, 0, 7, 1, double_type, &neues->v.x[2], repetition,optional))
 | 
|---|
| 757 |                                                 neues->v.x[2] = 0.;
 | 
|---|
| 758 |                                         // here we don't care if forces are present (last in trajectories is always equal to current position)
 | 
|---|
| 759 |                                         neues->type = elementhash[i]; // find element type
 | 
|---|
| 760 |                                         mol->AddAtom(neues);
 | 
|---|
| 761 |                                 }
 | 
|---|
| 762 |                         }
 | 
|---|
| 763 |                 }
 | 
|---|
| 764 |         }
 | 
|---|
| 765 |         file->close();
 | 
|---|
| 766 |         delete(file);
 | 
|---|
| 767 | };
 | 
|---|
| 768 | 
 | 
|---|
| 769 | /** Initializes config file structure by loading elements from a give file with old pcp syntax.
 | 
|---|
| 770 |  * \param *file input file stream being the opened config file with old pcp syntax
 | 
|---|
| 771 |  * \param *periode pointer to a periodentafel class with all elements
 | 
|---|
| 772 |  * \param *mol pointer to molecule containing all atoms of the molecule
 | 
|---|
| 773 |  */
 | 
|---|
| 774 | void config::LoadOld(char *filename, periodentafel *periode, molecule *mol)
 | 
|---|
| 775 | {
 | 
|---|
| 776 |         ifstream *file = new ifstream(filename);
 | 
|---|
| 777 |         if (file == NULL) {
 | 
|---|
| 778 |                 cerr << "ERROR: config file " << filename << " missing!" << endl;
 | 
|---|
| 779 |                 return;
 | 
|---|
| 780 |         }
 | 
|---|
| 781 |         RetrieveConfigPathAndName(filename);
 | 
|---|
| 782 |         // ParseParameters
 | 
|---|
| 783 | 
 | 
|---|
| 784 |         /* Oeffne Hauptparameterdatei */
 | 
|---|
| 785 |         int l, i, di;
 | 
|---|
| 786 |         double a,b;
 | 
|---|
| 787 |         double BoxLength[9];
 | 
|---|
| 788 |         string zeile;
 | 
|---|
| 789 |         string dummy;
 | 
|---|
| 790 |         element *elementhash[128];
 | 
|---|
| 791 |         int Z, No, AtomNo, found;
 | 
|---|
| 792 |         int verbose = 0;
 | 
|---|
| 793 | 
 | 
|---|
| 794 |         /* Namen einlesen */
 | 
|---|
| 795 | 
 | 
|---|
| 796 |         ParseForParameter(verbose,file, "mainname", 0, 1, 1, string_type, (config::mainname), 1, critical);
 | 
|---|
| 797 |         ParseForParameter(verbose,file, "defaultpath", 0, 1, 1, string_type, (config::defaultpath), 1, critical);
 | 
|---|
| 798 |         ParseForParameter(verbose,file, "pseudopotpath", 0, 1, 1, string_type, (config::pseudopotpath), 1, critical);
 | 
|---|
| 799 |         ParseForParameter(verbose,file, "ProcsGammaPsi", 0, 1, 1, int_type, &(config::ProcPEGamma), 1, critical);
 | 
|---|
| 800 |         ParseForParameter(verbose,file, "ProcsGammaPsi", 0, 2, 1, int_type, &(config::ProcPEPsi), 1, critical);
 | 
|---|
| 801 |         config::Seed = 1;
 | 
|---|
| 802 |         config::DoOutOrbitals = 0;
 | 
|---|
| 803 |         ParseForParameter(verbose,file,"DoOutVis", 0, 1, 1, int_type, &(config::DoOutVis), 1, critical);
 | 
|---|
| 804 |         if (config::DoOutVis < 0) config::DoOutVis = 0;
 | 
|---|
| 805 |         if (config::DoOutVis > 1) config::DoOutVis = 1;
 | 
|---|
| 806 |         config::VectorPlane = -1;
 | 
|---|
| 807 |         config::VectorCut = 0.;
 | 
|---|
| 808 |         ParseForParameter(verbose,file,"DoOutMes", 0, 1, 1, int_type, &(config::DoOutMes), 1, critical);
 | 
|---|
| 809 |         if (config::DoOutMes < 0) config::DoOutMes = 0;
 | 
|---|
| 810 |         if (config::DoOutMes > 1) config::DoOutMes = 1;
 | 
|---|
| 811 |         config::DoOutCurrent = 0;
 | 
|---|
| 812 |         ParseForParameter(verbose,file,"AddGramSch", 0, 1, 1, int_type, &(config::UseAddGramSch), 1, critical);
 | 
|---|
| 813 |         if (config::UseAddGramSch < 0) config::UseAddGramSch = 0;
 | 
|---|
| 814 |         if (config::UseAddGramSch > 2) config::UseAddGramSch = 2;
 | 
|---|
| 815 |         config::CommonWannier = 0;
 | 
|---|
| 816 |         config::SawtoothStart = 0.01;
 | 
|---|
| 817 | 
 | 
|---|
| 818 |         ParseForParameter(verbose,file,"MaxOuterStep", 0, 1, 1, double_type, &(config::MaxOuterStep), 1, critical);
 | 
|---|
| 819 |         ParseForParameter(verbose,file,"Deltat", 0, 1, 1, double_type, &(config::Deltat), 1, optional);
 | 
|---|
| 820 |         ParseForParameter(verbose,file,"VisOuterStep", 0, 1, 1, int_type, &(config::OutVisStep), 1, optional);
 | 
|---|
| 821 |         ParseForParameter(verbose,file,"VisSrcOuterStep", 0, 1, 1, int_type, &(config::OutSrcStep), 1, optional);
 | 
|---|
| 822 |         ParseForParameter(verbose,file,"TargetTemp", 0, 1, 1, double_type, &(config::TargetTemp), 1, optional);
 | 
|---|
| 823 |         ParseForParameter(verbose,file,"ScaleTempStep", 0, 1, 1, int_type, &(config::ScaleTempStep), 1, optional);
 | 
|---|
| 824 |         config::EpsWannier = 1e-8;
 | 
|---|
| 825 | 
 | 
|---|
| 826 |         // stop conditions
 | 
|---|
| 827 |         //if (config::MaxOuterStep <= 0) config::MaxOuterStep = 1;
 | 
|---|
| 828 |         ParseForParameter(verbose,file,"MaxPsiStep", 0, 1, 1, int_type, &(config::MaxPsiStep), 1, critical);
 | 
|---|
| 829 |         if (config::MaxPsiStep <= 0) config::MaxPsiStep = 3;
 | 
|---|
| 830 | 
 | 
|---|
| 831 |         ParseForParameter(verbose,file,"MaxMinStep", 0, 1, 1, int_type, &(config::MaxMinStep), 1, critical);
 | 
|---|
| 832 |         ParseForParameter(verbose,file,"MaxMinStep", 0, 2, 1, double_type, &(config::RelEpsTotalEnergy), 1, critical);
 | 
|---|
| 833 |         ParseForParameter(verbose,file,"MaxMinStep", 0, 3, 1, double_type, &(config::RelEpsKineticEnergy), 1, critical);
 | 
|---|
| 834 |         ParseForParameter(verbose,file,"MaxMinStep", 0, 4, 1, int_type, &(config::MaxMinStopStep), 1, critical);
 | 
|---|
| 835 |         if (config::MaxMinStep <= 0) config::MaxMinStep = config::MaxPsiStep;
 | 
|---|
| 836 |         if (config::MaxMinStopStep < 1) config::MaxMinStopStep = 1;
 | 
|---|
| 837 |         config::MaxMinGapStopStep = 1;
 | 
|---|
| 838 | 
 | 
|---|
| 839 |         ParseForParameter(verbose,file,"MaxInitMinStep", 0, 1, 1, int_type, &(config::MaxInitMinStep), 1, critical);
 | 
|---|
| 840 |         ParseForParameter(verbose,file,"MaxInitMinStep", 0, 2, 1, double_type, &(config::InitRelEpsTotalEnergy), 1, critical);
 | 
|---|
| 841 |         ParseForParameter(verbose,file,"MaxInitMinStep", 0, 3, 1, double_type, &(config::InitRelEpsKineticEnergy), 1, critical);
 | 
|---|
| 842 |         ParseForParameter(verbose,file,"MaxInitMinStep", 0, 4, 1, int_type, &(config::InitMaxMinStopStep), 1, critical);
 | 
|---|
| 843 |         if (config::MaxInitMinStep <= 0) config::MaxInitMinStep = config::MaxPsiStep;
 | 
|---|
| 844 |         if (config::InitMaxMinStopStep < 1) config::InitMaxMinStopStep = 1;
 | 
|---|
| 845 |         config::InitMaxMinGapStopStep = 1;
 | 
|---|
| 846 | 
 | 
|---|
| 847 |         ParseForParameter(verbose,file, "BoxLength", 0, 3, 3, lower_trigrid, BoxLength, 1, critical); /* Lattice->RealBasis */
 | 
|---|
| 848 |         mol->cell_size[0] = BoxLength[0];
 | 
|---|
| 849 |         mol->cell_size[1] = BoxLength[3];
 | 
|---|
| 850 |         mol->cell_size[2] = BoxLength[4];
 | 
|---|
| 851 |         mol->cell_size[3] = BoxLength[6];
 | 
|---|
| 852 |         mol->cell_size[4] = BoxLength[7];
 | 
|---|
| 853 |         mol->cell_size[5] = BoxLength[8];
 | 
|---|
| 854 |         if (1) fprintf(stderr,"\n");
 | 
|---|
| 855 |         config::DoPerturbation = 0;
 | 
|---|
| 856 |         config::DoFullCurrent = 0;
 | 
|---|
| 857 | 
 | 
|---|
| 858 |         ParseForParameter(verbose,file,"ECut", 0, 1, 1, double_type, &(config::ECut), 1, critical);
 | 
|---|
| 859 |         ParseForParameter(verbose,file,"MaxLevel", 0, 1, 1, int_type, &(config::MaxLevel), 1, critical);
 | 
|---|
| 860 |         ParseForParameter(verbose,file,"Level0Factor", 0, 1, 1, int_type, &(config::Lev0Factor), 1, critical);
 | 
|---|
| 861 |         if (config::Lev0Factor < 2) {
 | 
|---|
| 862 |                 config::Lev0Factor = 2;
 | 
|---|
| 863 |         }
 | 
|---|
| 864 |         ParseForParameter(verbose,file,"RiemannTensor", 0, 1, 1, int_type, &di, 1, critical);
 | 
|---|
| 865 |         if (di >= 0 && di < 2) {
 | 
|---|
| 866 |                 config::RiemannTensor = di;
 | 
|---|
| 867 |         } else {
 | 
|---|
| 868 |                 fprintf(stderr, "0 <= RiemanTensor < 2: 0 UseNotRT, 1 UseRT");
 | 
|---|
| 869 |                 exit(1);
 | 
|---|
| 870 |         }
 | 
|---|
| 871 |         switch (config::RiemannTensor) {
 | 
|---|
| 872 |                 case 0: //UseNoRT
 | 
|---|
| 873 |                         if (config::MaxLevel < 2) {
 | 
|---|
| 874 |                                 config::MaxLevel = 2;
 | 
|---|
| 875 |                         }
 | 
|---|
| 876 |                         config::LevRFactor = 2;
 | 
|---|
| 877 |                         config::RTActualUse = 0;
 | 
|---|
| 878 |                         break;
 | 
|---|
| 879 |                 case 1: // UseRT
 | 
|---|
| 880 |                         if (config::MaxLevel < 3) {
 | 
|---|
| 881 |                                 config::MaxLevel = 3;
 | 
|---|
| 882 |                         }
 | 
|---|
| 883 |                         ParseForParameter(verbose,file,"RiemannLevel", 0, 1, 1, int_type, &(config::RiemannLevel), 1, critical);
 | 
|---|
| 884 |                         if (config::RiemannLevel < 2) {
 | 
|---|
| 885 |                                 config::RiemannLevel = 2;
 | 
|---|
| 886 |                         }
 | 
|---|
| 887 |                         if (config::RiemannLevel > config::MaxLevel-1) {
 | 
|---|
| 888 |                                 config::RiemannLevel = config::MaxLevel-1;
 | 
|---|
| 889 |                         }
 | 
|---|
| 890 |                         ParseForParameter(verbose,file,"LevRFactor", 0, 1, 1, int_type, &(config::LevRFactor), 1, critical);
 | 
|---|
| 891 |                         if (config::LevRFactor < 2) {
 | 
|---|
| 892 |                                 config::LevRFactor = 2;
 | 
|---|
| 893 |                         }
 | 
|---|
| 894 |                         config::Lev0Factor = 2;
 | 
|---|
| 895 |                         config::RTActualUse = 2;
 | 
|---|
| 896 |                         break;
 | 
|---|
| 897 |         }
 | 
|---|
| 898 |         ParseForParameter(verbose,file,"PsiType", 0, 1, 1, int_type, &di, 1, critical);
 | 
|---|
| 899 |         if (di >= 0 && di < 2) {
 | 
|---|
| 900 |                 config::PsiType = di;
 | 
|---|
| 901 |         } else {
 | 
|---|
| 902 |                 fprintf(stderr, "0 <= PsiType < 2: 0 UseSpinDouble, 1 UseSpinUpDown");
 | 
|---|
| 903 |                 exit(1);
 | 
|---|
| 904 |         }
 | 
|---|
| 905 |         switch (config::PsiType) {
 | 
|---|
| 906 |         case 0: // SpinDouble
 | 
|---|
| 907 |                 ParseForParameter(verbose,file,"MaxPsiDouble", 0, 1, 1, int_type, &(config::MaxPsiDouble), 1, critical);
 | 
|---|
| 908 |                 config::AddPsis = 0;
 | 
|---|
| 909 |                 break;
 | 
|---|
| 910 |         case 1: // SpinUpDown
 | 
|---|
| 911 |                 if (config::ProcPEGamma % 2) config::ProcPEGamma*=2;
 | 
|---|
| 912 |                 ParseForParameter(verbose,file,"MaxPsiUp", 0, 1, 1, int_type, &(config::PsiMaxNoUp), 1, critical);
 | 
|---|
| 913 |                 ParseForParameter(verbose,file,"MaxPsiDown", 0, 1, 1, int_type, &(config::PsiMaxNoDown), 1, critical);
 | 
|---|
| 914 |                 config::AddPsis = 0;
 | 
|---|
| 915 |                 break;
 | 
|---|
| 916 |         }
 | 
|---|
| 917 | 
 | 
|---|
| 918 |         // IonsInitRead
 | 
|---|
| 919 | 
 | 
|---|
| 920 |         ParseForParameter(verbose,file,"RCut", 0, 1, 1, double_type, &(config::RCut), 1, critical);
 | 
|---|
| 921 |         ParseForParameter(verbose,file,"IsAngstroem", 0, 1, 1, int_type, &(config::IsAngstroem), 1, critical);
 | 
|---|
| 922 |         config::RelativeCoord = 0;
 | 
|---|
| 923 |         config::StructOpt = 0;
 | 
|---|
| 924 | 
 | 
|---|
| 925 |         // Routine from builder.cpp
 | 
|---|
| 926 | 
 | 
|---|
| 927 | 
 | 
|---|
| 928 |         for (i=MAX_ELEMENTS;i--;)
 | 
|---|
| 929 |                 elementhash[i] = NULL;
 | 
|---|
| 930 |         cout << Verbose(0) << "Parsing Ions ..." << endl;
 | 
|---|
| 931 |         No=0;
 | 
|---|
| 932 |         found = 0;
 | 
|---|
| 933 |         while (getline(*file,zeile,'\n')) {
 | 
|---|
| 934 |                 if (zeile.find("Ions_Data") == 0) {
 | 
|---|
| 935 |                         cout << Verbose(1) << "found Ions_Data...begin parsing" << endl;
 | 
|---|
| 936 |                         found ++;
 | 
|---|
| 937 |                 }
 | 
|---|
| 938 |                 if (found > 0) {
 | 
|---|
| 939 |                         if (zeile.find("Ions_Data") == 0)
 | 
|---|
| 940 |                                 getline(*file,zeile,'\n'); // read next line and parse this one
 | 
|---|
| 941 |                         istringstream input(zeile);
 | 
|---|
| 942 |                         input >> AtomNo;        // number of atoms
 | 
|---|
| 943 |                         input >> Z;                      // atomic number
 | 
|---|
| 944 |                         input >> a;
 | 
|---|
| 945 |                         input >> l;
 | 
|---|
| 946 |                         input >> l;
 | 
|---|
| 947 |                         input >> b;              // element mass
 | 
|---|
| 948 |                         elementhash[No] = periode->FindElement(Z);
 | 
|---|
| 949 |                         cout << Verbose(1) << "AtomNo: " << AtomNo << "\tZ: " << Z << "\ta:" << a << "\tl:"     << l << "\b:" << b << "\tElement:" << elementhash[No] << "\t:" << endl;
 | 
|---|
| 950 |                         for(i=0;i<AtomNo;i++) {
 | 
|---|
| 951 |                                 if (!getline(*file,zeile,'\n')) {// parse on and on
 | 
|---|
| 952 |                                         cout << Verbose(2) << "Error: Too few items in ionic list of element" << elementhash[No] << "." << endl << "Exiting." << endl;
 | 
|---|
| 953 |                                         // return 1;
 | 
|---|
| 954 |                                 } else {
 | 
|---|
| 955 |                                         //cout << Verbose(2) << "Reading line: " << zeile << endl;
 | 
|---|
| 956 |                                 }
 | 
|---|
| 957 |                                 istringstream input2(zeile);
 | 
|---|
| 958 |                                 atom *neues = new atom();
 | 
|---|
| 959 |                                 input2 >> neues->x.x[0]; // x
 | 
|---|
| 960 |                                 input2 >> neues->x.x[1]; // y
 | 
|---|
| 961 |                                 input2 >> neues->x.x[2]; // z
 | 
|---|
| 962 |                                 input2 >> l;
 | 
|---|
| 963 |                                 neues->type = elementhash[No]; // find element type
 | 
|---|
| 964 |                                 mol->AddAtom(neues);
 | 
|---|
| 965 |                         }
 | 
|---|
| 966 |                         No++;
 | 
|---|
| 967 |                 }
 | 
|---|
| 968 |         }
 | 
|---|
| 969 |         file->close();
 | 
|---|
| 970 |         delete(file);
 | 
|---|
| 971 | };
 | 
|---|
| 972 | 
 | 
|---|
| 973 | /** Stores all elements of config structure from which they can be re-read.
 | 
|---|
| 974 |  * \param *filename name of file
 | 
|---|
| 975 |  * \param *periode pointer to a periodentafel class with all elements
 | 
|---|
| 976 |  * \param *mol pointer to molecule containing all atoms of the molecule
 | 
|---|
| 977 |  */
 | 
|---|
| 978 | bool config::Save(const char *filename, periodentafel *periode, molecule *mol) const
 | 
|---|
| 979 | {
 | 
|---|
| 980 |         bool result = true;
 | 
|---|
| 981 |         // bring MaxTypes up to date
 | 
|---|
| 982 |         mol->CountElements();
 | 
|---|
| 983 |         ofstream *output = NULL;
 | 
|---|
| 984 |         output = new ofstream(filename, ios::out);
 | 
|---|
| 985 |         if (output != NULL) {
 | 
|---|
| 986 |                 *output << "# ParallelCarParinello - main configuration file - created with molecuilder" << endl;
 | 
|---|
| 987 |                 *output << endl;
 | 
|---|
| 988 |                 *output << "mainname\t" << config::mainname << "\t# programm name (for runtime files)" << endl;
 | 
|---|
| 989 |                 *output << "defaultpath\t" << config::defaultpath << "\t# where to put files during runtime" << endl;
 | 
|---|
| 990 |                 *output << "pseudopotpath\t" << config::pseudopotpath << "\t# where to find pseudopotentials" << endl;
 | 
|---|
| 991 |                 *output << endl;
 | 
|---|
| 992 |                 *output << "ProcPEGamma\t" << config::ProcPEGamma << "\t# for parallel computing: share constants" << endl;
 | 
|---|
| 993 |                 *output << "ProcPEPsi\t" << config::ProcPEPsi << "\t# for parallel computing: share wave functions" << endl;
 | 
|---|
| 994 |                 *output << "DoOutVis\t" << config::DoOutVis << "\t# Output data for OpenDX" << endl;
 | 
|---|
| 995 |                 *output << "DoOutMes\t" << config::DoOutMes << "\t# Output data for measurements" << endl;
 | 
|---|
| 996 |                 *output << "DoOutOrbitals\t" << config::DoOutOrbitals << "\t# Output all Orbitals" << endl;
 | 
|---|
| 997 |                 *output << "DoOutCurr\t" << config::DoOutCurrent << "\t# Ouput current density for OpenDx" << endl;
 | 
|---|
| 998 |                 *output << "DoOutNICS\t" << config::DoOutNICS << "\t# Output Nucleus independent current shieldings" << endl;
 | 
|---|
| 999 |                 *output << "DoPerturbation\t" << config::DoPerturbation << "\t# Do perturbation calculate and determine susceptibility and shielding" << endl;
 | 
|---|
| 1000 |                 *output << "DoFullCurrent\t" << config::DoFullCurrent << "\t# Do full perturbation" << endl;
 | 
|---|
| 1001 |                 *output << "CommonWannier\t" << config::CommonWannier << "\t# Put virtual centers at indivual orbits, all common, merged by variance, to grid point, to cell center" << endl;
 | 
|---|
| 1002 |                 *output << "SawtoothStart\t" << config::SawtoothStart << "\t# Absolute value for smooth transition at cell border " << endl;
 | 
|---|
| 1003 |                 *output << "VectorPlane\t" << config::VectorPlane << "\t# Cut plane axis (x, y or z: 0,1,2) for two-dim current vector plot" << endl;
 | 
|---|
| 1004 |                 *output << "VectorCut\t" << config::VectorCut << "\t# Cut plane axis value" << endl;
 | 
|---|
| 1005 |                 *output << "AddGramSch\t" << config::UseAddGramSch << "\t# Additional GramSchmidtOrtogonalization to be safe" << endl;
 | 
|---|
| 1006 |                 *output << "Seed\t\t" << config::Seed << "\t# initial value for random seed for Psi coefficients" << endl;
 | 
|---|
| 1007 |                 *output << endl;
 | 
|---|
| 1008 |                 *output << "MaxOuterStep\t" << config::MaxOuterStep << "\t# number of MolecularDynamics/Structure optimization steps" << endl;
 | 
|---|
| 1009 |                 *output << "Deltat\t" << config::Deltat << "\t# time per MD step" << endl;
 | 
|---|
| 1010 |                 *output << "OutVisStep\t" << config::OutVisStep << "\t# Output visual data every ...th step" << endl;
 | 
|---|
| 1011 |                 *output << "OutSrcStep\t" << config::OutSrcStep << "\t# Output \"restart\" data every ..th step" << endl;
 | 
|---|
| 1012 |                 *output << "TargetTemp\t" << config::TargetTemp << "\t# Target temperature" << endl;
 | 
|---|
| 1013 |                 *output << "MaxPsiStep\t" << config::MaxPsiStep << "\t# number of Minimisation steps per state (0 - default)" << endl;
 | 
|---|
| 1014 |                 *output << "EpsWannier\t" << config::EpsWannier << "\t# tolerance value for spread minimisation of orbitals" << endl;
 | 
|---|
| 1015 |                 *output << endl;
 | 
|---|
| 1016 |                 *output << "# Values specifying when to stop" << endl;
 | 
|---|
| 1017 |                 *output << "MaxMinStep\t" << config::MaxMinStep << "\t# Maximum number of steps" << endl;
 | 
|---|
| 1018 |                 *output << "RelEpsTotalE\t" << config::RelEpsTotalEnergy << "\t# relative change in total energy" << endl;
 | 
|---|
| 1019 |                 *output << "RelEpsKineticE\t" << config::RelEpsKineticEnergy << "\t# relative change in kinetic energy" << endl;
 | 
|---|
| 1020 |                 *output << "MaxMinStopStep\t" << config::MaxMinStopStep << "\t# check every ..th steps" << endl;
 | 
|---|
| 1021 |                 *output << "MaxMinGapStopStep\t" << config::MaxMinGapStopStep << "\t# check every ..th steps" << endl;
 | 
|---|
| 1022 |                 *output << endl;
 | 
|---|
| 1023 |                 *output << "# Values specifying when to stop for INIT, otherwise same as above" << endl;
 | 
|---|
| 1024 |                 *output << "MaxInitMinStep\t" << config::MaxInitMinStep << "\t# Maximum number of steps" << endl;
 | 
|---|
| 1025 |                 *output << "InitRelEpsTotalE\t" << config::InitRelEpsTotalEnergy << "\t# relative change in total energy" << endl;
 | 
|---|
| 1026 |                 *output << "InitRelEpsKineticE\t" << config::InitRelEpsKineticEnergy << "\t# relative change in kinetic energy" << endl;
 | 
|---|
| 1027 |                 *output << "InitMaxMinStopStep\t" << config::InitMaxMinStopStep << "\t# check every ..th steps" << endl;
 | 
|---|
| 1028 |                 *output << "InitMaxMinGapStopStep\t" << config::InitMaxMinGapStopStep << "\t# check every ..th steps" << endl;
 | 
|---|
| 1029 |                 *output << endl;
 | 
|---|
| 1030 |                 *output << "BoxLength\t\t\t# (Length of a unit cell)" << endl;
 | 
|---|
| 1031 |                 *output << mol->cell_size[0] << "\t" << endl;
 | 
|---|
| 1032 |                 *output << mol->cell_size[1] << "\t" << mol->cell_size[2] << "\t" << endl;
 | 
|---|
| 1033 |                 *output << mol->cell_size[3] << "\t" << mol->cell_size[4] << "\t" << mol->cell_size[5] << "\t" << endl;
 | 
|---|
| 1034 |                 // FIXME
 | 
|---|
| 1035 |                 *output << endl;
 | 
|---|
| 1036 |                 *output << "ECut\t\t" << config::ECut << "\t# energy cutoff for discretization in Hartrees" << endl;
 | 
|---|
| 1037 |                 *output << "MaxLevel\t" << config::MaxLevel << "\t# number of different levels in the code, >=2" << endl;
 | 
|---|
| 1038 |                 *output << "Level0Factor\t" << config::Lev0Factor << "\t# factor by which node number increases from S to 0 level" << endl;
 | 
|---|
| 1039 |                 *output << "RiemannTensor\t" << config::RiemannTensor << "\t# (Use metric)" << endl;
 | 
|---|
| 1040 |                 switch (config::RiemannTensor) {
 | 
|---|
| 1041 |                         case 0: //UseNoRT
 | 
|---|
| 1042 |                                 break;
 | 
|---|
| 1043 |                         case 1: // UseRT
 | 
|---|
| 1044 |                                 *output << "RiemannLevel\t" << config::RiemannLevel << "\t# Number of Riemann Levels" << endl;
 | 
|---|
| 1045 |                                 *output << "LevRFactor\t" << config::LevRFactor << "\t# factor by which node number increases from 0 to R level from" << endl;
 | 
|---|
| 1046 |                                 break;
 | 
|---|
| 1047 |                 }
 | 
|---|
| 1048 |                 *output << "PsiType\t\t" << config::PsiType << "\t# 0 - doubly occupied, 1 - SpinUp,SpinDown" << endl;
 | 
|---|
| 1049 |                 // write out both types for easier changing afterwards
 | 
|---|
| 1050 |         //      switch (PsiType) {
 | 
|---|
| 1051 |         //              case 0:
 | 
|---|
| 1052 |                                 *output << "MaxPsiDouble\t" << config::MaxPsiDouble << "\t# here: specifying both maximum number of SpinUp- and -Down-states" << endl;
 | 
|---|
| 1053 |         //                      break;
 | 
|---|
| 1054 |         //              case 1:
 | 
|---|
| 1055 |                                 *output << "PsiMaxNoUp\t" << config::PsiMaxNoUp << "\t# here: specifying maximum number of SpinUp-states" << endl;
 | 
|---|
| 1056 |                                 *output << "PsiMaxNoDown\t" << config::PsiMaxNoDown << "\t# here: specifying maximum number of SpinDown-states" << endl;
 | 
|---|
| 1057 |         //                      break;
 | 
|---|
| 1058 |         //      }
 | 
|---|
| 1059 |                 *output << "AddPsis\t\t" << config::AddPsis << "\t# Additional unoccupied Psis for bandgap determination" << endl;
 | 
|---|
| 1060 |                 *output << endl;
 | 
|---|
| 1061 |                 *output << "RCut\t\t" << config::RCut << "\t# R-cut for the ewald summation" << endl;
 | 
|---|
| 1062 |                 *output << "StructOpt\t" << config::StructOpt << "\t# Do structure optimization beforehand" << endl;
 | 
|---|
| 1063 |                 *output << "IsAngstroem\t" << config::IsAngstroem << "\t# 0 - Bohr, 1 - Angstroem" << endl;
 | 
|---|
| 1064 |                 *output << "RelativeCoord\t" << config::RelativeCoord << "\t# whether ion coordinates are relative (1) or absolute (0)" << endl;
 | 
|---|
| 1065 |                 *output << "MaxTypes\t" << mol->ElementCount << "\t# maximum number of different ion types" << endl;
 | 
|---|
| 1066 |                 *output << endl;
 | 
|---|
| 1067 |                 result = result && mol->Checkout(output);
 | 
|---|
| 1068 |                 if (mol->MDSteps <=1 )
 | 
|---|
| 1069 |                         result = result && mol->Output(output);
 | 
|---|
| 1070 |                 else
 | 
|---|
| 1071 |                         result = result && mol->OutputTrajectories(output);
 | 
|---|
| 1072 |                 output->close();
 | 
|---|
| 1073 |                 output->clear();
 | 
|---|
| 1074 |                 delete(output);
 | 
|---|
| 1075 |                 return result;
 | 
|---|
| 1076 |         } else
 | 
|---|
| 1077 |                 return false;
 | 
|---|
| 1078 | };
 | 
|---|
| 1079 | 
 | 
|---|
| 1080 | /** Stores all elements in a MPQC input file.
 | 
|---|
| 1081 |  * Note that this format cannot be parsed again.
 | 
|---|
| 1082 |  * \param *filename name of file (without ".in" suffix!)
 | 
|---|
| 1083 |  * \param *mol pointer to molecule containing all atoms of the molecule
 | 
|---|
| 1084 |  */
 | 
|---|
| 1085 | bool config::SaveMPQC(const char *filename, molecule *mol) const
 | 
|---|
| 1086 | {
 | 
|---|
| 1087 |         int ElementNo = 0;
 | 
|---|
| 1088 |         int AtomNo;
 | 
|---|
| 1089 |         atom *Walker = NULL;
 | 
|---|
| 1090 |         element *runner = NULL;
 | 
|---|
| 1091 |         Vector *center = NULL;
 | 
|---|
| 1092 |         ofstream *output = NULL;
 | 
|---|
| 1093 |         stringstream *fname = NULL;
 | 
|---|
| 1094 | 
 | 
|---|
| 1095 |         // first without hessian
 | 
|---|
| 1096 |         fname = new stringstream;
 | 
|---|
| 1097 |         *fname << filename << ".in";
 | 
|---|
| 1098 |         output = new ofstream(fname->str().c_str(), ios::out);
 | 
|---|
| 1099 |         *output << "% Created by MoleCuilder" << endl;
 | 
|---|
| 1100 |         *output << "mpqc: (" << endl;
 | 
|---|
| 1101 |         *output << "\tsavestate = no" << endl;
 | 
|---|
| 1102 |         *output << "\tdo_gradient = yes" << endl;
 | 
|---|
| 1103 |         *output << "\tmole<MBPT2>: (" << endl;
 | 
|---|
| 1104 |         *output << "\t\tmaxiter = 200" << endl;
 | 
|---|
| 1105 |         *output << "\t\tbasis = $:basis" << endl;
 | 
|---|
| 1106 |         *output << "\t\tmolecule = $:molecule" << endl;
 | 
|---|
| 1107 |         *output << "\t\treference<CLHF>: (" << endl;
 | 
|---|
| 1108 |         *output << "\t\t\tbasis = $:basis" << endl;
 | 
|---|
| 1109 |         *output << "\t\t\tmolecule = $:molecule" << endl;
 | 
|---|
| 1110 |         *output << "\t\t)" << endl;
 | 
|---|
| 1111 |         *output << "\t)" << endl;
 | 
|---|
| 1112 |         *output << ")" << endl;
 | 
|---|
| 1113 |         *output << "molecule<Molecule>: (" << endl;
 | 
|---|
| 1114 |         *output << "\tunit = " << (IsAngstroem ? "angstrom" : "bohr" ) << endl;
 | 
|---|
| 1115 |         *output << "\t{ atoms geometry } = {" << endl;
 | 
|---|
| 1116 |         center = mol->DetermineCenterOfAll(output);
 | 
|---|
| 1117 |         // output of atoms
 | 
|---|
| 1118 |         runner = mol->elemente->start;
 | 
|---|
| 1119 |         while (runner->next != mol->elemente->end) { // go through every element
 | 
|---|
| 1120 |                 runner = runner->next;
 | 
|---|
| 1121 |                 if (mol->ElementsInMolecule[runner->Z]) { // if this element got atoms
 | 
|---|
| 1122 |                         ElementNo++;
 | 
|---|
| 1123 |                         AtomNo = 0;
 | 
|---|
| 1124 |                         Walker = mol->start;
 | 
|---|
| 1125 |                         while (Walker->next != mol->end) { // go through every atom of this element
 | 
|---|
| 1126 |                                 Walker = Walker->next;
 | 
|---|
| 1127 |                                 if (Walker->type == runner) { // if this atom fits to element
 | 
|---|
| 1128 |                                         AtomNo++;
 | 
|---|
| 1129 |                                         *output << "\t\t" << Walker->type->symbol << " [ " << Walker->x.x[0]-center->x[0] << "\t" << Walker->x.x[1]-center->x[1] << "\t" << Walker->x.x[2]-center->x[2] << " ]" << endl;
 | 
|---|
| 1130 |                                 }
 | 
|---|
| 1131 |                         }
 | 
|---|
| 1132 |                 }
 | 
|---|
| 1133 |         }
 | 
|---|
| 1134 |         delete(center);
 | 
|---|
| 1135 |         *output << "\t}" << endl;
 | 
|---|
| 1136 |         *output << ")" << endl;
 | 
|---|
| 1137 |         *output << "basis<GaussianBasisSet>: (" << endl;
 | 
|---|
| 1138 |         *output << "\tname = \""<< basis << "\"" << endl;
 | 
|---|
| 1139 |         *output << "\tmolecule = $:molecule" << endl;
 | 
|---|
| 1140 |         *output << ")" << endl;
 | 
|---|
| 1141 |         output->close();
 | 
|---|
| 1142 |         delete(output);
 | 
|---|
| 1143 |         delete(fname);
 | 
|---|
| 1144 | 
 | 
|---|
| 1145 |         // second with hessian
 | 
|---|
| 1146 |         fname = new stringstream;
 | 
|---|
| 1147 |         *fname << filename << ".hess.in";
 | 
|---|
| 1148 |         output = new ofstream(fname->str().c_str(), ios::out);
 | 
|---|
| 1149 |         *output << "% Created by MoleCuilder" << endl;
 | 
|---|
| 1150 |         *output << "mpqc: (" << endl;
 | 
|---|
| 1151 |         *output << "\tsavestate = no" << endl;
 | 
|---|
| 1152 |         *output << "\tdo_gradient = yes" << endl;
 | 
|---|
| 1153 |         *output << "\tmole<CLHF>: (" << endl;
 | 
|---|
| 1154 |         *output << "\t\tmaxiter = 200" << endl;
 | 
|---|
| 1155 |         *output << "\t\tbasis = $:basis" << endl;
 | 
|---|
| 1156 |         *output << "\t\tmolecule = $:molecule" << endl;
 | 
|---|
| 1157 |         *output << "\t)" << endl;
 | 
|---|
| 1158 |         *output << "\tfreq<MolecularFrequencies>: (" << endl;
 | 
|---|
| 1159 |         *output << "\t\tmolecule=$:molecule" << endl;
 | 
|---|
| 1160 |         *output << "\t)" << endl;
 | 
|---|
| 1161 |         *output << ")" << endl;
 | 
|---|
| 1162 |         *output << "molecule<Molecule>: (" << endl;
 | 
|---|
| 1163 |         *output << "\tunit = " << (IsAngstroem ? "angstrom" : "bohr" ) << endl;
 | 
|---|
| 1164 |         *output << "\t{ atoms geometry } = {" << endl;
 | 
|---|
| 1165 |         center = mol->DetermineCenterOfAll(output);
 | 
|---|
| 1166 |         // output of atoms
 | 
|---|
| 1167 |         runner = mol->elemente->start;
 | 
|---|
| 1168 |         while (runner->next != mol->elemente->end) { // go through every element
 | 
|---|
| 1169 |                 runner = runner->next;
 | 
|---|
| 1170 |                 if (mol->ElementsInMolecule[runner->Z]) { // if this element got atoms
 | 
|---|
| 1171 |                         ElementNo++;
 | 
|---|
| 1172 |                         AtomNo = 0;
 | 
|---|
| 1173 |                         Walker = mol->start;
 | 
|---|
| 1174 |                         while (Walker->next != mol->end) { // go through every atom of this element
 | 
|---|
| 1175 |                                 Walker = Walker->next;
 | 
|---|
| 1176 |                                 if (Walker->type == runner) { // if this atom fits to element
 | 
|---|
| 1177 |                                         AtomNo++;
 | 
|---|
| 1178 |                                         *output << "\t\t" << Walker->type->symbol << " [ " << Walker->x.x[0]-center->x[0] << "\t" << Walker->x.x[1]-center->x[1] << "\t" << Walker->x.x[2]-center->x[2] << " ]" << endl;
 | 
|---|
| 1179 |                                 }
 | 
|---|
| 1180 |                         }
 | 
|---|
| 1181 |                 }
 | 
|---|
| 1182 |         }
 | 
|---|
| 1183 |         delete(center);
 | 
|---|
| 1184 |         *output << "\t}" << endl;
 | 
|---|
| 1185 |         *output << ")" << endl;
 | 
|---|
| 1186 |         *output << "basis<GaussianBasisSet>: (" << endl;
 | 
|---|
| 1187 |         *output << "\tname = \"3-21G\"" << endl;
 | 
|---|
| 1188 |         *output << "\tmolecule = $:molecule" << endl;
 | 
|---|
| 1189 |         *output << ")" << endl;
 | 
|---|
| 1190 |         output->close();
 | 
|---|
| 1191 |         delete(output);
 | 
|---|
| 1192 |         delete(fname);
 | 
|---|
| 1193 | 
 | 
|---|
| 1194 |         return true;
 | 
|---|
| 1195 | };
 | 
|---|
| 1196 | 
 | 
|---|
| 1197 | /** Reads parameter from a parsed file.
 | 
|---|
| 1198 |  * The file is either parsed for a certain keyword or if null is given for
 | 
|---|
| 1199 |  * the value in row yth and column xth. If the keyword was necessity#critical,
 | 
|---|
| 1200 |  * then an error is thrown and the programme aborted.
 | 
|---|
| 1201 |  * \warning value is modified (both in contents and position)!
 | 
|---|
| 1202 |  * \param verbose 1 - print found value to stderr, 0 - don't
 | 
|---|
| 1203 |  * \param file file to be parsed
 | 
|---|
| 1204 |  * \param name Name of value in file (at least 3 chars!)
 | 
|---|
| 1205 |  * \param sequential 1 - do not reset file pointer to begin of file, 0 - set to beginning
 | 
|---|
| 1206 |  *                              (if file is sequentially parsed this can be way faster! However, beware of multiple values per line, as whole line is read -
 | 
|---|
| 1207 |  *                               best approach: 0 0 0 1 (not resetted only on last value of line) - and of yth, which is now
 | 
|---|
| 1208 |  *                               counted from this unresetted position!)
 | 
|---|
| 1209 |  * \param xth Which among a number of parameters it is (in grid case it's row number as grid is read as a whole!)
 | 
|---|
| 1210 |  * \param yth In grid case specifying column number, otherwise the yth \a name matching line
 | 
|---|
| 1211 |  * \param type Type of the Parameter to be read
 | 
|---|
| 1212 |  * \param value address of the value to be read (must have been allocated)
 | 
|---|
| 1213 |  * \param repetition determines, if the keyword appears multiply in the config file, which repetition shall be parsed, i.e. 1 if not multiply
 | 
|---|
| 1214 |  * \param critical necessity of this keyword being specified (optional, critical)
 | 
|---|
| 1215 |  * \return 1 - found, 0 - not found
 | 
|---|
| 1216 |  * \note Routine is taken from the pcp project and hack-a-slack adapted to C++
 | 
|---|
| 1217 |  */
 | 
|---|
| 1218 | int config::ParseForParameter(int verbose, ifstream *file, const char *name, int sequential, int const xth, int const yth, int type, void *value, int repetition, int critical) {
 | 
|---|
| 1219 |         int i,j;        // loop variables
 | 
|---|
| 1220 |         int length = 0, maxlength = -1;
 | 
|---|
| 1221 |         long file_position = file->tellg(); // mark current position
 | 
|---|
| 1222 |         char *dummy1, *dummy, *free_dummy;              // pointers in the line that is read in per step
 | 
|---|
| 1223 |         dummy1 = free_dummy = (char *) Malloc(256 * sizeof(char), "config::ParseForParameter: *free_dummy");
 | 
|---|
| 1224 | 
 | 
|---|
| 1225 |         //fprintf(stderr,"Parsing for %s\n",name);
 | 
|---|
| 1226 |         if (repetition == 0)
 | 
|---|
| 1227 |                 //Error(SomeError, "ParseForParameter(): argument repetition must not be 0!");
 | 
|---|
| 1228 |                 return 0;
 | 
|---|
| 1229 | 
 | 
|---|
| 1230 |         int line = 0; // marks line where parameter was found
 | 
|---|
| 1231 |         int found = (type >= grid) ? 0 : (-yth + 1);    // marks if yth parameter name was found
 | 
|---|
| 1232 |         while((found != repetition)) {
 | 
|---|
| 1233 |                 dummy1 = dummy = free_dummy;
 | 
|---|
| 1234 |                 do {
 | 
|---|
| 1235 |                         file->getline(dummy1, 256); // Read the whole line
 | 
|---|
| 1236 |                         if (file->eof()) {
 | 
|---|
| 1237 |                                 if ((critical) && (found == 0)) {
 | 
|---|
| 1238 |                                         Free((void **)&free_dummy, "config::ParseForParameter: *free_dummy");
 | 
|---|
| 1239 |                                         //Error(InitReading, name);
 | 
|---|
| 1240 |                                         fprintf(stderr,"Error:InitReading, critical %s not found\n", name);
 | 
|---|
| 1241 |                                         exit(255);
 | 
|---|
| 1242 |                                 } else {
 | 
|---|
| 1243 |                                         //if (!sequential)
 | 
|---|
| 1244 |                                         file->clear();
 | 
|---|
| 1245 |                                         file->seekg(file_position, ios::beg);   // rewind to start position
 | 
|---|
| 1246 |                                         Free((void **)&free_dummy, "config::ParseForParameter: *free_dummy");
 | 
|---|
| 1247 |                                         return 0;
 | 
|---|
| 1248 |                                 }
 | 
|---|
| 1249 |                         }
 | 
|---|
| 1250 |                         line++;
 | 
|---|
| 1251 |                 } while (dummy != NULL && dummy1 != NULL && ((dummy1[0] == '#') || (dummy1[0] == '\0'))); // skip commentary and empty lines
 | 
|---|
| 1252 | 
 | 
|---|
| 1253 |                 // C++ getline removes newline at end, thus re-add
 | 
|---|
| 1254 |                 if ((dummy1 != NULL) && (strchr(dummy1,'\n') == NULL)) {
 | 
|---|
| 1255 |                         i = strlen(dummy1);
 | 
|---|
| 1256 |                         dummy1[i] = '\n';
 | 
|---|
| 1257 |                         dummy1[i+1] = '\0';
 | 
|---|
| 1258 |                 }
 | 
|---|
| 1259 |                 //fprintf(stderr,"line %i ends at %i, newline at %i\n", line, strlen(dummy1), strchr(dummy1,'\n')-free_dummy);
 | 
|---|
| 1260 | 
 | 
|---|
| 1261 |                 if (dummy1 == NULL) {
 | 
|---|
| 1262 |                         if (verbose) fprintf(stderr,"Error reading line %i\n",line);
 | 
|---|
| 1263 |                 } else {
 | 
|---|
| 1264 |                         //fprintf(stderr,"Now parsing the line %i: %s\n", line, dummy1);
 | 
|---|
| 1265 |                 }
 | 
|---|
| 1266 |                 // Seek for possible end of keyword on line if given ...
 | 
|---|
| 1267 |                 if (name != NULL) {
 | 
|---|
| 1268 |                         dummy = strchr(dummy1,'\t');    // set dummy on first tab or space which ever's nearer
 | 
|---|
| 1269 |                         if (dummy == NULL) {
 | 
|---|
| 1270 |                                 dummy = strchr(dummy1, ' ');    // if not found seek for space
 | 
|---|
| 1271 |                                 while ((dummy != NULL) && ((*dummy == '\t') || (*dummy == ' ')))                // skip some more tabs and spaces if necessary
 | 
|---|
| 1272 |                                         dummy++;
 | 
|---|
| 1273 |                         }
 | 
|---|
| 1274 |                         if (dummy == NULL) {
 | 
|---|
| 1275 |                                 dummy = strchr(dummy1, '\n'); // set on line end then (whole line = keyword)
 | 
|---|
| 1276 |                                 //fprintf(stderr,"Error: Cannot find tabs or spaces on line %i in search for %s\n", line, name);
 | 
|---|
| 1277 |                                 //Free((void **)&free_dummy);
 | 
|---|
| 1278 |                                 //Error(FileOpenParams, NULL);
 | 
|---|
| 1279 |                         } else {
 | 
|---|
| 1280 |                                 //fprintf(stderr,"found tab at %i\n",(char *)dummy-(char *)dummy1);
 | 
|---|
| 1281 |                         }
 | 
|---|
| 1282 |                 } else dummy = dummy1;
 | 
|---|
| 1283 |                 // ... and check if it is the keyword!
 | 
|---|
| 1284 |                 //fprintf(stderr,"name %p, dummy %i/%c, dummy1 %i/%c, strlen(name) %i\n", &name, dummy, *dummy, dummy1, *dummy1, strlen(name));
 | 
|---|
| 1285 |                 if ((name == NULL) || (((dummy-dummy1 >= 3) && (strncmp(dummy1, name, strlen(name)) == 0)) && ((unsigned int)(dummy-dummy1) == strlen(name)))) {
 | 
|---|
| 1286 |                         found++; // found the parameter!
 | 
|---|
| 1287 |                         //fprintf(stderr,"found %s at line %i between %i and %i\n", name, line, dummy1, dummy);
 | 
|---|
| 1288 | 
 | 
|---|
| 1289 |                         if (found == repetition) {
 | 
|---|
| 1290 |                                 for (i=0;i<xth;i++) { // i = rows
 | 
|---|
| 1291 |                                         if (type >= grid) {
 | 
|---|
| 1292 |                                                 // grid structure means that grid starts on the next line, not right after keyword
 | 
|---|
| 1293 |                                                 dummy1 = dummy = free_dummy;
 | 
|---|
| 1294 |                                                 do {
 | 
|---|
| 1295 |                                                         file->getline(dummy1, 256); // Read the whole line, skip commentary and empty ones
 | 
|---|
| 1296 |                                                         if (file->eof()) {
 | 
|---|
| 1297 |                                                                 if ((critical) && (found == 0)) {
 | 
|---|
| 1298 |                                                                         Free((void **)&free_dummy, "config::ParseForParameter: *free_dummy");
 | 
|---|
| 1299 |                                                                         //Error(InitReading, name);
 | 
|---|
| 1300 |                                                                         fprintf(stderr,"Error:InitReading, critical %s not found\n", name);
 | 
|---|
| 1301 |                                                                         exit(255);
 | 
|---|
| 1302 |                                                                 } else {
 | 
|---|
| 1303 |                                                                         //if (!sequential)
 | 
|---|
| 1304 |                                                                         file->clear();
 | 
|---|
| 1305 |                                                                         file->seekg(file_position, ios::beg);   // rewind to start position
 | 
|---|
| 1306 |                                                                         Free((void **)&free_dummy, "config::ParseForParameter: *free_dummy");
 | 
|---|
| 1307 |                                                                         return 0;
 | 
|---|
| 1308 |                                                                 }
 | 
|---|
| 1309 |                                                         }
 | 
|---|
| 1310 |                                                         line++;
 | 
|---|
| 1311 |                                                 } while ((dummy1[0] == '#') || (dummy1[0] == '\n'));
 | 
|---|
| 1312 |                                                 if (dummy1 == NULL){
 | 
|---|
| 1313 |                                                         if (verbose) fprintf(stderr,"Error reading line %i\n", line);
 | 
|---|
| 1314 |                                                 } else {
 | 
|---|
| 1315 |                                                         //fprintf(stderr,"Reading next line %i: %s\n", line, dummy1);
 | 
|---|
| 1316 |                                                 }
 | 
|---|
| 1317 |                                         } else { // simple int, strings or doubles start in the same line
 | 
|---|
| 1318 |                                                 while ((*dummy == '\t') || (*dummy == ' '))      // skip interjacent tabs and spaces
 | 
|---|
| 1319 |                                                         dummy++;
 | 
|---|
| 1320 |                                         }
 | 
|---|
| 1321 |                                         // C++ getline removes newline at end, thus re-add
 | 
|---|
| 1322 |                                         if ((dummy1 != NULL) && (strchr(dummy1,'\n') == NULL)) {
 | 
|---|
| 1323 |                                                 j = strlen(dummy1);
 | 
|---|
| 1324 |                                                 dummy1[j] = '\n';
 | 
|---|
| 1325 |                                                 dummy1[j+1] = '\0';
 | 
|---|
| 1326 |                                         }
 | 
|---|
| 1327 | 
 | 
|---|
| 1328 |                                         int start = (type >= grid) ? 0 : yth-1 ;
 | 
|---|
| 1329 |                                         for (j=start;j<yth;j++) { // j = columns
 | 
|---|
| 1330 |                                                 // check for lower triangular area and upper triangular area
 | 
|---|
| 1331 |                                                 if ( ((i > j) && (type == upper_trigrid)) || ((j > i) && (type == lower_trigrid))) {
 | 
|---|
| 1332 |                                                         *((double *)value) = 0.0;
 | 
|---|
| 1333 |                                                         fprintf(stderr,"%f\t",*((double *)value));
 | 
|---|
| 1334 |                                                         value = (void *)((long)value + sizeof(double));
 | 
|---|
| 1335 |                                                         //value += sizeof(double);
 | 
|---|
| 1336 |                                                 } else {
 | 
|---|
| 1337 |                                                         // otherwise we must skip all interjacent tabs and spaces and find next value
 | 
|---|
| 1338 |                                                         dummy1 = dummy;
 | 
|---|
| 1339 |                                                         dummy = strchr(dummy1, '\t'); // seek for tab or space
 | 
|---|
| 1340 |                                                         if (dummy == NULL)
 | 
|---|
| 1341 |                                                                 dummy = strchr(dummy1, ' ');    // if not found seek for space
 | 
|---|
| 1342 |                                                         if (dummy == NULL) { // if still zero returned ...
 | 
|---|
| 1343 |                                                                 dummy = strchr(dummy1, '\n');   // ... at line end then
 | 
|---|
| 1344 |                                                                 if ((j < yth-1) && (type < 4)) {        // check if xth value or not yet
 | 
|---|
| 1345 |                                                                         if (critical) {
 | 
|---|
| 1346 |                                                                                 if (verbose) fprintf(stderr,"Error: EoL at %i and still missing %i value(s) for parameter %s\n", line, yth-j, name);
 | 
|---|
| 1347 |                                                                                 Free((void **)&free_dummy, "config::ParseForParameter: *free_dummy");
 | 
|---|
| 1348 |                                                                                 //return 0;
 | 
|---|
| 1349 |                                                                                 exit(255);
 | 
|---|
| 1350 |                                                                                 //Error(FileOpenParams, NULL);
 | 
|---|
| 1351 |                                                                         } else {
 | 
|---|
| 1352 |                                                                                 //if (!sequential)
 | 
|---|
| 1353 |                                                                                 file->clear();
 | 
|---|
| 1354 |                                                                                 file->seekg(file_position, ios::beg);   // rewind to start position
 | 
|---|
| 1355 |                                                                                 Free((void **)&free_dummy, "config::ParseForParameter: *free_dummy");
 | 
|---|
| 1356 |                                                                                 return 0;
 | 
|---|
| 1357 |                                                                         }
 | 
|---|
| 1358 |                                                                 }
 | 
|---|
| 1359 |                                                         } else {
 | 
|---|
| 1360 |                                                                 //fprintf(stderr,"found tab at %i\n",(char *)dummy-(char *)free_dummy);
 | 
|---|
| 1361 |                                                         }
 | 
|---|
| 1362 |                                                         if (*dummy1 == '#') {
 | 
|---|
| 1363 |                                                                 // found comment, skipping rest of line
 | 
|---|
| 1364 |                                                                 //if (verbose) fprintf(stderr,"Error: '#' at %i and still missing %i value(s) for parameter %s\n", line, yth-j, name);
 | 
|---|
| 1365 |                                                                 if (!sequential) { // here we need it!
 | 
|---|
| 1366 |                                                                         file->seekg(file_position, ios::beg);   // rewind to start position
 | 
|---|
| 1367 |                                                                 }
 | 
|---|
| 1368 |                                                                 Free((void **)&free_dummy, "config::ParseForParameter: *free_dummy");
 | 
|---|
| 1369 |                                                                 return 0;
 | 
|---|
| 1370 |                                                         }
 | 
|---|
| 1371 |                                                         //fprintf(stderr,"value from %i to %i\n",(char *)dummy1-(char *)free_dummy,(char *)dummy-(char *)free_dummy);
 | 
|---|
| 1372 |                                                         switch(type) {
 | 
|---|
| 1373 |                                                                 case (row_int):
 | 
|---|
| 1374 |                                                                         *((int *)value) = atoi(dummy1);
 | 
|---|
| 1375 |                                                                         if ((verbose) && (i==0) && (j==0)) fprintf(stderr,"%s = ", name);
 | 
|---|
| 1376 |                                                                         if (verbose) fprintf(stderr,"%i\t",*((int *)value));
 | 
|---|
| 1377 |                                                                                 value = (void *)((long)value + sizeof(int));
 | 
|---|
| 1378 |                                                                                 //value += sizeof(int);
 | 
|---|
| 1379 |                                                                         break;
 | 
|---|
| 1380 |                                                                 case(row_double):
 | 
|---|
| 1381 |                                                                 case(grid):
 | 
|---|
| 1382 |                                                                 case(lower_trigrid):
 | 
|---|
| 1383 |                                                                 case(upper_trigrid):
 | 
|---|
| 1384 |                                                                         *((double *)value) = atof(dummy1);
 | 
|---|
| 1385 |                                                                         if ((verbose) && (i==0) && (j==0)) fprintf(stderr,"%s = ", name);
 | 
|---|
| 1386 |                                                                         if (verbose) fprintf(stderr,"%lg\t",*((double *)value));
 | 
|---|
| 1387 |                                                                         value = (void *)((long)value + sizeof(double));
 | 
|---|
| 1388 |                                                                         //value += sizeof(double);
 | 
|---|
| 1389 |                                                                         break;
 | 
|---|
| 1390 |                                                                 case(double_type):
 | 
|---|
| 1391 |                                                                         *((double *)value) = atof(dummy1);
 | 
|---|
| 1392 |                                                                         if ((verbose) && (i == xth-1)) fprintf(stderr,"%s = %lg\n", name, *((double *) value));
 | 
|---|
| 1393 |                                                                         //value += sizeof(double);
 | 
|---|
| 1394 |                                                                         break;
 | 
|---|
| 1395 |                                                                 case(int_type):
 | 
|---|
| 1396 |                                                                         *((int *)value) = atoi(dummy1);
 | 
|---|
| 1397 |                                                                         if ((verbose) && (i == xth-1)) fprintf(stderr,"%s = %i\n", name, *((int *) value));
 | 
|---|
| 1398 |                                                                         //value += sizeof(int);
 | 
|---|
| 1399 |                                                                         break;
 | 
|---|
| 1400 |                                                                 default:
 | 
|---|
| 1401 |                                                                 case(string_type):
 | 
|---|
| 1402 |                                                                         if (value != NULL) {
 | 
|---|
| 1403 |                                                                                 //if (maxlength == -1) maxlength = strlen((char *)value); // get maximum size of string array
 | 
|---|
| 1404 |                                                                                 maxlength = MAXSTRINGSIZE;
 | 
|---|
| 1405 |                                                                                 length = maxlength > (dummy-dummy1) ? (dummy-dummy1) : maxlength; // cap at maximum
 | 
|---|
| 1406 |                                                                                 strncpy((char *)value, dummy1, length); // copy as much
 | 
|---|
| 1407 |                                                                                 ((char *)value)[length] = '\0'; // and set end marker
 | 
|---|
| 1408 |                                                                                 if ((verbose) && (i == xth-1)) fprintf(stderr,"%s is '%s' (%i chars)\n",name,((char *) value), length);
 | 
|---|
| 1409 |                                                                                 //value += sizeof(char);
 | 
|---|
| 1410 |                                                                         } else {
 | 
|---|
| 1411 |                                                                         }
 | 
|---|
| 1412 |                                                                 break;
 | 
|---|
| 1413 |                                                         }
 | 
|---|
| 1414 |                                                 }
 | 
|---|
| 1415 |                                                 while (*dummy == '\t')
 | 
|---|
| 1416 |                                                         dummy++;
 | 
|---|
| 1417 |                                         }
 | 
|---|
| 1418 |                                 }
 | 
|---|
| 1419 |                         }
 | 
|---|
| 1420 |                 }
 | 
|---|
| 1421 |         }
 | 
|---|
| 1422 |         if ((type >= row_int) && (verbose)) fprintf(stderr,"\n");
 | 
|---|
| 1423 |         Free((void **)&free_dummy, "config::ParseForParameter: *free_dummy");
 | 
|---|
| 1424 |         if (!sequential) {
 | 
|---|
| 1425 |                 file->clear();
 | 
|---|
| 1426 |                 file->seekg(file_position, ios::beg);   // rewind to start position
 | 
|---|
| 1427 |         }
 | 
|---|
| 1428 |         //fprintf(stderr, "End of Parsing\n\n");
 | 
|---|
| 1429 | 
 | 
|---|
| 1430 |         return (found); // true if found, false if not
 | 
|---|
| 1431 | }
 | 
|---|