00001
00007
00008
00009
00010
00011
00012
00013
00020 #include "InfoStream.h"
00021 #include "ComputeNonbondedUtil.h"
00022 #include "ConfigList.h"
00023 #include "SimParameters.h"
00024 #include "ParseOptions.h"
00025 #include "structures.h"
00026 #include "Communicate.h"
00027 #include "MStream.h"
00028 #include "Output.h"
00029 #include <stdio.h>
00030 #include <time.h>
00031 #ifdef NAMD_FFTW
00032 #ifdef NAMD_FFTW_3
00033 #include <fftw3.h>
00034 #else
00035
00036 #define fftwf_malloc fftw_malloc
00037 #define fftwf_free fftw_free
00038 #ifdef NAMD_FFTW_NO_TYPE_PREFIX
00039 #include <fftw.h>
00040 #include <rfftw.h>
00041 #else
00042 #include <sfftw.h>
00043 #include <srfftw.h>
00044 #endif
00045 #endif
00046 #endif
00047 #if defined(WIN32) && !defined(__CYGWIN__)
00048 #include <direct.h>
00049 #define CHDIR _chdir
00050 #define MKDIR(X) mkdir(X)
00051 #define PATHSEP '\\'
00052 #define PATHSEPSTR "\\"
00053 #else
00054 #include <unistd.h>
00055 #define CHDIR chdir
00056 #define MKDIR(X) mkdir(X,0777)
00057 #define PATHSEP '/'
00058 #define PATHSEPSTR "/"
00059 #endif
00060 #include <fcntl.h>
00061 #include <sys/stat.h>
00062 #ifdef WIN32
00063 #include <io.h>
00064 #endif
00065 #include <fstream>
00066 using namespace std;
00067
00068 #ifdef WIN32
00069 extern "C" {
00070 double erfc(double);
00071 }
00072 #endif
00073
00074 #include "strlib.h"
00075
00076
00077 #include "Debug.h"
00078
00079 #define XXXBIGREAL 1.0e32
00080
00081
00082
00083
00084
00085
00086
00087
00088
00089
00090
00091
00092
00093 void SimParameters::initialize_config_data(ConfigList *config, char *&cwd)
00094
00095 {
00096
00097 ParseOptions opts;
00098
00099 config_parser(opts);
00100
00102 if (!opts.check_consistency())
00103 {
00104 NAMD_die("Internal error in configuration file parser");
00105 }
00106
00107
00108
00109 if (!opts.set(*config))
00110 {
00111 NAMD_die("ERROR(S) IN THE CONFIGURATION FILE");
00112 }
00113
00115
00116 check_config(opts,config,cwd);
00117
00118 print_config(opts,config,cwd);
00119
00120 }
00121
00122
00123
00124
00125
00126
00127
00128 int atobool(const char *s) {
00129 return ( (! strncasecmp(s,"yes",8)) ||
00130 (! strncasecmp(s,"on",8)) ||
00131 (! strncasecmp(s,"true",8)) );
00132 };
00133
00134 void SimParameters::scriptSet(const char *param, const char *value) {
00135
00136 if ( CkMyRank() ) return;
00137
00138 #define MAX_SCRIPT_PARAM_SIZE 128
00139 #define SCRIPT_PARSE_BOOL(NAME,VAR) { if ( ! strncasecmp(param,(NAME),MAX_SCRIPT_PARAM_SIZE) ) { (VAR) = atobool(value); return; } }
00140 #define SCRIPT_PARSE_INT(NAME,VAR) { if ( ! strncasecmp(param,(NAME),MAX_SCRIPT_PARAM_SIZE) ) { (VAR) = atoi(value); return; } }
00141 #define SCRIPT_PARSE_FLOAT(NAME,VAR) { if ( ! strncasecmp(param,(NAME),MAX_SCRIPT_PARAM_SIZE) ) { (VAR) = atof(value); return; } }
00142 #define SCRIPT_PARSE_MOD_FLOAT(NAME,VAR,MOD) { if ( ! strncasecmp(param,(NAME),MAX_SCRIPT_PARAM_SIZE) ) { (VAR) = atof(value) MOD; return; } }
00143 #define SCRIPT_PARSE_VECTOR(NAME,VAR) { if ( ! strncasecmp(param,(NAME),MAX_SCRIPT_PARAM_SIZE) ) { (VAR).set(value); return; } }
00144 #define SCRIPT_PARSE_STRING(NAME,VAR) { if ( ! strncasecmp(param,(NAME),MAX_SCRIPT_PARAM_SIZE) ) { strcpy(VAR,value); return; } }
00145
00146 SCRIPT_PARSE_FLOAT("scriptArg1",scriptArg1)
00147 SCRIPT_PARSE_FLOAT("scriptArg2",scriptArg2)
00148 SCRIPT_PARSE_FLOAT("scriptArg3",scriptArg3)
00149 SCRIPT_PARSE_FLOAT("scriptArg4",scriptArg4)
00150 SCRIPT_PARSE_FLOAT("scriptArg5",scriptArg5)
00151 SCRIPT_PARSE_INT("numsteps",N)
00152 SCRIPT_PARSE_INT("firsttimestep",firstTimestep)
00153 SCRIPT_PARSE_FLOAT("reassignTemp",reassignTemp)
00154 SCRIPT_PARSE_FLOAT("rescaleTemp",rescaleTemp)
00155
00156 SCRIPT_PARSE_FLOAT("langevinTemp",langevinTemp)
00157 SCRIPT_PARSE_FLOAT("loweAndersenTemp",loweAndersenTemp)
00158 SCRIPT_PARSE_FLOAT("initialTemp",initialTemp)
00159 SCRIPT_PARSE_BOOL("useGroupPressure",useGroupPressure)
00160 SCRIPT_PARSE_BOOL("useFlexibleCell",useFlexibleCell)
00161 SCRIPT_PARSE_BOOL("useConstantArea",useConstantArea)
00162 SCRIPT_PARSE_BOOL("fixCellDims",fixCellDims)
00163 SCRIPT_PARSE_BOOL("fixCellDimX",fixCellDimX)
00164 SCRIPT_PARSE_BOOL("fixCellDimY",fixCellDimY)
00165 SCRIPT_PARSE_BOOL("fixCellDimZ",fixCellDimZ)
00166 SCRIPT_PARSE_BOOL("useConstantRatio",useConstantRatio)
00167 SCRIPT_PARSE_BOOL("LangevinPiston",langevinPistonOn)
00168 SCRIPT_PARSE_MOD_FLOAT("LangevinPistonTarget",
00169 langevinPistonTarget,/PRESSUREFACTOR)
00170 SCRIPT_PARSE_FLOAT("LangevinPistonPeriod",langevinPistonPeriod)
00171 SCRIPT_PARSE_FLOAT("LangevinPistonDecay",langevinPistonDecay)
00172 SCRIPT_PARSE_FLOAT("LangevinPistonTemp",langevinPistonTemp)
00173 SCRIPT_PARSE_MOD_FLOAT("SurfaceTensionTarget",
00174 surfaceTensionTarget,*(100.0/PRESSUREFACTOR))
00175 SCRIPT_PARSE_BOOL("BerendsenPressure",berendsenPressureOn)
00176 SCRIPT_PARSE_MOD_FLOAT("BerendsenPressureTarget",
00177 berendsenPressureTarget,/PRESSUREFACTOR)
00178 SCRIPT_PARSE_MOD_FLOAT("BerendsenPressureCompressibility",
00179 berendsenPressureCompressibility,*PRESSUREFACTOR)
00180 SCRIPT_PARSE_FLOAT("BerendsenPressureRelaxationTime",
00181 berendsenPressureRelaxationTime)
00182 SCRIPT_PARSE_FLOAT("constraintScaling",constraintScaling)
00183 SCRIPT_PARSE_FLOAT("consForceScaling",consForceScaling)
00184 SCRIPT_PARSE_STRING("outputname",outputFilename)
00185 SCRIPT_PARSE_STRING("tclBCArgs",tclBCArgs)
00186 SCRIPT_PARSE_VECTOR("eField",eField)
00187 SCRIPT_PARSE_FLOAT("eFieldFreq",eFieldFreq)
00188 SCRIPT_PARSE_FLOAT("eFieldPhase",eFieldPhase)
00189 SCRIPT_PARSE_FLOAT("accelMDE",accelMDE)
00190 SCRIPT_PARSE_FLOAT("accelMDalpha",accelMDalpha)
00191 SCRIPT_PARSE_FLOAT("accelMDTE",accelMDTE)
00192 SCRIPT_PARSE_FLOAT("accelMDTalpha",accelMDTalpha)
00193 SCRIPT_PARSE_VECTOR("stirAxis",stirAxis)
00194 SCRIPT_PARSE_VECTOR("stirPivot",stirPivot)
00195 if ( ! strncasecmp(param,"mgridforcescale",MAX_SCRIPT_PARAM_SIZE) ) {
00196 NAMD_die("Can't yet modify mgridforcescale in a script");
00197 return;
00198 }
00199 if ( ! strncasecmp(param,"mgridforcevoff",MAX_SCRIPT_PARAM_SIZE) ) {
00200 NAMD_die("Can't yet modify mgridforcevoff in a script");
00201 return;
00202 }
00203
00204 if ( ! strncasecmp(param,"fixedatoms",MAX_SCRIPT_PARAM_SIZE) ) {
00205 if ( ! fixedAtomsOn )
00206 NAMD_die("FixedAtoms may not be enabled in a script.");
00207 if ( ! fixedAtomsForces )
00208 NAMD_die("To use fixedAtoms in script first use fixedAtomsForces yes.");
00209 fixedAtomsOn = atobool(value);
00210 return;
00211 }
00212
00213
00214 SCRIPT_PARSE_INT("alchEquilSteps",alchEquilSteps)
00215
00216 if ( ! strncasecmp(param,"alchFepWCArcut1",MAX_SCRIPT_PARAM_SIZE) ) {
00217 alchFepWCArcut1 = atof(value);
00218 ComputeNonbondedUtil::select();
00219 return;
00220 }
00221
00222 if ( ! strncasecmp(param,"alchFepWCArcut2",MAX_SCRIPT_PARAM_SIZE) ) {
00223 alchFepWCArcut2 = atof(value);
00224 ComputeNonbondedUtil::select();
00225 return;
00226 }
00227
00228 if ( ! strncasecmp(param,"alchLambda",MAX_SCRIPT_PARAM_SIZE) ) {
00229 alchLambda = atof(value);
00230 ComputeNonbondedUtil::select();
00231 return;
00232 }
00233
00234 if ( ! strncasecmp(param,"alchLambda2",MAX_SCRIPT_PARAM_SIZE) ) {
00235 alchLambda2 = atof(value);
00236 ComputeNonbondedUtil::select();
00237 return;
00238 }
00239
00240
00241
00242
00243
00244
00245
00246
00247
00248
00249 if ( ! strncasecmp(param,"nonbondedScaling",MAX_SCRIPT_PARAM_SIZE) ) {
00250 nonbondedScaling = atof(value);
00251 ComputeNonbondedUtil::select();
00252 return;
00253 }
00254 if ( ! strncasecmp(param,"commOnly",MAX_SCRIPT_PARAM_SIZE) ) {
00255 commOnly = atobool(value);
00256 ComputeNonbondedUtil::select();
00257 return;
00258 }
00259
00260 char *error = new char[2 * MAX_SCRIPT_PARAM_SIZE + 100];
00261 sprintf(error,"Setting parameter %s from script failed!\n",param);
00262 NAMD_die(error);
00263
00264 }
00265
00266
00267
00268
00269
00270
00271
00272 void SimParameters::config_parser(ParseOptions &opts) {
00273
00274
00275
00276
00277
00278 HydrogenBonds = FALSE;
00279 useAntecedent = TRUE;
00280 aaAngleExp = 2;
00281 haAngleExp = 4;
00282 distAttExp = 4;
00283 distRepExp = 6;
00284 dhaCutoffAngle = 100.0;
00285 dhaOnAngle = 60.0;
00286 dhaOffAngle = 80.0;
00287 daCutoffDist = 7.5;
00288 daOnDist = 5.5;
00289 daOffDist = 6.5;
00290
00291 config_parser_basic(opts);
00292 config_parser_fileio(opts);
00293 config_parser_fullelect(opts);
00294 config_parser_methods(opts);
00295 config_parser_constraints(opts);
00296 #ifdef OPENATOM_VERSION
00297 config_parser_openatom(opts);
00298 #endif // OPENATOM_VERSION
00299
00300 config_parser_gridforce(opts);
00301 config_parser_mgridforce(opts);
00302 config_parser_movdrag(opts);
00303 config_parser_rotdrag(opts);
00304 config_parser_constorque(opts);
00305 config_parser_boundary(opts);
00306 config_parser_misc(opts);
00307
00308 }
00309
00310 void SimParameters::config_parser_basic(ParseOptions &opts) {
00311
00312
00313
00314
00315 opts.optional("main", "obsolete", "used to flag obsolete options",
00316 PARSE_STRING);
00317
00319 opts.require("main", "timestep", "size of the timestep, in fs",
00320 &dt, 1.0);
00321 opts.range("timestep", NOT_NEGATIVE);
00322 opts.units("timestep", N_FSEC);
00323
00324 opts.optional("main", "numsteps", "number of timesteps to perform",
00325 &N,0);
00326 opts.range("numsteps", NOT_NEGATIVE);
00327
00328 opts.optional("main", "stepspercycle",
00329 "Number of steps between atom migrations",
00330 &stepsPerCycle, 20);
00331 opts.range("stepspercycle", POSITIVE);
00332
00333 opts.require("main", "cutoff", "local electrostatic and Vdw distance",
00334 &cutoff);
00335 opts.range("cutoff", POSITIVE);
00336 opts.units("cutoff", N_ANGSTROM);
00337
00338 opts.optional("main", "nonbondedScaling", "nonbonded scaling factor",
00339 &nonbondedScaling, 1.0);
00340 opts.range("nonbondedScaling", NOT_NEGATIVE);
00341
00342 opts.optional("main", "limitDist", "limit nonbonded below this distance",
00343 &limitDist, 0.0);
00344 opts.range("limitDist", NOT_NEGATIVE);
00345
00346 opts.require("main", "exclude", "Electrostatic and VDW exclusion policy",
00347 PARSE_STRING);
00348
00349 opts.optional("exclude", "1-4scaling", "1-4 electrostatic scaling factor",
00350 &scale14, 1.0);
00351 opts.range("1-4scaling", POSITIVE);
00352
00353 opts.optionalB("main", "switching",
00354 "Should a smoothing function be used?", &switchingActive, TRUE);
00355
00356 opts.optionalB("switching", "vdwForceSwitching",
00357 "Use force switching for vdw?", &vdwForceSwitching, FALSE);
00358
00359 opts.optional("switching", "switchdist",
00360 "Distance for switching function activation",
00361 &switchingDist);
00362 opts.range("switchdist", POSITIVE);
00363 opts.units("switchdist", N_ANGSTROM);
00364
00365 opts.optionalB("main", "martiniSwitching",
00366 "Use Martini residue-based coarse-grain switching?", &martiniSwitching, FALSE);
00367 opts.optionalB("main", "martiniDielAllow",
00368 "Allow use of dielectric != 15.0 when using Martini", &martiniDielAllow, FALSE);
00369
00370 opts.optional("main", "pairlistdist", "Pairlist inclusion distance",
00371 &pairlistDist);
00372 opts.range("pairlistdist", POSITIVE);
00373 opts.units("pairlistdist", N_ANGSTROM);
00374
00375 opts.optional("main", "pairlistMinProcs", "Min procs for pairlists",
00376 &pairlistMinProcs,1);
00377 opts.range("pairlistMinProcs", POSITIVE);
00378
00379 opts.optional("main", "pairlistsPerCycle", "regenerate x times per cycle",
00380 &pairlistsPerCycle,2);
00381 opts.range("pairlistsPerCycle", POSITIVE);
00382
00383 opts.optional("main", "outputPairlists", "how often to print warnings",
00384 &outputPairlists, 0);
00385 opts.range("outputPairlists", NOT_NEGATIVE);
00386
00387 opts.optional("main", "pairlistShrink", "tol *= (1 - x) on regeneration",
00388 &pairlistShrink,0.01);
00389 opts.range("pairlistShrink", NOT_NEGATIVE);
00390
00391 opts.optional("main", "pairlistGrow", "tol *= (1 + x) on trigger",
00392 &pairlistGrow, 0.01);
00393 opts.range("pairlistGrow", NOT_NEGATIVE);
00394
00395 opts.optional("main", "pairlistTrigger", "trigger is atom > (1 - x) * tol",
00396 &pairlistTrigger, 0.3);
00397 opts.range("pairlistTrigger", NOT_NEGATIVE);
00398
00399 opts.optional("main", "temperature", "initial temperature",
00400 &initialTemp);
00401 opts.range("temperature", NOT_NEGATIVE);
00402 opts.units("temperature", N_KELVIN);
00403
00404 opts.optionalB("main", "COMmotion", "allow initial center of mass movement",
00405 &comMove, FALSE);
00406
00407 opts.optionalB("main", "zeroMomentum", "constrain center of mass",
00408 &zeroMomentum, FALSE);
00409 opts.optionalB("zeroMomentum", "zeroMomentumAlt", "constrain center of mass",
00410 &zeroMomentumAlt, FALSE);
00411
00412 opts.optionalB("main", "wrapWater", "wrap waters around periodic boundaries on output",
00413 &wrapWater, FALSE);
00414 opts.optionalB("main", "wrapAll", "wrap all clusters around periodic boundaries on output",
00415 &wrapAll, FALSE);
00416 opts.optionalB("main", "wrapNearest", "wrap to nearest image to cell origin",
00417 &wrapNearest, FALSE);
00418
00419 opts.optional("main", "dielectric", "dielectric constant",
00420 &dielectric, 1.0);
00421 opts.range("dielectric", POSITIVE);
00422
00423 opts.optional("main", "margin", "Patch width margin", &margin, XXXBIGREAL);
00424 opts.range("margin", NOT_NEGATIVE);
00425 opts.units("margin", N_ANGSTROM);
00426
00427 opts.optional("main", "seed", "Initial random number seed", &randomSeed);
00428 opts.range("seed", POSITIVE);
00429
00430 opts.optional("main", "outputEnergies", "How often to print energies in timesteps",
00431 &outputEnergies, 1);
00432 opts.range("outputEnergies", POSITIVE);
00433
00434 opts.optional("main", "outputMomenta", "How often to print linear and angular momenta in timesteps",
00435 &outputMomenta, 0);
00436 opts.range("outputMomenta", NOT_NEGATIVE);
00437
00438 opts.optional("main", "outputTiming", "How often to print timing data in timesteps",
00439 &outputTiming);
00440 opts.range("outputTiming", NOT_NEGATIVE);
00441
00442 opts.optional("main", "outputCudaTiming", "How often to print CUDA timing data in timesteps",
00443 &outputCudaTiming, 0);
00444 opts.range("outputCudaTiming", NOT_NEGATIVE);
00445
00446 opts.optional("main", "outputPressure", "How often to print pressure data in timesteps",
00447 &outputPressure, 0);
00448 opts.range("outputPressure", NOT_NEGATIVE);
00449
00450 opts.optionalB("main", "mergeCrossterms", "merge crossterm energy with dihedral when printing?",
00451 &mergeCrossterms, TRUE);
00452
00453 opts.optional("main", "MTSAlgorithm", "Multiple timestep algorithm",
00454 PARSE_STRING);
00455
00456 opts.optional("main", "longSplitting", "Long range force splitting option",
00457 PARSE_STRING);
00458
00459 opts.optional("main", "splitPatch", "Atom into patch splitting option",
00460 PARSE_STRING);
00461 opts.optional("main", "hgroupCutoff", "Hydrogen margin", &hgroupCutoff, 2.5);
00462
00463 opts.optional("main", "extendedSystem",
00464 "Initial configuration of extended system variables and periodic cell",
00465 PARSE_STRING);
00466
00467 opts.optional("main", "cellBasisVector1", "Basis vector for periodic cell",
00468 &cellBasisVector1);
00469 opts.optional("main", "cellBasisVector2", "Basis vector for periodic cell",
00470 &cellBasisVector2);
00471 opts.optional("main", "cellBasisVector3", "Basis vector for periodic cell",
00472 &cellBasisVector3);
00473 opts.optional("main", "cellOrigin", "Fixed center of periodic cell",
00474 &cellOrigin);
00475
00476 opts.optionalB("main", "molly", "Rigid bonds to hydrogen",&mollyOn,FALSE);
00477 opts.optional("main", "mollyTolerance", "Error tolerance for MOLLY",
00478 &mollyTol, 0.00001);
00479 opts.optional("main", "mollyIterations",
00480 "Max number of iterations for MOLLY", &mollyIter, 100);
00481
00482 opts.optional("main", "rigidBonds", "Rigid bonds to hydrogen",PARSE_STRING);
00483 opts.optional("main", "rigidTolerance",
00484 "Error tolerance for rigid bonds to hydrogen",
00485 &rigidTol, 1.0e-8);
00486 opts.optional("main", "rigidIterations",
00487 "Max number of SHAKE iterations for rigid bonds to hydrogen",
00488 &rigidIter, 100);
00489 opts.optionalB("main", "rigidDieOnError",
00490 "Die if rigidTolerance is not achieved after rigidIterations",
00491 &rigidDie, TRUE);
00492 opts.optionalB("main", "useSettle",
00493 "Use the SETTLE algorithm for rigid waters",
00494 &useSettle, TRUE);
00495
00496 opts.optional("main", "nonbondedFreq", "Nonbonded evaluation frequency",
00497 &nonbondedFrequency, 1);
00498 opts.range("nonbondedFreq", POSITIVE);
00499
00500 opts.optionalB("main", "outputPatchDetails", "print number of atoms in each patch",
00501 &outputPatchDetails, FALSE);
00502 opts.optionalB("main", "staticAtomAssignment", "never migrate atoms",
00503 &staticAtomAssignment, FALSE);
00504 opts.optionalB("main", "lonePairs", "Enable lone pairs", &lonepairs, FALSE);
00505 opts.optional("main", "waterModel", "Water model to use", PARSE_STRING);
00506 opts.optionalB("main", "LJcorrection", "Apply analytical tail corrections for energy and virial", &LJcorrection, FALSE);
00507 }
00508
00509 void SimParameters::config_parser_fileio(ParseOptions &opts) {
00510
00512
00513 opts.optional("main", "cwd", "current working directory", PARSE_STRING);
00514
00515
00516
00517
00518
00519
00520
00521 opts.optional("main", "coordinates", "initial PDB coordinate file",
00522 PARSE_STRING);
00523
00524 opts.optional("main", "velocities",
00525 "initial velocities, given as a PDB file", PARSE_STRING);
00526 opts.optional("main", "binvelocities",
00527 "initial velocities, given as a binary restart", PARSE_STRING);
00528 opts.optional("main", "bincoordinates",
00529 "initial coordinates in a binary restart file", PARSE_STRING);
00530
00531
00532
00533 opts.optional("main", "structure", "initial PSF structure file",
00534 PARSE_STRING);
00535
00536
00537
00538
00539 opts.optional("main", "parameters",
00540 "CHARMm 19 or CHARMm 22 compatable force field file (multiple "
00541 "inputs allowed)", PARSE_MULTIPLES);
00542
00543
00544
00546 opts.optionalB("parameters", "paraTypeXplor", "Parameter file in Xplor format?", ¶TypeXplorOn, FALSE);
00547 opts.optionalB("parameters", "paraTypeCharmm", "Parameter file in Charmm format?", ¶TypeCharmmOn, FALSE);
00548
00549
00550
00551 opts.optionalB("main", "GoForcesOn", "Go forces will be calculated", &goForcesOn, FALSE);
00552 opts.require("GoForcesOn", "GoParameters", "Go parameter file", goParameters);
00553 opts.require("GoForcesOn", "GoCoordinates", "target coordinates for Go forces", goCoordinates);
00554
00555
00556
00557
00558 opts.require("GoForcesOn", "GoMethod", "Which type of matrix should be used to store Go contacts?", PARSE_STRING);
00559
00560
00561 opts.require("main", "outputname",
00562 "prefix for the final PDB position and velocity filenames",
00563 outputFilename);
00564
00565 opts.optional("main", "auxFile", "Filename for data stream output",
00566 auxFilename);
00567
00568 opts.optional("main", "numinputprocs", "Number of pes to use for parallel input",
00569 &numinputprocs, 0);
00570 opts.range("numinputprocs", NOT_NEGATIVE);
00571
00572 opts.optional("main", "numoutputprocs", "Number of pes to use for parallel output",
00573 &numoutputprocs, 0);
00574 opts.range("numoutputprocs", NOT_NEGATIVE);
00575 opts.optional("main", "numoutputwriters", "Number of output processors that simultaneously write to an output file",
00576 &numoutputwrts, 1);
00577 opts.range("numoutputwriters", NOT_NEGATIVE);
00578
00579 opts.optional("main", "DCDfreq", "Frequency of DCD trajectory output, in "
00580 "timesteps", &dcdFrequency, 0);
00581 opts.range("DCDfreq", NOT_NEGATIVE);
00582 opts.optional("DCDfreq", "DCDfile", "DCD trajectory output file name",
00583 dcdFilename);
00584 opts.optionalB("DCDfreq", "DCDunitcell", "Store unit cell in dcd timesteps?",
00585 &dcdUnitCell);
00586
00587 opts.optional("main", "velDCDfreq", "Frequency of velocity "
00588 "DCD output, in timesteps", &velDcdFrequency, 0);
00589 opts.range("velDCDfreq", NOT_NEGATIVE);
00590 opts.optional("velDCDfreq", "velDCDfile", "velocity DCD output file name",
00591 velDcdFilename);
00592
00593 opts.optional("main", "forceDCDfreq", "Frequency of force"
00594 "DCD output, in timesteps", &forceDcdFrequency, 0);
00595 opts.range("forceDCDfreq", NOT_NEGATIVE);
00596 opts.optional("forceDCDfreq", "forceDCDfile", "force DCD output file name",
00597 forceDcdFilename);
00598
00599 opts.optional("main", "XSTfreq", "Frequency of XST trajectory output, in "
00600 "timesteps", &xstFrequency, 0);
00601 opts.range("XSTfreq", NOT_NEGATIVE);
00602 opts.optional("XSTfreq", "XSTfile", "Extended sytem trajectory output "
00603 "file name", xstFilename);
00604
00605 opts.optional("main", "restartfreq", "Frequency of restart file "
00606 "generation", &restartFrequency, 0);
00607 opts.range("restartfreq", NOT_NEGATIVE);
00608 opts.optional("restartfreq", "restartname", "Prefix for the position and "
00609 "velocity PDB files used for restarting", restartFilename);
00610 opts.optionalB("restartfreq", "restartsave", "Save restart files with "
00611 "unique filenames rather than overwriting", &restartSave, FALSE);
00612
00613 opts.optionalB("restartfreq", "binaryrestart", "Specify use of binary restart files ",
00614 &binaryRestart, TRUE);
00615
00616 opts.optionalB("outputname", "binaryoutput", "Specify use of binary output files ",
00617 &binaryOutput, TRUE);
00618
00619 opts.optionalB("main", "amber", "Is it AMBER force field?",
00620 &amberOn, FALSE);
00621 opts.optionalB("amber", "readexclusions", "Read exclusions from parm file?",
00622 &readExclusions, TRUE);
00623 opts.require("amber", "scnb", "1-4 VDW interactions are divided by scnb",
00624 &vdwscale14, 2.0);
00625 opts.require("amber", "parmfile", "AMBER parm file", PARSE_STRING);
00626 opts.optional("amber", "ambercoor", "AMBER coordinate file", PARSE_STRING);
00627
00628
00629
00630 opts.optional("alch", "alchoutfreq", "Frequency of alchemical energy"
00631 "output in timesteps", &alchOutFreq, 5);
00632 opts.range("alchoutfreq", NOT_NEGATIVE);
00633 opts.optional("alchoutfreq", "alchoutfile", "Alchemical energy output"
00634 "filename", alchOutFile);
00635
00636
00637
00638
00639
00640
00641
00642
00643
00644
00645
00646
00647 opts.optionalB("main", "gromacs", "Use GROMACS-like force field?",
00648 &gromacsOn, FALSE);
00649 opts.require("gromacs", "grotopfile", "GROMACS topology file",
00650 PARSE_STRING);
00651 opts.optional("gromacs", "grocoorfile","GROMACS coordinate file",
00652 PARSE_STRING);
00653
00654
00655 opts.optionalB("main", "vdwGeometricSigma",
00656 "Use geometric mean to combine L-J sigmas, as for OPLS",
00657 &vdwGeometricSigma, FALSE);
00658
00659
00660 opts.optional("main", "computeMapFile", "Filename for computeMap",
00661 computeMapFilename);
00662 opts.optionalB("main", "storeComputeMap", "store computeMap?",
00663 &storeComputeMap, FALSE);
00664 opts.optionalB("main", "loadComputeMap", "load computeMap?",
00665 &loadComputeMap, FALSE);
00666 }
00667
00668
00669 void SimParameters::config_parser_fullelect(ParseOptions &opts) {
00670
00672 #ifdef DPMTA
00673 DebugM(1,"DPMTA setup start\n");
00674
00675 opts.optionalB("main", "FMA", "Should FMA be used?", &FMAOn, FALSE);
00676 opts.optional("FMA", "FMALevels", "Tree levels to use in FMA", &FMALevels,
00677 5);
00678 opts.range("FMALevels", POSITIVE);
00679 opts.optional("FMA", "FMAMp", "Number of FMA multipoles", &FMAMp, 8);
00680 opts.range("FMAMp", POSITIVE);
00681 opts.optionalB("FMA", "FMAFFT", "Use FFT enhancement in FMA?", &FMAFFTOn, TRUE);
00682 opts.optional("FMAFFT", "FMAFFTBlock", "FFT blocking factor",
00683 &FMAFFTBlock, 4);
00684 opts.range("FMAFFTBlock", POSITIVE);
00685 DebugM(1,"DPMTA setup end\n");
00686 #else
00687
00688 FMAOn = FALSE;
00689 FMALevels = 0;
00690 FMAMp = 0;
00691 FMAFFTOn = FALSE;
00692 FMAFFTBlock = 0;
00693 #endif
00694
00695 opts.optional("main", "fullElectFrequency",
00696 "Number of steps between full electrostatic executions",
00697 &fullElectFrequency);
00698 opts.range("fullElectFrequency", POSITIVE);
00699
00700
00701 opts.optional("main", "fmaFrequency",
00702 "Number of steps between full electrostatic executions",
00703 &fmaFrequency);
00704 opts.range("fmaFrequency", POSITIVE);
00705
00706 opts.optional("main", "fmaTheta",
00707 "FMA theta parameter value",
00708 &fmaTheta,0.715);
00709 opts.range("fmaTheta", POSITIVE);
00710
00711 opts.optionalB("main", "FullDirect", "Should direct calculations of full electrostatics be performed?",
00712 &fullDirectOn, FALSE);
00713
00714
00716
00717 opts.optionalB("main", "MSM",
00718 "Use multilevel summation method for electrostatics?",
00719 &MSMOn, FALSE);
00720 opts.optional("MSM", "MSMApprox", "MSM approximation",
00721 &MSMApprox, 0);
00722 opts.optional("MSM", "MSMSplit", "MSM splitting",
00723 &MSMSplit, 0);
00724 opts.optional("MSM", "MSMLevels", "MSM maximum number of levels",
00725 &MSMLevels, 0);
00726 opts.optional("MSM", "MSMGridSpacing", "MSM grid spacing (Angstroms)",
00727 &MSMGridSpacing, 2.5);
00728 opts.optional("MSM", "MSMPadding", "MSM padding (Angstroms)",
00729 &MSMPadding, 2.5);
00730 opts.optional("MSM", "MSMBlockSizeX",
00731 "MSM grid block size along X direction (for decomposing parallel work)",
00732 &MSMBlockSizeX, 4);
00733 opts.optional("MSM", "MSMBlockSizeY",
00734 "MSM grid block size along Y direction (for decomposing parallel work)",
00735 &MSMBlockSizeY, 4);
00736 opts.optional("MSM", "MSMBlockSizeZ",
00737 "MSM grid block size along Z direction (for decomposing parallel work)",
00738 &MSMBlockSizeZ, 4);
00739
00740 opts.optionalB("MSM", "MsmSerial",
00741 "Use MSM serial version for long-range calculation?",
00742 &MsmSerialOn, FALSE);
00743
00744
00746
00747 opts.optionalB("main", "PME", "Use particle mesh Ewald for electrostatics?",
00748 &PMEOn, FALSE);
00749 opts.optional("PME", "PMETolerance", "PME direct space tolerance",
00750 &PMETolerance, 1.e-6);
00751 opts.optional("PME", "PMEInterpOrder", "PME interpolation order",
00752 &PMEInterpOrder, 4);
00753 opts.optional("PME", "PMEGridSizeX", "PME grid in x dimension",
00754 &PMEGridSizeX, 0);
00755 opts.optional("PME", "PMEGridSizeY", "PME grid in y dimension",
00756 &PMEGridSizeY, 0);
00757 opts.optional("PME", "PMEGridSizeZ", "PME grid in z dimension",
00758 &PMEGridSizeZ, 0);
00759 opts.optional("PME", "PMEGridSpacing", "Maximum PME grid spacing (Angstroms)",
00760 &PMEGridSpacing, 0.);
00761 opts.range("PMEGridSpacing", NOT_NEGATIVE);
00762 opts.optional("PME", "PMEProcessors",
00763 "PME FFT and reciprocal sum processor count", &PMEProcessors, 0);
00764 opts.optional("PME", "PMEMinSlices",
00765 "minimum thickness of PME reciprocal sum slab", &PMEMinSlices, 2);
00766 opts.range("PMEMinSlices", NOT_NEGATIVE);
00767 opts.optional("PME", "PMEPencils",
00768 "PME FFT and reciprocal sum pencil grid size", &PMEPencils, -1);
00769 opts.optional("PME", "PMEMinPoints",
00770 "minimum points per PME reciprocal sum pencil", &PMEMinPoints, 10000);
00771 opts.range("PMEMinPoints", NOT_NEGATIVE);
00772 opts.optionalB("main", "PMEBarrier", "Use barrier in PME?",
00773 &PMEBarrier, FALSE);
00774
00775 opts.optionalB("PME", "useOptPME", "Use the new scalable PME optimization", &useOptPME, FALSE);
00776 opts.optionalB("PME", "useManyToMany", "Use the many-to-many PME optimization", &useManyToMany, FALSE);
00777 if (PMEOn && !useOptPME)
00778 useManyToMany = false;
00779
00780 #ifdef DPME
00781 opts.optionalB("PME", "useDPME", "Use old DPME code?", &useDPME, FALSE);
00782 #else
00783 useDPME = 0;
00784 #endif
00785 opts.optionalB("main", "FFTWPatient", "Use intensive plan creation to optimize FFTW?",
00786 #ifdef WIN32
00787 &FFTWPatient, TRUE);
00788 #else
00789 &FFTWPatient, FALSE);
00790 #endif
00791
00792 opts.optionalB("main", "FFTWEstimate", "Use estimates to optimize FFTW?",
00793 #ifdef WIN32
00794 &FFTWEstimate, TRUE);
00795 #else
00796 &FFTWEstimate, FALSE);
00797 #endif
00798 opts.optionalB("main", "FFTWUseWisdom", "Read/save wisdom file for FFTW?",
00799 #ifdef WIN32
00800 &FFTWUseWisdom, FALSE);
00801 #else
00802 &FFTWUseWisdom, TRUE);
00803 #endif
00804 opts.optional("FFTWUseWisdom", "FFTWWisdomFile", "File for FFTW wisdom",
00805 FFTWWisdomFile);
00806
00807 }
00808
00809 void SimParameters::config_parser_methods(ParseOptions &opts) {
00810
00812 opts.optionalB("main", "minimization", "Should minimization be performed?",
00813 &minimizeCGOn, FALSE);
00814 opts.optionalB("main", "minVerbose", "Print extra minimization diagnostics?",
00815 &minVerbose, FALSE);
00816 opts.optional("main", "minTinyStep", "very first minimization steps",
00817 &minTinyStep, 1.0e-6);
00818 opts.range("minTinyStep", POSITIVE);
00819 opts.optional("main", "minBabyStep", "initial minimization steps",
00820 &minBabyStep, 1.0e-2);
00821 opts.range("minBabyStep", POSITIVE);
00822 opts.optional("main", "minLineGoal", "line minimization gradient reduction",
00823 &minLineGoal, 1.0e-3);
00824 opts.range("minLineGoal", POSITIVE);
00825
00826 opts.optionalB("main", "velocityQuenching",
00827 "Should old-style minimization be performed?", &minimizeOn, FALSE);
00828
00829 opts.optional("main", "maximumMove", "Maximum atom movement per step", &maximumMove, 0.0);
00830 opts.range("maximumMove", NOT_NEGATIVE);
00831 opts.units("maximumMove", N_ANGSTROM);
00832
00833 opts.optionalB("main", "Langevin", "Should Langevin dynamics be performed?",
00834 &langevinOn, FALSE);
00835 opts.require("Langevin", "langevinTemp", "Temperature for heat bath in Langevin "
00836 "dynamics", &langevinTemp);
00837 opts.range("langevinTemp", NOT_NEGATIVE);
00838 opts.units("langevinTemp", N_KELVIN);
00839 opts.optional("Langevin", "langevinDamping", "Damping coefficient (1/ps)",
00840 &langevinDamping);
00841 opts.range("langevinDamping", POSITIVE);
00842 opts.optionalB("Langevin", "langevinHydrogen", "Should Langevin dynamics be applied to hydrogen atoms?",
00843 &langevinHydrogen);
00844 opts.optional("Langevin", "langevinFile", "PDB file with temperature "
00845 "coupling terms (B(i)) (default is the PDB input file)",
00846 PARSE_STRING);
00847 opts.optional("Langevin", "langevinCol", "Column in the langevinFile "
00848 "containing the temperature coupling term B(i);\n"
00849 "default is 'O'", PARSE_STRING);
00850
00851
00852 opts.optionalB("main", "LoweAndersen", "Should Lowe-Andersen dynamics be performed?",
00853 &loweAndersenOn, FALSE);
00854 opts.require("LoweAndersen", "loweAndersenTemp", "Temperature for heat bath in Lowe-Andersen "
00855 "dynamics", &loweAndersenTemp);
00856 opts.range("loweAndersenTemp", NOT_NEGATIVE);
00857 opts.units("loweAndersenTemp", N_KELVIN);
00858 opts.optional("LoweAndersen", "loweAndersenRate", "Collision rate (1/ps)",
00859 &loweAndersenRate, 50);
00860 opts.range("loweAndersenRate", POSITIVE);
00861 opts.optional("LoweAndersen", "loweAndersenCutoff", "Cutoff radius",
00862 &loweAndersenCutoff, 2.7);
00863 opts.range("loweAndersenCutoff", POSITIVE);
00864 opts.units("loweAndersenCutoff", N_ANGSTROM);
00865
00866
00867
00868
00869 opts.optionalB("main", "alch", "Is achemical simulation being performed?",
00870 &alchOn, FALSE);
00871 opts.optional("alch", "alchType", "Which alchemical method to use?",
00872 PARSE_STRING);
00873 opts.optionalB("alch", "alchFepWCARepuOn", "WCA decomposition repu interaction in use?",
00874 &alchFepWCARepuOn, FALSE);
00875 opts.optionalB("alch", "alchFepWCADispOn", "WCA decomposition disp interaction in use?",
00876 &alchFepWCADispOn, FALSE);
00877 opts.optionalB("alch", "alchEnsembleAvg", "Ensemble Average in use?",
00878 &alchEnsembleAvg, TRUE);
00879 opts.optional("alch", "alchFepWCArcut1", "WCA repulsion Coeff1 used for generating"
00880 "the altered alchemical vDW interactions", &alchFepWCArcut1, 0.0);
00881 opts.optional("alch", "alchFepWCArcut2", "WCA repulsion Coeff2 used for generating"
00882 "the altered alchemical vDW interactions", &alchFepWCArcut2, 1.0);
00883 opts.range("alchFepWCArcut1", NOT_NEGATIVE);
00884 opts.range("alchFepWCArcut2", NOT_NEGATIVE);
00885 opts.require("alch", "alchLambda", "Coupling parameter value",
00886 &alchLambda);
00887 opts.require("alch", "alchLambda2", "Coupling comparison value",
00888 &alchLambda2);
00889 opts.optional("alch", "alchFile", "PDB file with perturbation flags "
00890 "default is the input PDB file", PARSE_STRING);
00891 opts.optional("alch", "alchCol", "Column in the alchFile with the "
00892 "perturbation flag", PARSE_STRING);
00893 opts.optional("alch", "alchEquilSteps", "Equilibration steps, before "
00894 "data collection in the alchemical window", &alchEquilSteps, 0);
00895 opts.range("alchEquilSteps", NOT_NEGATIVE);
00896 opts.optional("alch", "alchVdwShiftCoeff", "Coeff used for generating"
00897 "the altered alchemical vDW interactions", &alchVdwShiftCoeff, 5.);
00898 opts.range("alchVdwShiftCoeff", NOT_NEGATIVE);
00899
00900 opts.optional("alch", "alchElecLambdaStart", "Lambda at which electrostatic"
00901 "scaling of exnihilated particles begins", &alchElecLambdaStart, 0.5);
00902 opts.range("alchElecLambdaStart", NOT_NEGATIVE);
00903
00904 opts.optional("alch", "alchVdwLambdaEnd", "Lambda at which vdW"
00905 "scaling of exnihilated particles begins", &alchVdwLambdaEnd, 1.0);
00906 opts.range("alchVdwLambdaEnd", NOT_NEGATIVE);
00907
00908
00909
00910
00911
00912
00913
00914
00915
00916
00917
00918
00919
00920
00921
00922
00923
00924
00925
00926
00927
00928
00929
00930
00931
00932
00933
00934
00935
00936
00937
00938 opts.optionalB("main", "alchDecouple", "Enable alchemical decoupling?",
00939 &alchDecouple, FALSE);
00940
00941 opts.optionalB("main", "les", "Is locally enhanced sampling enabled?",
00942 &lesOn, FALSE);
00943 opts.require("les", "lesFactor", "Local enhancement factor", &lesFactor);
00944 opts.optional("les", "lesFile", "PDB file with enhancement flags "
00945 "default is the input PDB file", PARSE_STRING);
00946 opts.optional("les", "lesCol", "Column in the lesFile with the "
00947 "enhancement flag", PARSE_STRING);
00948 opts.optionalB("les", "lesReduceTemp", "Reduce enhanced atom temperature?",
00949 &lesReduceTemp, FALSE);
00950 opts.optionalB("les", "lesReduceMass", "Reduce enhanced atom mass?",
00951 &lesReduceMass, FALSE);
00952
00953
00954 opts.optionalB("main", "drude", "Perform integration of Drude oscillators?",
00955 &drudeOn, FALSE);
00956 opts.require("drude", "drudeTemp", "Temperature for freezing "
00957 "Drude oscillators", &drudeTemp);
00958 opts.range("drudeTemp", NOT_NEGATIVE);
00959 opts.units("drudeTemp", N_KELVIN);
00960 opts.optional("drude", "drudeDamping", "Damping coefficient (1/ps) for "
00961 "Drude oscillators", &drudeDamping);
00962 opts.range("drudeDamping", POSITIVE);
00963 opts.optional("drude", "drudeBondLen", "Drude oscillator bond length "
00964 "beyond which to apply restraint", &drudeBondLen);
00965 opts.range("drudeBondLen", POSITIVE);
00966 opts.optional("drude", "drudeBondConst", "Drude oscillator restraining "
00967 "force constant", &drudeBondConst);
00968 opts.range("drudeBondConst", POSITIVE);
00969 opts.optional("drude", "drudeNbtholeCut", "Nonbonded Thole interactions "
00970 "interaction radius", &drudeNbtholeCut);
00971 opts.range("drudeNbtholeCut", POSITIVE);
00972
00973
00974 opts.optionalB("main", "pairInteraction",
00975 "Are pair interactions calculated?", &pairInteractionOn, FALSE);
00976 opts.optional("pairInteraction", "pairInteractionFile",
00977 "PDB files with interaction flags " "default is the input PDB file",
00978 PARSE_STRING);
00979 opts.optional("pairInteraction", "pairInteractionCol",
00980 "Column in the pairInteractionFile with the interaction flags",
00981 PARSE_STRING);
00982 opts.require("pairInteraction", "pairInteractionGroup1",
00983 "Flag for interaction group 1", &pairInteractionGroup1);
00984 opts.optional("pairInteraction", "pairInteractionGroup2",
00985 "Flag for interaction group 2", &pairInteractionGroup2, -1);
00986 opts.optionalB("pairInteraction", "pairInteractionSelf",
00987 "Compute only within-group interactions?", &pairInteractionSelf,
00988 FALSE);
00989
00990 opts.optionalB("main", "cosAngles", "Are some angles cosine-based?", &cosAngles, FALSE);
00991
00992
00993
00994 opts.optionalB("main", "globalTest", "Should global integration (for development) be used?",
00995 &globalOn, FALSE);
00996 opts.optionalB("main", "dihedral", "Should dihedral angle dynamics be performed?",
00997 &dihedralOn, FALSE);
00998 COLDOn = FALSE;
00999 opts.optionalB("dihedral", "COLD", "Should overdamped Langevin dynamics be performed?",
01000 &COLDOn, FALSE);
01001 opts.require("COLD", "COLDTemp", "Temperature for heat bath in COLD",
01002 &COLDTemp);
01003 opts.range("COLDTemp", NOT_NEGATIVE);
01004 opts.units("COLDTemp", N_KELVIN);
01005 opts.require("COLD", "COLDRate", "Damping rate for COLD",
01006 &COLDRate, 3000.0);
01007 opts.range("COLDRate", NOT_NEGATIVE);
01008
01009
01010 opts.optionalB("main", "tcouple",
01011 "Should temperature coupling be performed?",
01012 &tCoupleOn, FALSE);
01013 opts.require("tcouple", "tCoupleTemp",
01014 "Temperature for temperature coupling", &tCoupleTemp);
01015 opts.range("tCoupleTemp", NOT_NEGATIVE);
01016 opts.units("tCoupleTemp", N_KELVIN);
01017 opts.optional("tCouple", "tCoupleFile", "PDB file with temperature "
01018 "coupling terms (B(i)) (default is the PDB input file)",
01019 PARSE_STRING);
01020 opts.optional("tCouple", "tCoupleCol", "Column in the tCoupleFile "
01021 "containing the temperature coupling term B(i);\n"
01022 "default is 'O'", PARSE_STRING);
01023
01024 opts.optional("main", "rescaleFreq", "Number of steps between "
01025 "velocity rescaling", &rescaleFreq);
01026 opts.range("rescaleFreq", POSITIVE);
01027 opts.optional("main", "rescaleTemp", "Target temperature for velocity rescaling",
01028 &rescaleTemp);
01029 opts.range("rescaleTemp", NOT_NEGATIVE);
01030 opts.units("rescaleTemp", N_KELVIN);
01031
01032 opts.optional("main", "reassignFreq", "Number of steps between "
01033 "velocity reassignment", &reassignFreq);
01034 opts.range("reassignFreq", POSITIVE);
01035 opts.optional("main", "reassignTemp", "Target temperature for velocity reassignment",
01036 &reassignTemp);
01037 opts.range("reassignTemp", NOT_NEGATIVE);
01038 opts.units("reassignTemp", N_KELVIN);
01039 opts.optional("main", "reassignIncr", "Temperature increment for velocity reassignment",
01040 &reassignIncr);
01041 opts.units("reassignIncr", N_KELVIN);
01042 opts.optional("main", "reassignHold", "Final holding temperature for velocity reassignment",
01043 &reassignHold);
01044 opts.range("reassignHold", NOT_NEGATIVE);
01045 opts.units("reassignHold", N_KELVIN);
01046
01048 opts.optionalB("main", "useGroupPressure",
01049 "Use group rather than atomic quantities for pressure control?",
01050 &useGroupPressure, FALSE);
01051
01053 opts.optionalB("main", "useFlexibleCell",
01054 "Use anisotropic cell fluctuation for pressure control?",
01055 &useFlexibleCell, FALSE);
01056
01057
01058 opts.optionalB("main", "fixCellDims",
01059 "Fix some cell dimensions?",
01060 &fixCellDims, FALSE);
01061
01062 opts.optionalB("fixCellDims", "fixCellDimX",
01063 "Fix the X dimension?",
01064 &fixCellDimX, FALSE);
01065 opts.optionalB("fixCellDims", "fixCellDimY",
01066 "Fix the Y dimension?",
01067 &fixCellDimY, FALSE);
01068 opts.optionalB("fixCellDims", "fixCellDimZ",
01069 "Fix the Z dimension?",
01070 &fixCellDimZ, FALSE);
01071
01073 opts.optionalB("main", "useConstantRatio",
01074 "Use constant X-Y ratio for pressure control?",
01075 &useConstantRatio, FALSE);
01076
01078 opts.optionalB("main", "useConstantArea",
01079 "Use constant area for pressure control?",
01080 &useConstantArea, FALSE);
01081
01083 opts.optionalB("main", "excludeFromPressure",
01084 "Should some atoms be excluded from pressure rescaling?",
01085 &excludeFromPressure, FALSE);
01086 opts.optional("excludeFromPressure", "excludeFromPressureFile",
01087 "PDB file for atoms to be excluded from pressure",
01088 PARSE_STRING);
01089 opts.optional("excludeFromPressure", "excludeFromPressureCol",
01090 "Column in the excludeFromPressureFile"
01091 "containing the flags (nonzero means excluded);\n"
01092 "default is 'O'", PARSE_STRING);
01093
01095 opts.optionalB("main", "BerendsenPressure",
01096 "Should Berendsen pressure bath coupling be performed?",
01097 &berendsenPressureOn, FALSE);
01098 opts.require("BerendsenPressure", "BerendsenPressureTarget",
01099 "Target pressure for pressure coupling",
01100 &berendsenPressureTarget);
01101
01102 opts.require("BerendsenPressure", "BerendsenPressureCompressibility",
01103 "Isothermal compressibility for pressure coupling",
01104 &berendsenPressureCompressibility);
01105
01106 opts.require("BerendsenPressure", "BerendsenPressureRelaxationTime",
01107 "Relaxation time for pressure coupling",
01108 &berendsenPressureRelaxationTime);
01109 opts.range("BerendsenPressureRelaxationTime", POSITIVE);
01110 opts.units("BerendsenPressureRelaxationTime", N_FSEC);
01111 opts.optional("BerendsenPressure", "BerendsenPressureFreq",
01112 "Number of steps between volume rescaling",
01113 &berendsenPressureFreq, 1);
01114 opts.range("BerendsenPressureFreq", POSITIVE);
01115
01117 opts.optionalB("main", "LangevinPiston",
01118 "Should Langevin piston pressure control be used?",
01119 &langevinPistonOn, FALSE);
01120 opts.require("LangevinPiston", "LangevinPistonTarget",
01121 "Target pressure for pressure control",
01122 &langevinPistonTarget);
01123 opts.require("LangevinPiston", "LangevinPistonPeriod",
01124 "Oscillation period for pressure control",
01125 &langevinPistonPeriod);
01126 opts.range("LangevinPistonPeriod", POSITIVE);
01127 opts.units("LangevinPistonPeriod", N_FSEC);
01128 opts.require("LangevinPiston", "LangevinPistonDecay",
01129 "Decay time for pressure control",
01130 &langevinPistonDecay);
01131 opts.range("LangevinPistonDecay", POSITIVE);
01132 opts.units("LangevinPistonDecay", N_FSEC);
01133 opts.require("LangevinPiston", "LangevinPistonTemp",
01134 "Temperature for pressure control piston",
01135 &langevinPistonTemp);
01136 opts.range("LangevinPistonTemp", POSITIVE);
01137 opts.units("LangevinPistonTemp", N_KELVIN);
01138 opts.optional("LangevinPiston", "StrainRate",
01139 "Initial strain rate for pressure control (x y z)",
01140 &strainRate);
01141
01143 opts.optional("main", "SurfaceTensionTarget",
01144 "Surface tension in the x-y plane",
01145 &surfaceTensionTarget, 0);
01146
01148 opts.optionalB("main", "pressureprofile", "Compute pressure profile?",
01149 &pressureProfileOn, FALSE);
01150 opts.require("pressureprofile", "pressureprofileslabs",
01151 "Number of pressure profile slabs", &pressureProfileSlabs, 10);
01152 opts.optional("pressureprofile", "pressureprofilefreq",
01153 "How often to store profile data", &pressureProfileFreq, 1);
01154 opts.optional("pressureprofile", "pressureProfileAtomTypes",
01155 "Number of pressure profile atom types", &pressureProfileAtomTypes, 1);
01156 opts.range("pressureProfileAtomTypes", POSITIVE);
01157 opts.optional("pressureProfile", "pressureProfileAtomTypesFile",
01158 "PDB files with pressure profile atom types" "default is the input PDB file",
01159 PARSE_STRING);
01160 opts.optional("pressureProfile", "pressureProfileAtomTypesCol",
01161 "Column in the pressureProfileAtomTypesFile with the atom types ",
01162 PARSE_STRING);
01163 opts.optionalB("pressureProfile", "pressureProfileEwald",
01164 "Compute Ewald contribution to pressure profile",
01165 &pressureProfileEwaldOn, FALSE);
01166 opts.optional("pressureProfile", "pressureProfileEwaldX",
01167 "Ewald grid size X", &pressureProfileEwaldX, 10);
01168 opts.range("pressureProfileEwaldX", POSITIVE);
01169 opts.optional("pressureProfile", "pressureProfileEwaldY",
01170 "Ewald grid size Y", &pressureProfileEwaldY, 10);
01171 opts.range("pressureProfileEwaldY", POSITIVE);
01172 opts.optional("pressureProfile", "pressureProfileEwaldZ",
01173 "Ewald grid size Z", &pressureProfileEwaldZ, 10);
01174 opts.range("pressureProfileEwaldZ", POSITIVE);
01175
01177 opts.optionalB("main", "accelMD", "Perform acclerated MD?", &accelMDOn, FALSE);
01178 opts.optional("accelMD", "accelMDFirstStep", "First accelMD step", &accelMDFirstStep, 0);
01179 opts.range("accelMDFirstStep", NOT_NEGATIVE);
01180 opts.optional("accelMD", "accelMDLastStep", "Last accelMD step", &accelMDLastStep, 0);
01181 opts.range("accelMDLastStep", NOT_NEGATIVE);
01182 opts.optional("accelMD", "accelMDOutFreq", "Frequency of accelMD output", &accelMDOutFreq, 1);
01183 opts.range("accelMDOutFreq", POSITIVE);
01184 opts.optionalB("accelMD", "accelMDdihe", "Apply boost to dihedral potential", &accelMDdihe, TRUE);
01185 opts.optionalB("accelMD", "accelMDDebugOn", "Debugging accelMD", &accelMDDebugOn, FALSE);
01186 opts.require("accelMD", "accelMDE","E for AMD", &accelMDE);
01187 opts.units("accelMDE", N_KCAL);
01188 opts.require("accelMD", "accelMDalpha","alpha for AMD", &accelMDalpha);
01189 opts.units("accelMDalpha", N_KCAL);
01190 opts.range("accelMDalpha", POSITIVE);
01191 opts.optionalB("accelMD", "accelMDdual", "Apply dual boost", &accelMDdual, FALSE);
01192 opts.require("accelMDdual", "accelMDTE","E for total potential under accelMDdual mode", &accelMDTE);
01193 opts.units("accelMDTE", N_KCAL);
01194 opts.require("accelMDdual", "accelMDTalpha","alpha for total potential under accelMDdual mode", &accelMDTalpha);
01195 opts.units("accelMDTalpha", N_KCAL);
01196 opts.range("accelMDTalpha", POSITIVE);
01197
01198
01199 opts.optionalB("main", "adaptTempMD", "Perform adaptive temperature sampling", &adaptTempOn, FALSE);
01200 opts.optional("adaptTempMD", "adaptTempFirstStep", "First adaptTemp step", &adaptTempFirstStep, 0);
01201 opts.range("adaptTempFirstStep", NOT_NEGATIVE);
01202 opts.optional("adaptTempMD", "adaptTempLastStep", "Last adaptTemp step", &adaptTempLastStep, 0);
01203 opts.range("adaptTempLastStep", NOT_NEGATIVE);
01204 opts.optional("adaptTempMD", "adaptTempOutFreq", "Frequency of adaptTemp output", &adaptTempOutFreq, 10);
01205 opts.range("adaptTempOutFreq", POSITIVE);
01206 opts.optional("adaptTempMD", "adaptTempFreq", "Frequency of writing average energies to adaptTempOutFile", &adaptTempFreq, 10);
01207 opts.range("adaptTempFreq", POSITIVE);
01208 opts.optionalB("adaptTempMD", "adaptTempDebug", "Print debug output for adaptTemp", &adaptTempDebug, FALSE);
01209 opts.optional("adaptTempMD", "adaptTempTmin","Minimun temperature for adaptTemp", &adaptTempTmin);
01210 opts.units("adaptTempTmin", N_KELVIN);
01211 opts.range("adaptTempTmin", POSITIVE);
01212 opts.optional("adaptTempMD", "adaptTempTmax","Maximum temperature for adaptTemp", &adaptTempTmax);
01213 opts.units("adaptTempTmax", N_KELVIN);
01214 opts.range("adaptTempTmax", POSITIVE);
01215 opts.optional("adaptTempMD", "adaptTempBins","Number of bins to store average energies", &adaptTempBins,0);
01216 opts.range("adaptTempBins", NOT_NEGATIVE);
01217 opts.optional("adaptTempMD", "adaptTempDt", "Integration timestep for Temp. updates", &adaptTempDt, 0.0001);
01218 opts.units("adaptTempDt", N_FSEC);
01219 opts.range("adaptTempDt", NOT_NEGATIVE);
01220 opts.optional("adaptTempMD", "adaptTempAutoDt", "Average temperature update in percent of temperature range", &adaptTempAutoDt, 0.0);
01221 opts.range("adaptTempAutoDt", NOT_NEGATIVE);
01222 opts.optional("adaptTempMD", "adaptTempCgamma", "Adaptive bin averaging constant", &adaptTempCgamma, 0.1);
01223 opts.range("adaptTempCgamma", NOT_NEGATIVE);
01224 opts.optionalB("adaptTempMD","adaptTempLangevin","Send adaptTemp temperature to langevin thermostat",&adaptTempLangevin,TRUE);
01225 opts.optionalB("adaptTempMD","adaptTempRescaling","Send adaptTemp temperature to velocity rescaling thermostat", &adaptTempRescale,TRUE);
01226 opts.optional("adaptTempMD", "adaptTempInFile", "File containing restart information for adaptTemp", adaptTempInFile);
01227 opts.optional("adaptTempMD", "adaptTempRestartFile", "File for writing adaptTemp restart information", adaptTempRestartFile);
01228 opts.require("adaptTempRestartFile","adaptTempRestartFreq", "Frequency of writing restart file", &adaptTempRestartFreq,0);
01229 opts.range("adaptTempRestartFreq",NOT_NEGATIVE);
01230 opts.optionalB("adaptTempMD", "adaptTempRandom", "Randomly assign a temperature if we step out of range", &adaptTempRandom, FALSE);
01231 }
01232
01233 void SimParameters::config_parser_constraints(ParseOptions &opts) {
01234
01236 opts.optionalB("main", "fixedatoms", "Are there fixed atoms?",
01237 &fixedAtomsOn, FALSE);
01238 opts.optionalB("fixedatoms", "fixedAtomsForces",
01239 "Calculate forces between fixed atoms? (Required to unfix during run.)",
01240 &fixedAtomsForces, FALSE);
01241 opts.optional("fixedatoms", "fixedAtomsFile", "PDB file with flags for "
01242 "fixed atoms (default is the PDB input file)",
01243 PARSE_STRING);
01244 opts.optional("fixedatoms", "fixedAtomsCol", "Column in the fixedAtomsFile "
01245 "containing the flags (nonzero means fixed);\n"
01246 "default is 'O'", PARSE_STRING);
01247 opts.optional("fixedatoms", "fixedAtomListFile", "the text input file for fixed atoms "
01248 "used for parallel input IO", PARSE_STRING);
01249
01251 opts.optionalB("main", "constraints", "Are harmonic constraints active?",
01252 &constraintsOn, FALSE);
01253 opts.require("constraints", "consexp", "Exponent for harmonic potential",
01254 &constraintExp, 2);
01255 opts.range("consexp", POSITIVE);
01256 opts.require("constraints", "consref", "PDB file containing reference "
01257 "positions",
01258 PARSE_STRING);
01259 opts.require("constraints", "conskfile", "PDB file containing force "
01260 "constaints in one of the columns", PARSE_STRING);
01261 opts.require("constraints", "conskcol", "Column of conskfile to use "
01262 "for the force constants", PARSE_STRING);
01263 opts.require("constraints", "constraintScaling", "constraint scaling factor",
01264 &constraintScaling, 1.0);
01265 opts.range("constraintScaling", NOT_NEGATIVE);
01266
01267
01268
01269
01270
01272 opts.optionalB("constraints", "selectConstraints",
01273 "Restrain only selected Cartesian components of the coordinates?",
01274 &selectConstraintsOn, FALSE);
01275 opts.optionalB("selectConstraints", "selectConstrX",
01276 "Restrain X components of coordinates ", &constrXOn, FALSE);
01277 opts.optionalB("selectConstraints", "selectConstrY",
01278 "Restrain Y components of coordinates ", &constrYOn, FALSE);
01279 opts.optionalB("selectConstraints", "selectConstrZ",
01280 "Restrain Z components of coordinates ", &constrZOn, FALSE);
01281
01282
01283
01284
01285
01287 opts.optionalB("constraints", "movingConstraints",
01288 "Are some of the constraints moving?",
01289 &movingConstraintsOn, FALSE);
01290 opts.require("movingConstraints", "movingConsVel",
01291 "Velocity of the movement, A/timestep", &movingConsVel);
01292
01293
01294
01295 opts.optionalB("constraints", "rotConstraints",
01296 "Are the constraints rotating?",
01297 &rotConstraintsOn, FALSE);
01298 opts.require("rotConstraints", "rotConsAxis",
01299 "Axis of rotation", &rotConsAxis);
01300 opts.require("rotConstraints", "rotConsPivot",
01301 "Pivot point of rotation",
01302 &rotConsPivot);
01303 opts.require("rotConstraints", "rotConsVel",
01304 "Velocity of rotation, deg/timestep", &rotConsVel);
01305
01306
01307
01308
01309 opts.optionalB("main", "extForces", "External command forces?",
01310 &extForcesOn, FALSE);
01311 opts.require("extForces", "extForcesCommand",
01312 "External forces command", extForcesCommand);
01313 opts.require("extForces", "extCoordFilename",
01314 "External forces coordinate filename", extCoordFilename);
01315 opts.require("extForces", "extForceFilename",
01316 "External forces force filename", extForceFilename);
01317
01318
01319 opts.optionalB("main", "printBadContacts", "Print atoms with huge forces?",
01320 &printBadContacts, FALSE);
01321
01322
01323
01324 opts.optionalB("main", "GBIS", "Use GB implicit solvent?",
01325 &GBISOn, FALSE);
01326 opts.optionalB("main", "GBISSer", "Use GB implicit solvent?",
01327 &GBISserOn, FALSE);
01328
01329 opts.optional("GBIS", "solventDielectric",
01330 "Solvent Dielectric", &solvent_dielectric, 78.5);
01331 opts.optional("GBIS", "intrinsicRadiusOffset",
01332 "Coulomb Radius Offset", &coulomb_radius_offset, 0.09);
01333 opts.optional("GBIS", "ionConcentration",
01334 "Ion Concentration", &ion_concentration, 0.2);
01335 opts.optional("GBIS", "GBISDelta",
01336 "delta from GBOBC", &gbis_delta, 1.0);
01337 opts.optional("GBIS", "GBISBeta",
01338 "beta from GBOBC", &gbis_beta, 0.8);
01339 opts.optional("GBIS", "GBISGamma",
01340 "gamma from GBOBC", &gbis_gamma, 4.85);
01341 opts.optional("GBIS", "alphaCutoff",
01342 "cutoff for calculating effective born radius", &alpha_cutoff, 15);
01343 opts.optional("GBIS", "alphaMax",
01344 "maximum allowable born radius", &alpha_max, 30);
01345 opts.optional("GBIS", "fsMax",
01346 "maximum screened intrinsic radius", &fsMax, 1.728);
01347
01348 opts.optionalB("main", "SASA", "Use Linear Combination of Pairwise Overlaps (LCPO) for calculating SASA",
01349 &LCPOOn, FALSE);
01350 opts.optional("SASA", "surfaceTension",
01351 "Surfce Tension for SASA (kcal/mol/Ang^2)", &surface_tension, 0.005);
01352
01353
01354
01355
01356 opts.optionalB("main", "SMD",
01357 "Do we use SMD option?",
01358 &SMDOn, FALSE);
01359 opts.require("SMD", "SMDVel",
01360 "Velocity of the movement, A/timestep", &SMDVel);
01361 opts.range("SMDVel", NOT_NEGATIVE);
01362 opts.require("SMD", "SMDDir",
01363 "Direction of movement", &SMDDir);
01364 opts.require("SMD", "SMDk",
01365 "Elastic constant for SMD", &SMDk);
01366 opts.optional("SMD", "SMDk2",
01367 "Transverse elastic constant for SMD", &SMDk2, 0);
01368 opts.range("SMDk", NOT_NEGATIVE);
01369 opts.range("SMDk2", NOT_NEGATIVE);
01370 opts.require("SMD", "SMDFile",
01371 "File for SMD information",
01372 SMDFile);
01373 opts.optional("SMD", "SMDOutputFreq",
01374 "Frequency of output",
01375 &SMDOutputFreq, 1);
01376 opts.range("SMDOutputFreq", POSITIVE);
01377
01378
01379
01380
01381 opts.optionalB("main", "tabulatedEnergies", "Do we get energies from a table?", &tabulatedEnergies, FALSE);
01382
01383 opts.require("tabulatedEnergies", "tabulatedEnergiesFile", "File containing energy table", tabulatedEnergiesFile);
01384 opts.require("tabulatedEnergies", "tableInterpType", "Cubic or linear interpolation", tableInterpType);
01385
01386
01387 opts.optionalB("main", "TMD", "Perform Targeted MD?", &TMDOn, FALSE);
01388 opts.optional("TMD", "TMDk", "Elastic constant for TMD", &TMDk, 0);
01389 opts.range("TMDk", NOT_NEGATIVE);
01390 opts.require("TMD", "TMDFile", "File for TMD information", TMDFile);
01391 opts.optionalB("TMD", "TMDDiffRMSD", "Restrain Difference between the RMSD from two structures", &TMDDiffRMSD, FALSE);
01392 opts.require("TMDDiffRMSD", "TMDFile2", "Second file for TMD information", TMDFile2);
01393
01394 opts.optional("TMD", "TMDOutputFreq", "Frequency of TMD output",
01395 &TMDOutputFreq, 1);
01396 opts.range("TMDOutputFreq", POSITIVE);
01397 opts.require("TMD", "TMDLastStep", "Last TMD timestep", &TMDLastStep);
01398 opts.range("TMDLastStep", POSITIVE);
01399 opts.optional("TMD", "TMDFirstStep", "First TMD step (default 0)", &TMDFirstStep, 0);
01400 opts.optional("TMD", "TMDInitialRMSD", "Target RMSD at first TMD step (default -1 to use initial coordinates)", &TMDInitialRMSD);
01401 TMDInitialRMSD = -1;
01402 opts.optional("TMD", "TMDFinalRMSD", "Target RMSD at last TMD step (default 0 )", &TMDFinalRMSD, 0);
01403 opts.range("TMDInitialRMSD", NOT_NEGATIVE);
01404
01405
01406
01407 opts.optionalB("main", "symmetryRestraints", "Enable symmetry restraints?", &symmetryOn, FALSE);
01408 opts.optional("symmetryRestraints", "symmetryk", "Elastic constant for symmetry restraints", &symmetryk, 0);
01409 opts.range("symmetryk", NOT_NEGATIVE);
01410 opts.optional("symmetryRestraints", "symmetrykfile", "PDB file specifying force contants on a per-atom basis", PARSE_MULTIPLES);
01411 opts.optionalB("symmetryRestraints", "symmetryScaleForces", "Scale applied forces over time?", &symmetryScaleForces, FALSE);
01412 opts.require("symmetryRestraints", "symmetryFile", "File for symmetry information", PARSE_MULTIPLES);
01413 opts.optional("symmetryRestraints", "symmetryMatrixFile", "File(s) for transfromation matrices", PARSE_MULTIPLES);
01414 opts.optional("symmetryRestraints", "symmetryLastStep", "Last symmetry timestep", &symmetryLastStep, -1);
01415 opts.optional("symmetryRestraints", "symmetryFirstStep", "First symmetry step (default 0)", &symmetryFirstStep, 0);
01416 opts.optional("symmetryRestraints", "symmetryLastFullStep", "Last full force symmetry timestep (default symmetryLastStep)", &symmetryLastFullStep, symmetryLastStep);
01417 opts.optional("symmetryRestraints", "symmetryFirstFullStep", "First full force symmetry step (default symmetryFirstStep)", &symmetryFirstFullStep, symmetryFirstStep);
01418
01419
01421 opts.optionalB("main", "tclForces", "Are Tcl global forces active?",
01422 &tclForcesOn, FALSE);
01423 opts.require("tclForces", "tclForcesScript",
01424 "Tcl script for global forces", PARSE_MULTIPLES);
01425
01427 opts.optionalB("main", "tclBC", "Are Tcl boundary forces active?",
01428 &tclBCOn, FALSE);
01429 opts.require("tclBC", "tclBCScript",
01430 "Tcl script defining calcforces for boundary forces", PARSE_STRING);
01431 tclBCScript = 0;
01432 opts.optional("tclBC", "tclBCArgs", "Extra args for calcforces command",
01433 tclBCArgs);
01434 tclBCArgs[0] = 0;
01435
01437 opts.optionalB("main", "miscForces", "Are misc global forces active?",
01438 &miscForcesOn, FALSE);
01439 opts.optional("miscForces", "miscForcesScript",
01440 "script for misc forces", PARSE_MULTIPLES);
01441
01443 opts.optionalB("main", "freeEnergy", "Perform free energy perturbation?",
01444 &freeEnergyOn, FALSE);
01445 opts.require("freeEnergy", "freeEnergyConfig",
01446 "Configuration file for free energy perturbation", PARSE_MULTIPLES);
01447
01449 opts.optionalB("main", "constantforce", "Apply constant force?",
01450 &consForceOn, FALSE);
01451 opts.optional("constantforce", "consForceFile",
01452 "Configuration file for constant forces", PARSE_STRING);
01453 opts.require("constantforce", "consForceScaling",
01454 "Scaling factor for constant forces", &consForceScaling, 1.0);
01455
01457 opts.optionalB("main", "colvars", "Is the colvars module enabled?",
01458 &colvarsOn, FALSE);
01459 opts.require("colvars", "colvarsConfig",
01460 "configuration for the collective variables", PARSE_STRING);
01461 opts.optional("colvars", "colvarsInput",
01462 "input restart file for the collective variables", PARSE_STRING);
01463
01464 }
01465
01466 #ifdef OPENATOM_VERSION
01467 void SimParameters::config_parser_openatom(ParseOptions &opts) {
01468 opts.optionalB("main", "openatom", "OpenAtom active?", &openatomOn, FALSE);
01469 opts.require("openatom", "openatomDriverFile", "What config file specifies openatom input parameters", PARSE_STRING);
01470 opts.require("openatom", "openatomPhysicsFile", "What structure file specifies openatom input system", PARSE_STRING);
01471 opts.require("openatom", "openatomPdbFile", "NAMD input file defining QM and MM regions", PARSE_STRING);
01472 opts.optional("openatom", "openatomCol", "Column in the openatomPdb with the QM/MM flag", PARSE_STRING);
01473 }
01474 #endif // OPENATOM_VERSION
01475
01476
01477 void SimParameters::config_parser_mgridforce(ParseOptions &opts) {
01479 opts.optionalB("main", "mgridforce", "Is Multiple gridforce active?",
01480 &mgridforceOn, FALSE);
01481 opts.optional("mgridforce", "mgridforcevolts", "Is Gridforce using Volts/eV as units?",
01482 PARSE_MULTIPLES);
01483 opts.require("mgridforce", "mgridforcescale", "Scale factor by which to multiply "
01484 "grid forces", PARSE_MULTIPLES);
01485 opts.require("mgridforce", "mgridforcefile", "PDB file containing force "
01486 "multipliers in one of the columns", PARSE_MULTIPLES);
01487 opts.require("mgridforce", "mgridforcecol", "Column of gridforcefile to "
01488 "use for force multiplier", PARSE_MULTIPLES);
01489 opts.optional("mgridforce", "mgridforcechargecol", "Column of gridforcefile to "
01490 "use for charge", PARSE_MULTIPLES);
01491 opts.require("mgridforce", "mgridforcepotfile", "Gridforce potential file",
01492 PARSE_MULTIPLES);
01493 opts.optional("mgridforce", "mgridforcecont1", "Use continuous grid "
01494 "in K1 direction?", PARSE_MULTIPLES);
01495 opts.optional("mgridforce", "mgridforcecont2", "Use continuous grid "
01496 "in K2 direction?", PARSE_MULTIPLES);
01497 opts.optional("mgridforce", "mgridforcecont3", "Use continuous grid "
01498 "in K3 direction?", PARSE_MULTIPLES);
01499 opts.optional("mgridforce", "mgridforcevoff", "Gridforce potential offsets",
01500 PARSE_MULTIPLES);
01501 opts.optional("mgridforce", "mgridforcelite", "Use Gridforce Lite?",
01502 PARSE_MULTIPLES);
01503 }
01504
01505 void SimParameters::config_parser_gridforce(ParseOptions &opts) {
01507 opts.optionalB("main", "gridforce", "Is Gridforce active?",
01508 &gridforceOn, FALSE);
01509 opts.optionalB("gridforce", "gridforcevolts", "Is Gridforce using Volts/eV as units?",
01510 &gridforceVolts, FALSE);
01511 opts.require("gridforce", "gridforcescale", "Scale factor by which to multiply "
01512 "grid forces", &gridforceScale);
01513 opts.require("gridforce", "gridforcefile", "PDB file containing force "
01514 "multipliers in one of the columns", PARSE_STRING);
01515 opts.require("gridforce", "gridforcecol", "Column of gridforcefile to "
01516 "use for force multiplier", PARSE_STRING);
01517 opts.optional("gridforce", "gridforcechargecol", "Column of gridforcefile to "
01518 "use for charge", PARSE_STRING);
01519 opts.require("gridforce", "gridforcepotfile", "Gridforce potential file",
01520 PARSE_STRING);
01521 opts.optionalB("gridforce", "gridforcecont1", "Use continuous grid "
01522 "in A1 direction?", &gridforceContA1, FALSE);
01523 opts.optionalB("gridforce", "gridforcecont2", "Use continuous grid "
01524 "in A2 direction?", &gridforceContA2, FALSE);
01525 opts.optionalB("gridforce", "gridforcecont3", "Use continuous grid "
01526 "in A3 direction?", &gridforceContA3, FALSE);
01527 opts.optional("gridforce", "gridforcevoff", "Gridforce potential offsets",
01528 &gridforceVOffset);
01529 opts.optionalB("gridforce", "gridforcelite", "Use Gridforce Lite?",
01530 &gridforceLite, FALSE);
01531 }
01532
01533
01534 void SimParameters::config_parser_movdrag(ParseOptions &opts) {
01536 opts.optionalB("main", "movDragOn", "Do we apply moving drag?",
01537 &movDragOn, FALSE);
01538 opts.require("movDragOn", "movDragFile",
01539 "Main moving drag PDB file", movDragFile);
01540 opts.require("movDragOn", "movDragCol",
01541 "Main moving drag PDB column", PARSE_STRING);
01542 opts.require("movDragOn", "movDragGlobVel",
01543 "Global moving drag velocity (A/step)", &movDragGlobVel);
01544 opts.require("movDragOn", "movDragVelFile",
01545 "Moving drag linear velocity file", movDragVelFile);
01546 }
01547
01548 void SimParameters::config_parser_rotdrag(ParseOptions &opts) {
01550 opts.optionalB("main", "rotDragOn", "Do we apply rotating drag?",
01551 &rotDragOn, FALSE);
01552 opts.require("rotDragOn", "rotDragFile",
01553 "Main rotating drag PDB file", rotDragFile);
01554 opts.require("rotDragOn", "rotDragCol",
01555 "Main rotating drag PDB column", PARSE_STRING);
01556 opts.require("rotDragOn", "rotDragAxisFile",
01557 "Rotating drag axis file", rotDragAxisFile);
01558 opts.require("rotDragOn", "rotDragPivotFile",
01559 "Rotating drag pivot point file", rotDragPivotFile);
01560 opts.require("rotDragOn", "rotDragGlobVel",
01561 "Global rotating drag angular velocity (deg/step)", &rotDragGlobVel);
01562 opts.require("rotDragOn", "rotDragVelFile",
01563 "Rotating drag angular velocity file", rotDragVelFile);
01564 opts.require("rotDragOn", "rotDragVelCol",
01565 "Rotating drag angular velocity column", PARSE_STRING);
01566 }
01567
01568 void SimParameters::config_parser_constorque(ParseOptions &opts) {
01570 opts.optionalB("main", "consTorqueOn", "Do we apply \"constant\" torque?",
01571 &consTorqueOn, FALSE);
01572 opts.require("consTorqueOn", "consTorqueFile",
01573 "Main \"constant\" torque PDB file", consTorqueFile);
01574 opts.require("consTorqueOn", "consTorqueCol",
01575 "Main \"constant\" torque PDB column", PARSE_STRING);
01576 opts.require("consTorqueOn", "consTorqueAxisFile",
01577 "\"Constant\" torque axis file", consTorqueAxisFile);
01578 opts.require("consTorqueOn", "consTorquePivotFile",
01579 "\"Constant\" torque pivot point file", consTorquePivotFile);
01580 opts.require("consTorqueOn", "consTorqueGlobVal",
01581 "Global \"constant\" torque value (Kcal/(mol*A^2))", &consTorqueGlobVal);
01582 opts.require("consTorqueOn", "consTorqueValFile",
01583 "\"constant\" torque factors file", consTorqueValFile);
01584 opts.require("consTorqueOn", "consTorqueValCol",
01585 "\"constant\" torque factors column", PARSE_STRING);
01586 }
01587
01588 void SimParameters::config_parser_boundary(ParseOptions &opts) {
01589
01591 opts.optionalB("main", "sphericalBC", "Are spherical boundary counditions "
01592 "active?", &sphericalBCOn, FALSE);
01593 opts.require("sphericalBC", "sphericalBCCenter",
01594 "Center of spherical boundaries", &sphericalCenter);
01595 opts.require("sphericalBC", "sphericalBCr1", "Radius for first sphere "
01596 "potential", &sphericalBCr1);
01597 opts.range("sphericalBCr1", POSITIVE);
01598 opts.units("sphericalBCr1", N_ANGSTROM);
01599 opts.require("sphericalBC", "sphericalBCk1", "Force constant for first "
01600 "sphere potential (+ is an inward force, - outward)",
01601 &sphericalBCk1);
01602 opts.units("sphericalBCk1", N_KCAL);
01603 opts.optional("sphericalBC", "sphericalBCexp1", "Exponent for first "
01604 "sphere potential", &sphericalBCexp1, 2);
01605 opts.range("sphericalBCexp1", POSITIVE);
01606
01607 opts.optional("sphericalBCr1", "sphericalBCr2", "Radius for second sphere "
01608 "potential", &sphericalBCr2);
01609 opts.range("sphericalBCr2", POSITIVE);
01610 opts.units("sphericalBCr2", N_ANGSTROM);
01611 opts.require("sphericalBCr2", "sphericalBCk2", "Force constant for second "
01612 "sphere potential (+ is an inward force, - outward)",
01613 &sphericalBCk2);
01614 opts.units("sphericalBCk2", N_KCAL);
01615 opts.optional("sphericalBCr2", "sphericalBCexp2", "Exponent for second "
01616 "sphere potential", &sphericalBCexp2, 2);
01617 opts.range("sphericalBCexp2", POSITIVE);
01618
01620 opts.optionalB("main", "cylindricalBC", "Are cylindrical boundary counditions "
01621 "active?", &cylindricalBCOn, FALSE);
01622 opts.require("cylindricalBC", "cylindricalBCr1", "Radius for first cylinder "
01623 "potential", &cylindricalBCr1);
01624 opts.range("cylindricalBCr1", POSITIVE);
01625 opts.units("cylindricalBCr1", N_ANGSTROM);
01626 opts.require("cylindricalBC", "cylindricalBCk1", "Force constant for first "
01627 "cylinder potential (+ is an inward force, - outward)",
01628 &cylindricalBCk1);
01629 opts.units("cylindricalBCk1", N_KCAL);
01630 opts.optional("cylindricalBC", "cylindricalBCexp1", "Exponent for first "
01631 "cylinder potential", &cylindricalBCexp1, 2);
01632 opts.range("cylindricalBCexp1", POSITIVE);
01633
01634
01635
01636 opts.optional("cylindricalBC", "cylindricalBCAxis", "Cylinder axis (defaults to x)",
01637 PARSE_STRING);
01638 opts.require("cylindricalBC", "cylindricalBCCenter",
01639 "Center of cylindrical boundaries", &cylindricalCenter);
01640 opts.require ("cylindricalBC", "cylindricalBCl1", "Length of first cylinder",
01641 &cylindricalBCl1);
01642 opts.range("cylindricalBCl1", POSITIVE);
01643 opts.units("cylindricalBCl1", N_ANGSTROM);
01644 opts.optional ("cylindricalBCl1", "cylindricalBCl2", "Length of second cylinder",
01645 &cylindricalBCl2);
01646 opts.range ("cylindricalBCl2", POSITIVE);
01647 opts.units ("cylindricalBCl2", N_ANGSTROM);
01648
01649
01650 opts.optional("cylindricalBCr1", "cylindricalBCr2", "Radius for second cylinder "
01651 "potential", &cylindricalBCr2);
01652 opts.range("cylindricalBCr2", POSITIVE);
01653 opts.units("cylindricalBCr2", N_ANGSTROM);
01654 opts.require("cylindricalBCr2", "cylindricalBCk2", "Force constant for second "
01655 "cylinder potential (+ is an inward force, - outward)",
01656 &cylindricalBCk2);
01657 opts.units("cylindricalBCk2", N_KCAL);
01658 opts.optional("cylindricalBCr2", "cylindricalBCexp2", "Exponent for second "
01659 "cylinder potential", &cylindricalBCexp2, 2);
01660 opts.range("cylindricalBCexp2", POSITIVE);
01661
01663 opts.optionalB("main", "eFieldOn", "Should an electric field be applied",
01664 &eFieldOn, FALSE);
01665 opts.require("eFieldOn", "eField", "Electric field vector", &eField);
01666 opts.optional("eFieldOn", "eFieldFreq", "Electric field frequency", &eFieldFreq);
01667 opts.optional("eFieldOn", "eFieldPhase", "Electric field phase", &eFieldPhase);
01668
01670 opts.optionalB("main", "stirOn", "Should stirring torque be applied",
01671 &stirOn, FALSE);
01672 opts.optional("stirOn", "stirFilename", "PDB file with flags for "
01673 "stirred atoms (default is the PDB input file)",
01674 PARSE_STRING);
01675 opts.optional("stirOn", "stirredAtomsCol", "Column in the stirredAtomsFile "
01676 "containing the flags (nonzero means fixed);\n"
01677 "default is 'O'", PARSE_STRING);
01678 opts.require("stirOn", "stirStartingTheta", "Stir starting theta offset", &stirStartingTheta);
01679 opts.require("stirOn", "stirK", "Stir force harmonic spring constant", &stirK);
01680
01681 opts.require("stirOn", "stirVel", "Stir angular velocity (deg/timestep)", &stirVel);
01682 opts.require("stirOn", "stirAxis", "Stir axis (direction vector)", &stirAxis);
01683 opts.require("stirOn", "stirPivot", "Stir pivot point (coordinate)", &stirPivot);
01684
01686 opts.optionalB("main", "extraBonds",
01687 "Should extra bonded forces be applied",
01688 &extraBondsOn, FALSE);
01689 opts.optional("extraBonds", "extraBondsFile",
01690 "file with list of extra bonds",
01691 PARSE_MULTIPLES);
01692
01693 }
01694
01695 void SimParameters::config_parser_misc(ParseOptions &opts) {
01696
01698 opts.optional("main", "ldBalancer", "Load balancer",
01699 loadBalancer);
01700 opts.optional("main", "ldbStrategy", "Load balancing strategy",
01701 loadStrategy);
01702 opts.optional("main", "ldbPeriod", "steps between load balancing",
01703 &ldbPeriod);
01704 opts.range("ldbPeriod", POSITIVE);
01705 opts.optional("main", "firstLdbStep", "when to start load balancing",
01706 &firstLdbStep);
01707 opts.range("firstLdbStep", POSITIVE);
01708 opts.optional("main", "lastLdbStep", "when to stop load balancing",
01709 &lastLdbStep);
01710 opts.range("lastLdbStep", POSITIVE);
01711 opts.optional("main", "hybridGroupSize", "Hybrid load balancing group size",
01712 &hybridGroupSize);
01713 opts.optional("main", "ldbBackgroundScaling",
01714 "background load scaling", &ldbBackgroundScaling);
01715 opts.range("ldbBackgroundScaling", NOT_NEGATIVE);
01716 opts.optional("main", "ldbPMEBackgroundScaling",
01717 "PME node background load scaling", &ldbPMEBackgroundScaling);
01718 opts.range("ldbPMEBackgroundScaling", NOT_NEGATIVE);
01719 opts.optional("main", "ldbHomeBackgroundScaling",
01720 "home node background load scaling", &ldbHomeBackgroundScaling);
01721 opts.range("ldbHomeBackgroundScaling", NOT_NEGATIVE);
01722 opts.optional("main", "ldbRelativeGrainsize",
01723 "fraction of average load per compute", &ldbRelativeGrainsize, 0.);
01724 opts.range("ldbRelativeGrainsize", NOT_NEGATIVE);
01725
01726 opts.optional("main", "traceStartStep", "when to start tracing", &traceStartStep);
01727 opts.range("traceStartStep", POSITIVE);
01728 opts.optional("main", "numTraceSteps", "the number of timesteps to be traced", &numTraceSteps);
01729 opts.range("numTraceSteps", POSITIVE);
01730
01731 #ifdef MEASURE_NAMD_WITH_PAPI
01732 opts.optionalB("main", "papiMeasure", "whether use PAPI to measure performacne", &papiMeasure, FALSE);
01733 opts.optional("main", "papiMeasureStartStep", "when to measure performacne using PAPI", &papiMeasureStartStep);
01734 opts.range("papiMeasureStartStep", POSITIVE);
01735 opts.optional("main", "numPapiMeasureSteps", "the number of timesteps to be measured using PAPI", &numPapiMeasureSteps);
01736 opts.range("numPapiMeasureSteps", POSITIVE);
01737 #endif
01738
01739 opts.optionalB("main", "outputMaps", "whether to dump compute map and patch map for analysis just before load balancing", &outputMaps, FALSE);
01740 opts.optionalB("main", "benchTimestep", "whether to do benchmarking timestep in which case final file output is disabled", &benchTimestep, FALSE);
01741 opts.optional("main", "useNodeHelper", "whether to use NodeHelper library to parallelize a loop in a function like OpenMP", &useNodeHelper, 0);
01742 opts.range("useNodeHelper", NOT_NEGATIVE);
01743
01744 opts.optionalB("main", "simulateInitialMapping", "whether to study the initial mapping scheme", &simulateInitialMapping, FALSE);
01745 opts.optional("main", "simulatedPEs", "the number of PEs to be used for studying initial mapping", &simulatedPEs);
01746 opts.range("simulatedPEs", POSITIVE);
01747 opts.optional("main", "simulatedNodeSize", "the node size to be used for studying initial mapping", &simulatedNodeSize);
01748 opts.range("simulatedNodeSize", POSITIVE);
01749
01750 opts.optionalB("main", "ldbUnloadPME", "no load on PME nodes",
01751 &ldbUnloadPME, FALSE);
01752 opts.optionalB("main", "ldbUnloadZero", "no load on pe zero",
01753 &ldbUnloadZero, FALSE);
01754 opts.optionalB("main", "ldbUnloadOne", "no load on pe one",
01755 &ldbUnloadOne, FALSE);
01756 opts.optionalB("main", "noPatchesOnZero", "no patches on pe zero",
01757 &noPatchesOnZero, FALSE);
01758 opts.optionalB("main", "noPatchesOnOne", "no patches on pe one",
01759 &noPatchesOnOne, FALSE);
01760 opts.optionalB("main", "useCompressedPsf", "The structure file psf is in the compressed format",
01761 &useCompressedPsf, FALSE);
01762 opts.optionalB("main", "genCompressedPsf", "Generate the compressed version of the psf file",
01763 &genCompressedPsf, FALSE);
01764 opts.optionalB("main", "usePluginIO", "Use the plugin I/O to load the molecule system",
01765 &usePluginIO, FALSE);
01766 opts.optionalB("main", "mallocTest", "test how much memory all PEs can allocate",
01767 &mallocTest, FALSE);
01768 opts.optionalB("main", "printExclusions", "print exclusion lists to stdout",
01769 &printExclusions, FALSE);
01770 opts.optional("main", "proxySendSpanningTree", "using spanning tree to send proxies",
01771 &proxySendSpanningTree, -1);
01772 opts.optional("main", "proxyRecvSpanningTree", "using spanning tree to receive proxies",
01773 &proxyRecvSpanningTree, -1);
01774 opts.optional("main", "proxyTreeBranchFactor", "the branch factor when building a spanning tree",
01775 &proxyTreeBranchFactor, 4);
01776 opts.optionalB("main", "twoAwayX", "half-size patches in 1st dimension",
01777 &twoAwayX, -1);
01778 opts.optionalB("main", "twoAwayY", "half-size patches in 2nd dimension",
01779 &twoAwayY, -1);
01780 opts.optionalB("main", "twoAwayZ", "half-size patches in 3rd dimension",
01781 &twoAwayZ, -1);
01782 opts.optional("main", "maxPatches", "maximum patch count", &maxPatches, -1);
01783
01785 opts.optional("main", "firsttimestep", "Timestep to start simulation at",
01786 &firstTimestep, 0);
01787 opts.range("firsttimestep", NOT_NEGATIVE);
01788
01790 opts.optionalB("main", "test", "Perform self-tests rather than simulation",
01791 &testOn, FALSE);
01792 opts.optionalB("main", "commOnly", "Do not evaluate forces or integrate",
01793 &commOnly, FALSE);
01794
01795 opts.optionalB("main", "statsOn", "counters in machine layer",
01796 &statsOn, FALSE);
01798 opts.optionalB("main", "hbonds", "Use explicit hydrogen bond term",
01799 &HydrogenBonds, FALSE);
01800 opts.optionalB("hbonds","hbAntecedents","Include Antecedent in hbond term",
01801 &useAntecedent, TRUE);
01802 opts.optional("hbonds","hbAAexp","Hbond AA-A-H angle cos exponential",
01803 &aaAngleExp, 2);
01804 opts.optional("hbonds","hbHAexp","Hbond D-H-A angle cos exponential",
01805 &haAngleExp, 4);
01806 opts.optional("hbonds","hbDistAexp","Hbond A-D dist attractive exponential",
01807 &distAttExp, 4);
01808 opts.optional("hbonds","hbDistRexp","Hbond A-D dist repulstive exponential",
01809 &distRepExp, 6);
01810 opts.optional("hbonds","hbCutoffAngle","Hbond D-H-A cutoff angle",
01811 &dhaCutoffAngle, 100.0);
01812 opts.range("hbCutoffAngle", NOT_NEGATIVE);
01813 opts.optional("hbonds","hbOnAngle","Hbond D-H-A switch function on angle",
01814 &dhaOnAngle, 60.0);
01815 opts.range("hbOnAngle", NOT_NEGATIVE);
01816 opts.optional("hbonds","hbOffAngle","Hbond D-H-A switch function off angle",
01817 &dhaOffAngle, 80.0);
01818 opts.range("hbOffAngle", NOT_NEGATIVE);
01819 opts.optional("hbonds","hbCutoffDist","Hbond A-D cutoff distance",
01820 &daCutoffDist, 7.5);
01821 opts.range("hbCutoffDist", POSITIVE);
01822 opts.units("hbCutoffDist", N_ANGSTROM);
01823 opts.optional("hbonds","hbOnDist","Hbond A-D switch function on distance",
01824 &daOnDist, 5.5);
01825 opts.range("hbOnDist", POSITIVE);
01826 opts.units("hbOnDist", N_ANGSTROM);
01827 opts.optional("hbonds","hbOffDist","Hbond A-D switch function off distance",
01828 &daOffDist, 6.5);
01829 opts.range("hbOffDist", POSITIVE);
01830 opts.units("hbOffDist", N_ANGSTROM);
01831
01832
01833 opts.optionalB("main","IMDon","Connect using IMD?",&IMDon, FALSE);
01834 opts.require("IMDon","IMDport", "Port to which to bind", &IMDport);
01835 opts.range("IMDport",POSITIVE);
01836 opts.require("IMDon","IMDfreq", "Frequency at which to report", &IMDfreq);
01837 opts.range("IMDfreq",POSITIVE);
01838 opts.optionalB("IMDon","IMDwait","Pause until IMD connection?",&IMDwait,
01839 FALSE);
01840 opts.optionalB("IMDon","IMDignore","Ignore forces, etc.?",&IMDignore,
01841 FALSE);
01842
01843
01844 opts.optional("ldBalancer", "maxSelfPart",
01845 "maximum number of self partitions in one patch", &maxSelfPart, 20);
01846 opts.range("maxSelfPart",POSITIVE);
01847 opts.optional("ldBalancer", "maxPairPart",
01848 "maximum number of pair partitions in one patch", &maxPairPart, 8);
01849 opts.range("maxPairPart",POSITIVE);
01850 opts.optional("ldBalancer", "numAtomsSelf",
01851 "maximum number of atoms in one self compute distribution",
01852 &numAtomsSelf, 154);
01853 opts.range("numAtomsSelf",NOT_NEGATIVE);
01854
01855 opts.optional("ldBalancer", "numAtomsSelf2",
01856 "maximum number of atoms in one self compute distribution",
01857 &numAtomsSelf2, 154);
01858 opts.range("numAtomsSelf2",NOT_NEGATIVE);
01859
01860 opts.optional("ldBalancer", "numAtomsPair",
01861 "maximum number of atoms in one pair compute distribution",
01862 &numAtomsPair, 318);
01863 opts.range("numAtomsPair",NOT_NEGATIVE);
01864 opts.optional("ldBalancer", "numAtomsPair2",
01865 "maximum number of atoms in one pair compute distribution",
01866 &numAtomsPair2, 637);
01867 opts.range("numAtomsPair2",NOT_NEGATIVE);
01868 opts.optional("main", "minAtomsPerPatch",
01869 "minimum average atoms per patch",
01870 &minAtomsPerPatch, 40);
01871 opts.range("minAtomsPerPatch",NOT_NEGATIVE);
01872
01873
01874 opts.optional("main", "maxExclusionFlags",
01875 "maximum number of exclusion flags per atom", &maxExclusionFlags, 256);
01876 opts.range("maxExclusionFlags",POSITIVE);
01877
01878 }
01879
01880 #ifdef MEM_OPT_VERSION
01881
01882 extern char *gWorkDir;
01883 #endif
01884
01885 void SimParameters::check_config(ParseOptions &opts, ConfigList *config, char *&cwd) {
01886
01887 int len;
01888 StringList *current;
01889
01890 #ifdef MEM_OPT_VERSION
01891 char *namdWorkDir = NULL;
01892 #endif
01893
01894 if ( opts.defined("obsolete") ) {
01895 iout << iWARN <<
01896 "\"obsolete\" defined, silently ignoring obsolete options\n" << endi;
01897 }
01898
01899
01900 if (opts.defined("cwd"))
01901 {
01902
01903 current = config->find("cwd");
01904
01905 len = strlen(current->data);
01906
01907 if ( CHDIR(current->data) )
01908 {
01909 NAMD_die("chdir() to given cwd failed!");
01910 } else {
01911 iout << iINFO << "Changed directory to " << current->data << "\n" << endi;
01912 }
01913
01914 if (current->data[len-1] != PATHSEP)
01915 len++;
01916
01917 cwd = new char[len+1];
01918
01919 strcpy(cwd, current->data);
01920
01921 if (current->data[strlen(current->data)-1] != PATHSEP)
01922 strcat(cwd, PATHSEPSTR);
01923 }
01924
01925 #ifdef MEM_OPT_VERSION
01926 if(cwd!=NULL)namdWorkDir = cwd;
01927 else namdWorkDir = gWorkDir;
01928 int dirlen = strlen(namdWorkDir);
01929
01930 char *tmpDir;
01931 if(namdWorkDir[dirlen-1]=='/'){
01932 tmpDir = new char[dirlen+1];
01933 tmpDir[dirlen] = 0;
01934 }else{
01935 tmpDir = new char[dirlen+2];
01936 tmpDir[dirlen]='/';
01937 tmpDir[dirlen+1]=0;
01938 }
01939 memcpy(tmpDir, namdWorkDir, dirlen);
01940 namdWorkDir = tmpDir;
01941
01942 delete [] gWorkDir;
01943 #endif
01944
01945
01946
01947 if ( usePluginIO && opts.defined("coordinates") ) {
01948 NAMD_die("Separate coordinates file not allowed with plugin IO, coordinates will be taken from structure file.");
01949 }
01950
01951
01952
01953 if (!amberOn && !gromacsOn) {
01954 #ifndef MEM_OPT_VERSION
01955 if (useCompressedPsf)
01956 NAMD_die("useCompressedPsf requires memory-optimized build!");
01957 if (!usePluginIO && !genCompressedPsf && !opts.defined("coordinates"))
01958 NAMD_die("coordinates not found in the configuration file!");
01959 #else
01960 if(!usePluginIO && !opts.defined("bincoordinates")) {
01961 NAMD_die("bincoordinates not found in the configuration file for the memory optimized version!");
01962 }
01963 if(!usePluginIO && opts.defined("coordinates")) {
01964 NAMD_die("coordinates not allowed in the configuration file for the memory optimized version!");
01965 }
01966 #endif
01967 if (!opts.defined("structure"))
01968 NAMD_die("structure not found in the configuration file!");
01969 if (!opts.defined("parameters"))
01970 NAMD_die("parameters not found in the configuration file!");
01971 }
01972
01973
01974
01975 if (opts.defined("coordinates") && opts.defined("ambercoor"))
01976 NAMD_die("Cannot specify both coordinates and ambercoor!");
01977 #ifndef MEM_OPT_VERSION
01978 if (!genCompressedPsf && !opts.defined("coordinates") && !opts.defined("ambercoor")
01979 && !opts.defined("grocoorfile") && !usePluginIO)
01980 NAMD_die("Coordinate file not found!");
01981 #endif
01982
01983
01984
01985 if (opts.defined("temperature") &&
01986 (opts.defined("velocities") || opts.defined("binvelocities")) )
01987 {
01988 NAMD_die("Cannot specify both an initial temperature and a velocity file");
01989 }
01990
01991 #ifdef MEM_OPT_VERSION
01992
01993 binAtomFile = NULL;
01994 binCoorFile = NULL;
01995 binVelFile = NULL;
01996
01997 char *curfile = NULL;
01998 dirlen = strlen(namdWorkDir);
01999 current = config->find("structure");;
02000 curfile = current->data;
02001 int filelen = strlen(curfile);
02002 if(*curfile == '/' || *curfile=='~') {
02003
02004
02005
02006
02007
02008 binAtomFile = new char[filelen+5];
02009 memcpy(binAtomFile, curfile, filelen);
02010 memcpy(binAtomFile+filelen, ".bin", 4);
02011 binAtomFile[filelen+4] = 0;
02012 }else{
02013 binAtomFile = new char[dirlen+filelen+5];
02014 memcpy(binAtomFile, namdWorkDir, dirlen);
02015 memcpy(binAtomFile+dirlen, curfile, filelen);
02016 memcpy(binAtomFile+dirlen+filelen, ".bin", 4);
02017 binAtomFile[dirlen+filelen+4] = 0;
02018 }
02019
02020 current = config->find("bincoordinates");
02021 curfile = current->data;
02022 filelen = strlen(curfile);
02023 if(*curfile == '/' || *curfile=='~') {
02024 binCoorFile = new char[filelen+1];
02025 memcpy(binCoorFile, curfile, filelen);
02026 binCoorFile[filelen] = 0;
02027 }else{
02028 binCoorFile = new char[dirlen+filelen+1];
02029 memcpy(binCoorFile, namdWorkDir, dirlen);
02030 memcpy(binCoorFile+dirlen, curfile, filelen);
02031 binCoorFile[dirlen+filelen] = 0;
02032 }
02033
02034 if(opts.defined("binvelocities")){
02035 current = config->find("binvelocities");
02036 curfile = current->data;
02037 filelen = strlen(curfile);
02038 if(*curfile == '/' || *curfile=='~') {
02039 binVelFile = new char[filelen+1];
02040 memcpy(binVelFile, curfile, filelen);
02041 binVelFile[filelen] = 0;
02042 }else{
02043 binVelFile = new char[dirlen+filelen+1];
02044 memcpy(binVelFile, namdWorkDir, dirlen);
02045 memcpy(binVelFile+dirlen, curfile, filelen);
02046 binVelFile[dirlen+filelen] = 0;
02047 }
02048 }
02049
02050
02051 if(outputFilename[0] != '/' && outputFilename[0]!='~') {
02052 filelen = strlen(outputFilename);
02053 char *tmpout = new char[filelen];
02054 memcpy(tmpout, outputFilename, filelen);
02055 CmiAssert(filelen+dirlen <= 120);
02056 memcpy(outputFilename, namdWorkDir, dirlen);
02057 memcpy(outputFilename+dirlen, tmpout, filelen);
02058 outputFilename[filelen+dirlen] = 0;
02059 delete [] tmpout;
02060 }
02061
02062 if ( dcdFrequency && opts.defined("dcdfile") &&
02063 dcdFilename[0] != '/' && dcdFilename[0]!='~' ) {
02064 filelen = strlen(dcdFilename);
02065 char *tmpout = new char[filelen];
02066 memcpy(tmpout, dcdFilename, filelen);
02067 CmiAssert(filelen+dirlen <= 120);
02068 memcpy(dcdFilename, namdWorkDir, dirlen);
02069 memcpy(dcdFilename+dirlen, tmpout, filelen);
02070 dcdFilename[filelen+dirlen] = 0;
02071 delete [] tmpout;
02072 }
02073
02074 if ( velDcdFrequency && opts.defined("veldcdfile") &&
02075 velDcdFilename[0] != '/' && velDcdFilename[0]!='~' ) {
02076 filelen = strlen(velDcdFilename);
02077 char *tmpout = new char[filelen];
02078 memcpy(tmpout, velDcdFilename, filelen);
02079 CmiAssert(filelen+dirlen <= 120);
02080 memcpy(velDcdFilename, namdWorkDir, dirlen);
02081 memcpy(velDcdFilename+dirlen, tmpout, filelen);
02082 velDcdFilename[filelen+dirlen] = 0;
02083 delete [] tmpout;
02084 }
02085
02086 if ( forceDcdFrequency && opts.defined("forcedcdfile") &&
02087 forceDcdFilename[0] != '/' && forceDcdFilename[0]!='~' ) {
02088 filelen = strlen(forceDcdFilename);
02089 char *tmpout = new char[filelen];
02090 memcpy(tmpout, forceDcdFilename, filelen);
02091 CmiAssert(filelen+dirlen <= 120);
02092 memcpy(forceDcdFilename, namdWorkDir, dirlen);
02093 memcpy(forceDcdFilename+dirlen, tmpout, filelen);
02094 forceDcdFilename[filelen+dirlen] = 0;
02095 delete [] tmpout;
02096 }
02097
02098 if ( restartFrequency && opts.defined("restartname") &&
02099 restartFilename[0] != '/' && restartFilename[0]!='~' ) {
02100 filelen = strlen(restartFilename);
02101 char *tmpout = new char[filelen];
02102 memcpy(tmpout, restartFilename, filelen);
02103 CmiAssert(filelen+dirlen <= 120);
02104 memcpy(restartFilename, namdWorkDir, dirlen);
02105 memcpy(restartFilename+dirlen, tmpout, filelen);
02106 restartFilename[filelen+dirlen] = 0;
02107 delete [] tmpout;
02108 }
02109
02110 delete [] namdWorkDir;
02111
02112 if (opts.defined("numinputprocs")) {
02113 if(numinputprocs > CkNumPes()) {
02114 iout << iWARN << "The number of input processors exceeds the total number of processors. Resetting to half of the number of total processors.\n" << endi;
02115 numinputprocs = (CkNumPes()>>1)+(CkNumPes()&1);
02116 }
02117 }
02118
02119 if (opts.defined("numoutputprocs")) {
02120 if(numoutputprocs > CkNumPes()) {
02121 iout << iWARN << "The number of output processors exceeds the total number of processors. Resetting to half of the number of total processors.\n" << endi;
02122 numoutputprocs = (CkNumPes()>>1)+(CkNumPes()&1);
02123 }
02124 }
02125
02126 #if !OUTPUT_SINGLE_FILE
02127
02128 create_output_directories("coor");
02129 create_output_directories("vel");
02130 if(dcdFrequency) {
02131 create_output_directories("dcd");
02132 if(opts.defined("dcdfile")){
02133 iout << iWARN << "The dcd file output has been changed to directory: " << outputFilename << ".\n" << endi;
02134 }
02135 }
02136 if (velDcdFrequency) {
02137 create_output_directories("veldcd");
02138 if(opts.defined("veldcdfile")){
02139 iout << iWARN << "The veldcd file output has been changed to directory: " << outputFilename << ".\n" << endi;
02140 }
02141 }
02142 if (forceDcdFrequency) {
02143 create_output_directories("forcedcd");
02144 if(opts.defined("forcedcdfile")){
02145 iout << iWARN << "The forcedcd file output has been changed to directory: " << outputFilename << ".\n" << endi;
02146 }
02147 }
02148 #endif
02149 #endif
02150
02151 if (! opts.defined("auxFile")) {
02152 strcpy(auxFilename,outputFilename);
02153 strcat(auxFilename,".aux");
02154 }
02155
02156
02157 if (dcdFrequency) {
02158 if (! opts.defined("dcdfile")) {
02159 strcpy(dcdFilename,outputFilename);
02160 strcat(dcdFilename,".dcd");
02161 }
02162 } else {
02163 dcdFilename[0] = STRINGNULL;
02164 }
02165
02166 if (velDcdFrequency) {
02167 if (! opts.defined("veldcdfile")) {
02168 strcpy(velDcdFilename,outputFilename);
02169 strcat(velDcdFilename,".veldcd");
02170 }
02171 } else {
02172 velDcdFilename[0] = STRINGNULL;
02173 }
02174
02175 if (forceDcdFrequency) {
02176 if (! opts.defined("forcedcdfile")) {
02177 strcpy(forceDcdFilename,outputFilename);
02178 strcat(forceDcdFilename,".forcedcd");
02179 }
02180 } else {
02181 forceDcdFilename[0] = STRINGNULL;
02182 }
02183
02184 if (xstFrequency) {
02185 if (! opts.defined("xstfile")) {
02186 strcpy(xstFilename,outputFilename);
02187 strcat(xstFilename,".xst");
02188 }
02189 } else {
02190 xstFilename[0] = STRINGNULL;
02191 }
02192
02193 if (restartFrequency) {
02194 if (! opts.defined("restartname")) {
02195 strcpy(restartFilename,outputFilename);
02196 strcat(restartFilename,".restart");
02197 }
02198 } else {
02199 restartFilename[0] = STRINGNULL;
02200 restartSave = FALSE;
02201 binaryRestart = FALSE;
02202 }
02203
02204 if (storeComputeMap || loadComputeMap) {
02205 if (! opts.defined("computeMapFile")) {
02206 strcpy(computeMapFilename,"computeMapFile");
02207 strcat(computeMapFilename,".txt");
02208 }
02209 }
02210
02211
02212 if (!amberOn)
02213 {
02215 if (!paraTypeXplorOn && !paraTypeCharmmOn)
02216 {
02217 paraTypeXplorOn = TRUE;
02218 }
02220 if (paraTypeXplorOn && paraTypeCharmmOn)
02221 {
02222 NAMD_die("Please specify either XPLOR or CHARMM format for parameters!");
02223 }
02224
02225 }
02226
02227
02228
02229 if (!(minimizeOn||minimizeCGOn) && !opts.defined("temperature") &&
02230 !opts.defined("velocities") && !opts.defined("binvelocities") )
02231 {
02232 NAMD_die("Must have either an initial temperature or a velocity file");
02233 }
02234
02235 if (minimizeOn||minimizeCGOn) { initialTemp = 0.0; }
02236 if (opts.defined("velocities") || opts.defined("binvelocities") )
02237 {
02238 initialTemp = -1.0;
02239 }
02240
02242
02243 if ( opts.defined("extendedSystem") )
02244 {
02245 current = config->find("extendedSystem");
02246
02247 iout << iINFO << "EXTENDED SYSTEM FILE "
02248 << current->data << "\n" << endi;
02249
02250 ifstream xscFile(current->data);
02251 if ( ! xscFile ) NAMD_die("Unable to open extended system file.\n");
02252
02253 char labels[1024];
02254 do {
02255 if ( ! xscFile ) NAMD_die("Error reading extended system file.\n");
02256 xscFile.getline(labels,1023);
02257 } while ( strncmp(labels,"#$LABELS ",9) );
02258
02259 int a_x, a_y, a_z, b_x, b_y, b_z, c_x, c_y, c_z;
02260 a_x = a_y = a_z = b_x = b_y = b_z = c_x = c_y = c_z = -1;
02261 int o_x, o_y, o_z, s_u, s_v, s_w, s_x, s_y, s_z;
02262 o_x = o_y = o_z = s_u = s_v = s_w = s_x = s_y = s_z = -1;
02263
02264 int pos = 0;
02265 char *l_i = labels + 8;
02266 while ( *l_i ) {
02267 if ( *l_i == ' ' ) { ++l_i; continue; }
02268 char *l_i2;
02269 for ( l_i2 = l_i; *l_i2 && *l_i2 != ' '; ++l_i2 );
02270 if ( (l_i2 - l_i) == 3 && (l_i[1] == '_') ) {
02271 if (l_i[0] == 'a' && l_i[2] == 'x') a_x = pos;
02272 if (l_i[0] == 'a' && l_i[2] == 'y') a_y = pos;
02273 if (l_i[0] == 'a' && l_i[2] == 'z') a_z = pos;
02274 if (l_i[0] == 'b' && l_i[2] == 'x') b_x = pos;
02275 if (l_i[0] == 'b' && l_i[2] == 'y') b_y = pos;
02276 if (l_i[0] == 'b' && l_i[2] == 'z') b_z = pos;
02277 if (l_i[0] == 'c' && l_i[2] == 'x') c_x = pos;
02278 if (l_i[0] == 'c' && l_i[2] == 'y') c_y = pos;
02279 if (l_i[0] == 'c' && l_i[2] == 'z') c_z = pos;
02280 if (l_i[0] == 'o' && l_i[2] == 'x') o_x = pos;
02281 if (l_i[0] == 'o' && l_i[2] == 'y') o_y = pos;
02282 if (l_i[0] == 'o' && l_i[2] == 'z') o_z = pos;
02283 if (l_i[0] == 's' && l_i[2] == 'u') s_u = pos;
02284 if (l_i[0] == 's' && l_i[2] == 'v') s_v = pos;
02285 if (l_i[0] == 's' && l_i[2] == 'w') s_w = pos;
02286 if (l_i[0] == 's' && l_i[2] == 'x') s_x = pos;
02287 if (l_i[0] == 's' && l_i[2] == 'y') s_y = pos;
02288 if (l_i[0] == 's' && l_i[2] == 'z') s_z = pos;
02289 }
02290 ++pos;
02291 l_i = l_i2;
02292 }
02293 int numpos = pos;
02294
02295 for ( pos = 0; pos < numpos; ++pos ) {
02296 double tmp;
02297 xscFile >> tmp;
02298 if ( ! xscFile ) NAMD_die("Error reading extended system file.\n");
02299 if ( pos == a_x ) cellBasisVector1.x = tmp;
02300 if ( pos == a_y ) cellBasisVector1.y = tmp;
02301 if ( pos == a_z ) cellBasisVector1.z = tmp;
02302 if ( pos == b_x ) cellBasisVector2.x = tmp;
02303 if ( pos == b_y ) cellBasisVector2.y = tmp;
02304 if ( pos == b_z ) cellBasisVector2.z = tmp;
02305 if ( pos == c_x ) cellBasisVector3.x = tmp;
02306 if ( pos == c_y ) cellBasisVector3.y = tmp;
02307 if ( pos == c_z ) cellBasisVector3.z = tmp;
02308 if ( pos == o_x ) cellOrigin.x = tmp;
02309 if ( pos == o_y ) cellOrigin.y = tmp;
02310 if ( pos == o_z ) cellOrigin.z = tmp;
02311 if ( pos == s_u ) strainRate2.x = tmp;
02312 if ( pos == s_v ) strainRate2.y = tmp;
02313 if ( pos == s_w ) strainRate2.z = tmp;
02314 if ( pos == s_x ) strainRate.x = tmp;
02315 if ( pos == s_y ) strainRate.y = tmp;
02316 if ( pos == s_z ) strainRate.z = tmp;
02317 }
02318
02319 }
02320
02321 if ( LJcorrection && ! cellBasisVector3.length2() ) {
02322 NAMD_die("Can't use LJ tail corrections without periodic boundary conditions!");
02323 }
02324
02325 if ( cellBasisVector3.length2() && ! cellBasisVector2.length2() ) {
02326 NAMD_die("Used cellBasisVector3 without cellBasisVector2!");
02327 }
02328
02329 if ( cellBasisVector2.length2() && ! cellBasisVector1.length2() ) {
02330 NAMD_die("Used cellBasisVector2 without cellBasisVector1!");
02331 }
02332
02333 if ( cellOrigin.length2() && ! cellBasisVector1.length2() ) {
02334 NAMD_die("Used cellOrigin without cellBasisVector1!");
02335 }
02336
02337 lattice.set(cellBasisVector1,cellBasisVector2,cellBasisVector3,cellOrigin);
02338
02339 if (! opts.defined("DCDunitcell")) {
02340 dcdUnitCell = lattice.a_p() && lattice.b_p() && lattice.c_p();
02341 }
02342
02343 char s[129];
02344
02346 if ( ! opts.defined("cylindricalBCAxis") )
02347 {
02348 cylindricalBCAxis = 'x';
02349 }
02350 else
02351 {
02352 opts.get("cylindricalBCAxis", s);
02353
02354 if (!strcasecmp(s, "x"))
02355 {
02356 cylindricalBCAxis = 'x';
02357 }
02358 else if (!strcasecmp(s, "y"))
02359 {
02360 cylindricalBCAxis = 'y';
02361 }
02362 else if (!strcasecmp(s, "z"))
02363 {
02364 cylindricalBCAxis = 'z';
02365 }
02366 else
02367 {
02368 char err_msg[128];
02369
02370 sprintf(err_msg, "Illegal value '%s' for 'cylindricalBCAxis' in configuration file", s);
02371 NAMD_die(err_msg);
02372 }
02373 }
02374
02375 if (!opts.defined("splitPatch"))
02376 {
02377 splitPatch = SPLIT_PATCH_HYDROGEN;
02378 }
02379 else
02380 {
02381 opts.get("splitPatch", s);
02382 if (!strcasecmp(s, "position"))
02383 splitPatch = SPLIT_PATCH_POSITION;
02384 else if (!strcasecmp(s,"hydrogen"))
02385 splitPatch = SPLIT_PATCH_HYDROGEN;
02386 else
02387 {
02388 char err_msg[129];
02389 sprintf(err_msg,
02390 "Illegal value '%s' for 'splitPatch' in configuration file",
02391 s);
02392 NAMD_die(err_msg);
02393 }
02394 }
02395
02397 opts.get("exclude", s);
02398
02399 if (!strcasecmp(s, "none"))
02400 {
02401 exclude = NONE;
02402 splitPatch = SPLIT_PATCH_POSITION;
02403 }
02404 else if (!strcasecmp(s, "1-2"))
02405 {
02406 exclude = ONETWO;
02407 splitPatch = SPLIT_PATCH_POSITION;
02408 }
02409 else if (!strcasecmp(s, "1-3"))
02410 {
02411 exclude = ONETHREE;
02412 }
02413 else if (!strcasecmp(s, "1-4"))
02414 {
02415 exclude = ONEFOUR;
02416 }
02417 else if (!strcasecmp(s, "scaled1-4"))
02418 {
02419 exclude = SCALED14;
02420 }
02421 else
02422 {
02423 char err_msg[128];
02424
02425 sprintf(err_msg, "Illegal value '%s' for 'exclude' in configuration file",
02426 s);
02427 NAMD_die(err_msg);
02428 }
02429
02430 if (scale14 != 1.0 && exclude != SCALED14)
02431 {
02432 iout << iWARN << "Exclude is not scaled1-4; 1-4scaling ignored.\n" << endi;
02433 }
02434
02435
02436 if (!opts.defined("waterModel")) {
02437 watmodel = WAT_TIP3;
02438 } else {
02439 opts.get("waterModel", s);
02440 if (!strncasecmp(s, "tip4", 4)) {
02441 iout << iINFO << "Using TIP4P water model.\n" << endi;
02442 watmodel = WAT_TIP4;
02443 } else if (!strncasecmp(s, "tip3", 4)) {
02444 iout << iINFO << "Using TIP3P water model.\n" << endi;
02445 } else if (!strncasecmp(s, "swm4", 4)) {
02446 iout << iINFO << "Using SWM4-DP water model.\n" << endi;
02447 } else {
02448 char err_msg[128];
02449 sprintf(err_msg,
02450 "Illegal value %s for 'waterModel' in configuration file", s);
02451 NAMD_die(err_msg);
02452 }
02453 }
02454 if (watmodel == WAT_SWM4 && !drudeOn) {
02455 NAMD_die("Must have 'drudeOn' enabled to use SWM4-DP water model.");
02456 }
02457 if (drudeOn && watmodel != WAT_SWM4) {
02458 watmodel = WAT_SWM4;
02459 iout << iWARN
02460 << "Setting water model to 'swm4' (SWM4-DP) for Drude polarization.\n"
02461 << endi;
02462 }
02463
02464
02465 if (watmodel == WAT_SWM4) {
02466 lonepairs = TRUE;
02467 }
02468
02469
02470 if (goForcesOn) {
02471 iout << iINFO << "Go forces are on\n" << endi;
02472
02473 int * gomethod = &goMethod;
02474 if (!opts.defined("GoMethod")) {
02475 *gomethod = 0;
02476
02477 } else {
02478 opts.get("GoMethod",s);
02479
02480 *gomethod = atoi(s);
02481 }
02482 if (!strcasecmp(s, "matrix")) {
02483 goMethod = 1;
02484
02485 } else if (!strcasecmp(s, "sparse")) {
02486 goMethod = 2;
02487
02488 } else if (!strcasecmp(s, "lowmem")) {
02489 goMethod = 3;
02490
02491
02492
02493
02494
02495
02496
02497 }
02498 else {
02499 char err_msg[129];
02500 sprintf(err_msg,
02501 "Illegal value '%s' for 'GoMethod' in configuration file",
02502 s);
02503 NAMD_die(err_msg);
02504 }
02505 }
02506
02507
02508
02509 if (!opts.defined("MTSAlgorithm"))
02510 {
02511 MTSAlgorithm = VERLETI;
02512 }
02513 else
02514 {
02515 opts.get("MTSAlgorithm", s);
02516
02517 if (!strcasecmp(s, "naive"))
02518 {
02519 MTSAlgorithm = NAIVE;
02520 }
02521 else if (!strcasecmp(s, "constant"))
02522 {
02523 MTSAlgorithm = NAIVE;
02524 }
02525 else if (!strcasecmp(s, "impulse"))
02526 {
02527 MTSAlgorithm = VERLETI;
02528 }
02529 else if (!strcasecmp(s, "verleti"))
02530 {
02531 MTSAlgorithm = VERLETI;
02532 }
02533 else
02534 {
02535 char err_msg[129];
02536
02537 sprintf(err_msg,
02538 "Illegal value '%s' for 'MTSAlgorithm' in configuration file",
02539 s);
02540 NAMD_die(err_msg);
02541 }
02542 }
02543
02544
02545 if (!opts.defined("longSplitting"))
02546 {
02547 longSplitting = C1;
02548 }
02549 else
02550 {
02551 opts.get("longSplitting", s);
02552 if (!strcasecmp(s, "sharp"))
02553 longSplitting = SHARP;
02554 else if (!strcasecmp(s, "xplor"))
02555 longSplitting = XPLOR;
02556 else if (!strcasecmp(s, "c1"))
02557 longSplitting = C1;
02558 else if (!strcasecmp(s, "c2"))
02559 longSplitting = C2;
02560 else
02561 {
02562 char err_msg[129];
02563
02564 sprintf(err_msg,
02565 "Illegal value '%s' for 'longSplitting' in configuration file",
02566 s);
02567 NAMD_die(err_msg);
02568 }
02569 }
02570
02571
02572 if (!opts.defined("rigidBonds"))
02573 {
02574 rigidBonds = RIGID_NONE;
02575 }
02576 else
02577 {
02578 opts.get("rigidBonds", s);
02579 if (!strcasecmp(s, "all"))
02580 {
02581 rigidBonds = RIGID_ALL;
02582 }
02583 else if (!strcasecmp(s, "water"))
02584 {
02585 rigidBonds = RIGID_WATER;
02586 }
02587 else if (!strcasecmp(s, "none"))
02588 {
02589 rigidBonds = RIGID_NONE;
02590 }
02591 else
02592 {
02593 char err_msg[256];
02594 sprintf(err_msg,
02595 "Illegal value '%s' for 'rigidBonds' in configuration file", s);
02596 NAMD_die(err_msg);
02597 }
02598 }
02599
02600
02601 if (switchingActive)
02602 {
02603
02604 if (!opts.defined("switchDist")) {
02605 NAMD_die("switchDist must be defined when switching is enabled");
02606 }
02607
02608 if ( (switchingDist>cutoff) || (switchingDist<0) )
02609 {
02610 char err_msg[129];
02611
02612 sprintf(err_msg,
02613 "switchDist muct be between 0 and cutoff, which is %f", cutoff);
02614 NAMD_die(err_msg);
02615 }
02616
02617 }
02618
02619 if ( martiniSwitching )
02620 {
02621 if ( ! switchingActive )
02622 {
02623 NAMD_die("martiniSwitching requires switching");
02624 }
02625 if ( vdwForceSwitching )
02626 {
02627 NAMD_die("martiniSwitching and vdwForceSwitching are exclusive to one another. Select only one.");
02628 }
02629 if ( dielectric != 15.0 && ! martiniDielAllow )
02630 {
02631 iout << iWARN << "USE DIELECTRIC OF 15.0 WITH MARTINI.\n";
02632 iout << iWARN << "SETTING dielectric 15.0\n";
02633 iout << iWARN << "FOR NON-STANDARD DIELECTRIC WITH MARTINI, SET: martiniDielAllow on\n";
02634 dielectric = 15.0;
02635 }
02636 if ( ! cosAngles )
02637 {
02638 iout << iWARN << "USE COSINE BASED ANGLES WITH MARTINI.\n";
02639 iout << iWARN << "SETTING cosAngles on\n";
02640 cosAngles = TRUE;
02641 }
02642 if ( PMEOn )
02643 {
02644 NAMD_die("Do not use Particle Mesh Ewald with Martini. Set: PME off");
02645 }
02646 if ( MSMOn )
02647 {
02648 NAMD_die("Do not use Multilevel Summation Method with Martini. Set: MSM off");
02649 }
02650 }
02651
02652
02653 if (!opts.defined("pairlistDist"))
02654 {
02655 pairlistDist = cutoff;
02656 }
02657 else if (pairlistDist < cutoff)
02658 {
02659 NAMD_die("pairlistDist must be >= cutoff distance");
02660 }
02661
02662 patchDimension = pairlistDist;
02663
02664 if ( splitPatch == SPLIT_PATCH_HYDROGEN ) {
02665 patchDimension += hgroupCutoff;
02666 }
02667
02668 BigReal defaultMargin = 0.0;
02669 if (berendsenPressureOn || langevinPistonOn) {
02670 defaultMargin = ( useFlexibleCell ? 0.06 : 0.03 ) * patchDimension;
02671 }
02672 if ( margin == XXXBIGREAL ) {
02673 margin = defaultMargin;
02674 }
02675 if ( defaultMargin != 0.0 && margin == 0.0 ) {
02676 margin = defaultMargin;
02677 iout << iWARN << "ALWAYS USE NON-ZERO MARGIN WITH CONSTANT PRESSURE!\n";
02678 iout << iWARN << "CHANGING MARGIN FROM 0 to " << margin << "\n" << endi;
02679 }
02680
02681 patchDimension += margin;
02682
02683
02684 if (GBISOn) {
02685
02686 if (fullDirectOn) {
02687 NAMD_die("GBIS not compatible with FullDirect");
02688 }
02689 if (PMEOn) {
02690 NAMD_die("GBIS not compatible with PME");
02691 }
02692 if (MSMOn) {
02693 NAMD_die("GBIS not compatible with MSM");
02694 }
02695 if (alchOn) {
02696 NAMD_die("GBIS not compatible with Alchemical Transformations");
02697 }
02698 if (lesOn) {
02699 NAMD_die("GBIS not compatible with Locally Enhanced Sampling");
02700 }
02701 if (FMAOn) {
02702 NAMD_die("GBIS not compatible with FMA");
02703 }
02704 if (drudeOn) {
02705 NAMD_die("GBIS not compatible with Drude Polarization");
02706 }
02707 if (fixedAtomsOn) {
02708 #ifdef NAMD_CUDA
02709 NAMD_die("GBIS CUDA not yet compatible with fixed atoms");
02710 #endif
02711 }
02712
02713 if (alpha_cutoff > patchDimension) {
02714 patchDimension = alpha_cutoff;
02715 }
02716
02717 BigReal tmp = (initialTemp > 0) ? initialTemp : 300;
02718 kappa = 50.29216*sqrt(ion_concentration/solvent_dielectric/tmp);
02719
02720 }
02721
02722 if (LCPOOn) {
02723 #ifdef MEM_OPT_VERSION
02724 NAMD_die("SASA not yet available for memory optimized builds");
02725 #endif
02726 if ( lattice.volume() > 0 ) {
02727 NAMD_die("SASA does not yet support periodic boundary conditions.");
02728 }
02729
02730
02731 }
02732
02733
02734
02735 if ( dihedralOn ) globalOn = TRUE;
02736
02737
02738 if ((minimizeOn||minimizeCGOn) && langevinOn)
02739 {
02740 NAMD_die("Minimization and Langevin dynamics are mutually exclusive dynamics modes");
02741 }
02742
02743 if ((minimizeOn||minimizeCGOn) && tCoupleOn)
02744 {
02745 NAMD_die("Minimization and temperature coupling are mutually exclusive dynamics modes");
02746 }
02747
02748
02749 if ((minimizeOn||minimizeCGOn) && loweAndersenOn)
02750 {
02751 NAMD_die("Minimization and Lowe-Andersen dynamics are mutually exclusive dynamics modes");
02752 }
02753 #ifdef NAMD_CUDA
02754 if (loweAndersenOn) {
02755 NAMD_die("Lowe-Andersen dynamics not compatible with CUDA at this time");
02756 }
02757 #endif
02758
02759
02760 if (langevinOn && tCoupleOn)
02761 {
02762 NAMD_die("Langevin dynamics and temperature coupling are mutually exclusive dynamics modes");
02763 }
02764
02765
02766 if (loweAndersenOn && (langevinOn || tCoupleOn))
02767 {
02768 NAMD_die("Lowe-Andersen dynamics, Langevin dynamics and temperature coupling are mutually exclusive dynamics modes");
02769 }
02770
02771
02772 if (tCoupleOn && opts.defined("rescaleFreq") )
02773 {
02774 NAMD_die("Temperature coupling and temperature rescaling are mutually exclusive");
02775 }
02776
02777 if (globalOn && CkNumPes() > 1)
02778 {
02779 NAMD_die("Global integration does not run in parallel (yet).");
02780 }
02781
02782 if (COLDOn && langevinOn)
02783 {
02784 NAMD_die("COLD and Langevin dynamics are mutually exclusive dynamics modes");
02785 }
02786 if (COLDOn && minimizeOn)
02787 {
02788 NAMD_die("COLD and minimization are mutually exclusive dynamics modes");
02789 }
02790 if (COLDOn && tCoupleOn)
02791 {
02792 NAMD_die("COLD and temperature coupling are mutually exclusive dynamics modes");
02793 }
02794 if (COLDOn && opts.defined("rescaleFreq"))
02795 {
02796 NAMD_die("COLD and velocity rescaling are mutually exclusive dynamics modes");
02797 }
02798
02799 if (splitPatch == SPLIT_PATCH_POSITION && mollyOn )
02800 {
02801 NAMD_die("splitPatch hydrogen is required for MOLLY");
02802 }
02803
02804 if (splitPatch == SPLIT_PATCH_POSITION && rigidBonds != RIGID_NONE)
02805 {
02806 NAMD_die("splitPatch hydrogen is required for rigidBonds");
02807 }
02808
02809
02810
02811 if (minimizeOn && (maximumMove == 0.0))
02812 {
02813 maximumMove = 0.75 * pairlistDist/stepsPerCycle;
02814 }
02815 if (adaptTempOn) {
02816 if (!adaptTempRescale && !adaptTempLangevin)
02817 NAMD_die("Adaptive tempering needs to be coupled to either the Langevin thermostat or velocity rescaling.");
02818 if (opts.defined("adaptTempInFile") && (opts.defined("adaptTempTmin") ||
02819 opts.defined("adaptTempTmax") ||
02820 adaptTempBins != 0))
02821 NAMD_die("cannot simultaneously specify adaptTempInFile and any of {adaptTempTmin, adaptTempTmax,adaptTempBins} as these are read from the input file");
02822 if (!opts.defined("adaptTempInFile") && !(opts.defined("adaptTempTmin") &&
02823 opts.defined("adaptTempTmax") &&
02824 adaptTempBins != 0 ))
02825 NAMD_die("Need to specify either adaptTempInFile or all of {adaptTempTmin, adaptTempTmax,adaptTempBins} if adaptTempMD is on.");
02826 }
02827 if (langevinOn) {
02828 if ( ! opts.defined("langevinDamping") ) langevinDamping = 0.0;
02829 if ( ! opts.defined("langevinHydrogen") ) langevinHydrogen = TRUE;
02830 if ( (opts.defined("langevinDamping") || opts.defined("langevinHydrogen"))
02831 && (opts.defined("langevinFile") || opts.defined("langevinCol")) )
02832 NAMD_die("To specify Langevin dynamics parameters, use either langevinDamping and langevinHydrogen or langevinFile and langevinCol. Do not combine them.");
02833 if ( opts.defined("langevinHydrogen") && langevinDamping == 0.0 )
02834 NAMD_die("langevinHydrogen requires langevinDamping to be set.");
02835 }
02836
02837
02838 if (loweAndersenOn) {
02839 if (!opts.defined("loweAndersenRate")) loweAndersenRate = 100;
02840 if (!opts.defined("loweAndersenCutoff")) loweAndersenCutoff = 2.7;
02841 }
02842
02843
02844 if (opts.defined("rescaleFreq"))
02845 {
02846 if (!opts.defined("rescaleTemp"))
02847 {
02848 if (opts.defined("temperature"))
02849 {
02850 rescaleTemp = initialTemp;
02851 }
02852 else
02853 {
02854 NAMD_die("Must give a rescale temperature if rescaleFreq is defined");
02855 }
02856 }
02857 }
02858 else
02859 {
02860 rescaleFreq = -1;
02861 rescaleTemp = 0.0;
02862 }
02863
02864 if (opts.defined("rescaleTemp"))
02865 {
02866 if (!opts.defined("rescaleFreq"))
02867 {
02868 NAMD_die("Must give a rescale freqency if rescaleTemp is given");
02869 }
02870 }
02871
02872 if ((minimizeOn||minimizeCGOn) && rescaleFreq > 0)
02873 {
02874 NAMD_die("Minimization and temperature rescaling are mutually exclusive dynamics modes");
02875 }
02876
02877 if (opts.defined("reassignFreq"))
02878 {
02879 if (!opts.defined("reassignTemp"))
02880 {
02881 if (opts.defined("temperature"))
02882 {
02883 reassignTemp = initialTemp;
02884 }
02885 else
02886 {
02887 NAMD_die("Must give a reassign temperature if reassignFreq is defined");
02888 }
02889 }
02890 }
02891 else
02892 {
02893 reassignFreq = -1;
02894 reassignTemp = 0.0;
02895 }
02896
02897 if (opts.defined("reassignTemp"))
02898 {
02899 if (!opts.defined("reassignFreq"))
02900 {
02901 NAMD_die("Must give a reassignment freqency if reassignTemp is given");
02902 }
02903 }
02904
02905 if (opts.defined("reassignIncr"))
02906 {
02907 if (!opts.defined("reassignFreq"))
02908 {
02909 NAMD_die("Must give a reassignment freqency if reassignIncr is given");
02910 }
02911 }
02912 else
02913 {
02914 reassignIncr = 0.0;
02915 }
02916
02917 if (opts.defined("reassignHold"))
02918 {
02919 if (!opts.defined("reassignIncr"))
02920 {
02921 NAMD_die("Must give a reassignment increment if reassignHold is given");
02922 }
02923 }
02924 else
02925 {
02926 reassignHold = 0.0;
02927 }
02928
02929 if ((minimizeOn||minimizeCGOn) && reassignFreq > 0)
02930 {
02931 NAMD_die("Minimization and temperature reassignment are mutually exclusive dynamics modes");
02932 }
02933
02934 if (!opts.defined("seed"))
02935 {
02936 randomSeed = (unsigned int) time(NULL);
02937 }
02938
02939
02940
02941 alchFepOn = FALSE;
02942 alchThermIntOn = FALSE;
02943
02944 if (alchOn) {
02945
02946 if (!opts.defined("alchType"))
02947 {
02948 NAMD_die("Must define type of alchemical simulation: fep or ti\n");
02949 }
02950 else
02951 {
02952 opts.get("alchType",s);
02953 if (!strcasecmp(s, "fep"))
02954 {
02955 alchFepOn = TRUE;
02956 }
02957 else if (!strcasecmp(s, "ti"))
02958 {
02959 alchThermIntOn = TRUE;
02960 }
02961 else
02962 {
02963 NAMD_die("Unknown type of alchemical simulation; choices are fep or ti\n");
02964 }
02965 }
02966
02967 if (rescaleFreq > 0) alchTemp = rescaleTemp;
02968 else if (reassignFreq > 0) alchTemp = reassignTemp;
02969 else if (langevinOn) alchTemp = langevinTemp;
02970 else if (tCoupleOn) alchTemp = tCoupleTemp;
02971 else NAMD_die("Alchemical FEP can be performed only in constant temperature simulations\n");
02972
02973 if (reassignFreq > 0 && reassignIncr != 0)
02974 NAMD_die("reassignIncr cannot be used in alchemical simulations\n");
02975
02976 if (alchLambda < 0.0 || alchLambda > 1.0 || alchLambda2 < 0.0 || alchLambda2 > 1.0)
02977 NAMD_die("Alchemical lambda values should be in the range [0.0, 1.0]\n");
02978
02979 if ( alchOn && alchVdwLambdaEnd > 1.0)
02980 NAMD_die("Gosh tiny Elvis, you kicked soft-core in the van der Waals! alchVdwLambdaEnd should be in the range [0.0, 1.0]\n");
02981
02982 if ( alchOn && alchElecLambdaStart > 1.0)
02983 NAMD_die("alchElecLambdaStart should be in the range [0.0, 1.0]\n");
02984
02985 if (alchFepOn)
02986 {
02987 if (!opts.defined("alchoutfile")) {
02988 strcpy(alchOutFile, outputFilename);
02989 strcat(alchOutFile, ".fep");
02990 }
02991 if (alchFepWCARepuOn && alchFepWCADispOn)
02992 NAMD_die("With WCA decomposition, repulsion and dispersion can NOT be in the same FEP stage");
02993 if (alchFepWCARepuOn && (!opts.defined("alchFepWCArcut1")||!opts.defined("alchFepWCArcut2")))
02994 NAMD_die("When using WCA repulsion, alchFepWCArcut1 and alchFepWCArcut2 must be defined!");
02995 if (alchFepWCARepuOn && (alchFepWCArcut1 >= alchFepWCArcut2))
02996 NAMD_die("When using WCA repulsion, alchFepWCArcut2 must be larger than alchFEPWCArcut1!");
02997 if ((alchFepWCARepuOn || alchFepWCADispOn) && (alchElecLambdaStart < 1.0))
02998 NAMD_die("When using WCA decomposition, repulsion, dispersion and electrostatic must be in 3 different stages!");
02999 }
03000 else if (alchThermIntOn)
03001 {
03002 if (!opts.defined("alchoutfile")) {
03003 strcpy(alchOutFile, outputFilename);
03004 strcat(alchOutFile, ".ti");
03005 }
03006 }
03007
03008 } else {
03009 alchLambda = alchLambda2 = 0;
03010 alchElecLambdaStart = 0;
03011 alchOutFile[0] = STRINGNULL;
03012 }
03013
03014
03015
03016
03017 if ( alchOn && alchFepOn && alchThermIntOn )
03018 NAMD_die("Sorry, combined TI and FEP is not implemented.\n");
03019 if ( alchOn && lesOn )
03020 NAMD_die("Sorry, combined LES with FEP or TI is not implemented.\n");
03021 if ( alchOn && alchThermIntOn && lesOn )
03022 NAMD_die("Sorry, combined LES and TI is not implemented.\n");
03023 if ( alchDecouple && (! (alchFepOn || alchThermIntOn) ) ) {
03024 iout << iWARN << "Alchemical decoupling was requested but \
03025 alchemical free energy calculation is not active. Setting \
03026 alchDecouple to off.\n" << endi;
03027 alchDecouple = FALSE;
03028
03029 }
03030
03031 if ( lesOn && ( lesFactor < 1 || lesFactor > 255 ) ) {
03032 NAMD_die("lesFactor must be positive and less than 256");
03033 }
03034 if ((pairInteractionOn && alchFepOn) || (pairInteractionOn && lesOn) || (pairInteractionOn && alchThermIntOn) )
03035 NAMD_die("Sorry, pair interactions may not be calculated when LES, FEP or TI is enabled.");
03036
03037
03038 if (drudeOn) {
03039 if ( ! langevinOn ) {
03040 NAMD_die("Drude model requires use of Langevin thermostat.");
03041 }
03042 if ( ! opts.defined("drudeDamping")) {
03043 drudeDamping = langevinDamping;
03044 iout << iWARN << "Undefined 'drudeDamping' will be set to "
03045 "value of 'langevinDamping'\n" << endi;
03046 }
03047 if ( ! opts.defined("drudeBondConst")) {
03048 drudeBondConst = 0;
03049 if (opts.defined("drudeBondLen")) {
03050 iout << iWARN << "Resetting 'drudeBondLen' to 0 "
03051 "since 'drudeBondConst' is unset\n" << endi;
03052 }
03053 drudeBondLen = 0;
03054 }
03055 if ( ! opts.defined("drudeNbtholeCut")) {
03056 drudeNbtholeCut = 0;
03057 iout << iWARN << "Resetting 'drudeNbtholeCut' to 0 "
03058 "since 'drudeNbtholeCut' is unset\n" << endi;
03059 }
03060 }
03061
03062
03063 if (opts.defined("ldBalancer")) {
03064 if (strcasecmp(loadBalancer, "none") == 0)
03065 ldBalancer = LDBAL_NONE;
03066 else if (strcasecmp(loadBalancer, "hybrid") == 0)
03067 ldBalancer = LDBAL_HYBRID;
03068 else
03069 NAMD_die("Unknown ldBalancer selected");
03070 } else {
03071 ldBalancer = LDBAL_CENTRALIZED;
03072 }
03073
03074 if (opts.defined("ldbStrategy")) {
03075
03076 if (strcasecmp(loadStrategy, "comprehensive") == 0)
03077 ldbStrategy = LDBSTRAT_COMPREHENSIVE;
03078 else if (strcasecmp(loadStrategy, "refineonly") == 0)
03079 ldbStrategy = LDBSTRAT_REFINEONLY;
03080 else if (strcasecmp(loadStrategy, "old") == 0)
03081 ldbStrategy = LDBSTRAT_OLD;
03082 else
03083 NAMD_die("Unknown ldbStrategy selected");
03084 } else {
03085 ldbStrategy = LDBSTRAT_DEFAULT;
03086 }
03087
03088 if (!opts.defined("ldbPeriod")) {
03089 ldbPeriod=200*stepsPerCycle;
03090 }
03091
03092
03093 if (!opts.defined("firstLdbStep")) {
03094 firstLdbStep=5*stepsPerCycle;
03095 }
03096
03097 if (ldbPeriod <= firstLdbStep) {
03098 NAMD_die("ldbPeriod must greater than firstLdbStep.");
03099 }
03100
03101 if (!opts.defined("lastLdbStep")) {
03102 lastLdbStep = -1;
03103 }
03104
03105 if (!opts.defined("hybridGroupSize")) {
03106 hybridGroupSize = 512;
03107 }
03108
03109
03110
03111 bool specialTracing = traceAvailable() && (traceIsOn() == 0);
03112
03113 if(!opts.defined("traceStartStep")) {
03114 traceStartStep = 4 * firstLdbStep + 2 * ldbPeriod;
03115 }
03116 if(!opts.defined("numTraceSteps")) {
03117 numTraceSteps = 100;
03118 }
03119
03120 if(specialTracing) {
03121 if (!opts.defined("firstLdbStep")) firstLdbStep = 20;
03122 if (!opts.defined("ldbPeriod")) ldbPeriod = 100;
03123
03124 if(!opts.defined("traceStartStep")) {
03125 traceStartStep = 4 * firstLdbStep + 2 * ldbPeriod;
03126 }
03127
03128 if(!opts.defined("numTraceSteps")) {
03129 numTraceSteps = 80;
03130 }
03131 }
03132
03133 #ifdef MEASURE_NAMD_WITH_PAPI
03134 if(papiMeasure){
03135 if(!opts.defined("papiMeasureStartStep")) {
03136 papiMeasureStartStep = 3 * firstLdbStep;
03137 }
03138 if(!opts.defined("numPapiMeasureSteps")) {
03139 numPapiMeasureSteps = 8;
03140 }
03141 }
03142 #endif
03143
03144 if(simulateInitialMapping) {
03145 if(!opts.defined("simulatedPEs")){
03146 simulatedPEs = CkNumPes();
03147 }
03148 if(!opts.defined("simulatedNodeSize")){
03149 simulatedNodeSize = CkMyNodeSize();
03150 }
03151 }
03152
03153 #ifdef MEM_OPT_VERSION
03154
03155
03156
03157
03158
03159
03160
03161
03162 if(ldbPeriod%stepsPerCycle!=0 || firstLdbStep%stepsPerCycle!=0) {
03163 iout << iWARN << "In memory optimized version, the ldbPeriod parameter or firstLdbStep parameter is better set to be a multiple of stepsPerCycle parameter!\n";
03164 }
03165 #endif
03166
03167 if (N < firstTimestep) { N = firstTimestep; }
03168
03169 if ( (firstTimestep%stepsPerCycle) != 0)
03170 {
03171 NAMD_die("First timestep must be a multiple of stepsPerCycle!!");
03172 }
03173
03174
03175 {
03176 int i = 0;
03177 if ( FMAOn ) ++i;
03178 if ( PMEOn ) ++i;
03179 if ( MSMOn ) ++i;
03180 if ( fullDirectOn ) ++i;
03181 if ( i > 1 )
03182 NAMD_die("More than one full electrostatics algorithm selected!!!");
03183 }
03184
03185 if (!opts.defined("ldbBackgroundScaling")) {
03186 ldbBackgroundScaling = 1.0;
03187 }
03188 if (!opts.defined("ldbPMEBackgroundScaling")) {
03189 ldbPMEBackgroundScaling = ldbBackgroundScaling;
03190 }
03191 if (!opts.defined("ldbHomeBackgroundScaling")) {
03192 ldbHomeBackgroundScaling = ldbBackgroundScaling;
03193 }
03194
03195
03196 if (PMEOn) {
03197 if ( lattice.volume() == 0. ) {
03198 NAMD_die("PME requires periodic boundary conditions.");
03199 }
03200 if ( PMEGridSpacing == 0. ) {
03201 if ( PMEGridSizeX * PMEGridSizeY * PMEGridSizeZ == 0 )
03202 NAMD_die("Either PMEGridSpacing or PMEGridSizeX, PMEGridSizeY, and PMEGridSizeZ must be specified.");
03203 else PMEGridSpacing = 1.5;
03204 }
03205 #ifndef TEST_PME_GRID
03206 for ( int idim = 0; idim < 3; ++idim ) {
03207 int *gridSize;
03208 BigReal cellLength;
03209 const char *direction;
03210 switch ( idim ) {
03211 case 0: direction = "X";
03212 gridSize = &PMEGridSizeX; cellLength = lattice.a().length();
03213 break;
03214 case 1: direction = "Y";
03215 gridSize = &PMEGridSizeY; cellLength = lattice.b().length();
03216 break;
03217 case 2: direction = "Z";
03218 gridSize = &PMEGridSizeZ; cellLength = lattice.c().length();
03219 break;
03220 }
03221 int minSize = (int) ceil(cellLength/PMEGridSpacing);
03222 #else
03223 for ( int minSize = 1; minSize < 300; ++minSize ) {
03224 #endif
03225 int bestSize = 10 * (minSize + 10);
03226 int max2, max3, ts;
03227 for ( max2=2, ts=1; ts < minSize; ++max2 ) ts *= 2;
03228 for ( max3=2, ts=1; ts < minSize; ++max3 ) ts *= 3;
03229 int max5 = 2;
03230 int max7 = 1;
03231 int max11 = 1;
03232 for ( int i2 = 0; i2 <= max2; ++i2 ) {
03233 for ( int i3 = 0; i3 <= max3; ++i3 ) {
03234 for ( int i5 = 0; i5 <= max5; ++i5 ) {
03235 for ( int i7 = 0; i7 <= max7; ++i7 ) {
03236 for ( int i11 = 0; i11 <= max11; ++i11 ) {
03237 if ( i5 + i7 + i11 > i2 ) continue;
03238 int testSize = 2;
03239 for ( int j2 = 0; j2 < i2; ++j2 ) testSize *= 2;
03240 if ( testSize > bestSize ) continue;
03241 for ( int j3 = 0; j3 < i3; ++j3 ) testSize *= 3;
03242 if ( testSize > bestSize ) continue;
03243 for ( int j5 = 0; j5 < i5; ++j5 ) testSize *= 5;
03244 if ( testSize > bestSize ) continue;
03245 for ( int j7 = 0; j7 < i7; ++j7 ) testSize *= 7;
03246 if ( testSize > bestSize ) continue;
03247 for ( int j11 = 0; j11 < i11; ++j11 ) testSize *= 11;
03248 if ( testSize > bestSize ) continue;
03249 if ( testSize >= minSize ) bestSize = testSize;
03250 } } } } }
03251 #ifdef TEST_PME_GRID
03252 iout << minSize << " " << bestSize << "\n" << endi;
03253 #else
03254 if ( ! *gridSize ) {
03255 *gridSize = bestSize;
03256 }
03257 if ( *gridSize * PMEGridSpacing < cellLength ) {
03258 char errmsg[512];
03259 sprintf(errmsg, "PMEGridSize%s %d is too small for cell length %f and PMEGridSpacing %f\n",
03260 direction, *gridSize, cellLength, PMEGridSpacing);
03261 NAMD_die(errmsg);
03262 }
03263 #endif
03264 }
03265 if ( PMEGridSizeX < 5 ) {
03266 NAMD_die("PMEGridSizeX (number of grid points) is very small.");
03267 }
03268 if ( PMEGridSizeY < 5 ) {
03269 NAMD_die("PMEGridSizeY (number of grid points) is very small.");
03270 }
03271 if ( PMEGridSizeZ < 5 ) {
03272 NAMD_die("PMEGridSizeZ (number of grid points) is very small.");
03273 }
03274 BigReal tolerance = PMETolerance;
03275 BigReal ewaldcof = 1.0;
03276 while ( erfc(ewaldcof*cutoff)/cutoff >= tolerance ) ewaldcof *= 2.0;
03277 BigReal ewaldcof_lo = 0.;
03278 BigReal ewaldcof_hi = ewaldcof;
03279 for ( int i = 0; i < 100; ++i ) {
03280 ewaldcof = 0.5 * ( ewaldcof_lo + ewaldcof_hi );
03281 if ( erfc(ewaldcof*cutoff)/cutoff >= tolerance ) {
03282 ewaldcof_lo = ewaldcof;
03283 } else {
03284 ewaldcof_hi = ewaldcof;
03285 }
03286 }
03287 PMEEwaldCoefficient = ewaldcof;
03288 } else {
03289 useDPME = 0;
03290 PMEGridSizeX = 0;
03291 PMEGridSizeY = 0;
03292 PMEGridSizeZ = 0;
03293 PMEGridSpacing = 1000.;
03294 PMEEwaldCoefficient = 0;
03295 }
03296
03297
03298
03299
03300 if (!FMAOn)
03301 {
03302 FMALevels = 0;
03303 FMAMp = 0;
03304 FMAFFTOn = FALSE;
03305 FMAFFTBlock = 0;
03306 }
03307 else
03308 {
03309
03310
03311 if (FMAFFTBlock != 4)
03312 NAMD_die("FMAFFTBlock: Block length must be 4 for short FFT's");
03313 if (FMAMp % FMAFFTBlock != 0)
03314 NAMD_die("FMAMp: multipole term must be multiple of block length (FMAFFTBlock)");
03315 }
03316
03317 if ( (nonbondedFrequency > stepsPerCycle) || ( (stepsPerCycle % nonbondedFrequency) != 0) )
03318 {
03319 NAMD_die("stepsPerCycle must be a multiple of nonbondedFreq");
03320 }
03321
03322 if (!LCPOOn && !GBISOn && !GBISserOn && !FMAOn && !PMEOn && !MSMOn && !fullDirectOn)
03323 {
03324 fullElectFrequency = 0;
03325 }
03326 else
03327 {
03328 if (!opts.defined("fullElectFrequency"))
03329 {
03330 if (opts.defined("fmaFrequency")) {
03331 iout << iWARN << "The parameter fmaFrequency has been renamed fullElectFrequency.\n" << endi;
03332 fullElectFrequency = fmaFrequency;
03333 } else {
03334 iout << iWARN << "The parameter fullElectFrequency now defaults to nonbondedFreq (" << nonbondedFrequency << ") rather than stepsPerCycle.\n" << endi;
03335 fullElectFrequency = nonbondedFrequency;
03336 }
03337 }
03338 else
03339 {
03340 if (opts.defined("fmaFrequency")) {
03341 iout << iWARN << "Ignoring redundant parameter fmaFrequency in favor of fullElectFrequency.\n" << endi;
03342 }
03343 if ( (fullElectFrequency > stepsPerCycle) || ( (stepsPerCycle % fullElectFrequency) != 0) )
03344 {
03345 NAMD_die("stepsPerCycle must be a multiple of fullElectFrequency");
03346 }
03347 }
03348
03349 if ( (nonbondedFrequency > fullElectFrequency) || ( (fullElectFrequency % nonbondedFrequency) != 0) )
03350 {
03351 NAMD_die("fullElectFrequency must be a multiple of nonbondedFreq");
03352 }
03353
03354
03355 if (!opts.defined("fmaTheta"))
03356 fmaTheta=0.715;
03357 }
03358
03359 if ( lesOn && ( FMAOn || useDPME || fullDirectOn ) ) {
03360 NAMD_die("Sorry, LES is only implemented for PME full electrostatics.");
03361 }
03362 if ( alchFepOn && ( FMAOn || useDPME || fullDirectOn ) ) {
03363 NAMD_die("Sorry, FEP is only implemented for PME full electrostatics.");
03364 }
03365 if ( alchThermIntOn && ( FMAOn || useDPME || fullDirectOn ) ) {
03366 NAMD_die("Sorry, TI is only implemented for PME full electrostatics.");
03367 }
03368 if ( pairInteractionOn && FMAOn ) {
03369 NAMD_die("Sorry, pairInteraction not implemented for FMA.");
03370 }
03371 if ( pairInteractionOn && useDPME ) {
03372 NAMD_die("Sorry, pairInteraction not implemented for DPME.");
03373 }
03374 if ( pairInteractionOn && fullDirectOn ) {
03375 NAMD_die("Sorry, pairInteraction not implemented for full direct electrostatics.");
03376 }
03377 if ( ! pairInteractionOn ) {
03378 pairInteractionSelf = 0;
03379 }
03380 if ( pairInteractionOn && !pairInteractionSelf && !config->find("pairInteractionGroup2"))
03381 NAMD_die("pairInteractionGroup2 must be specified");
03382
03383 if ( ! fixedAtomsOn ) {
03384 fixedAtomsForces = 0;
03385 }
03386
03387 if ( gridforceOn || mgridforceOn ) {
03388 parse_mgrid_params(config);
03389 }
03390
03391 if (!opts.defined("constraints"))
03392 {
03393 constraintExp = 0;
03394 constraintScaling = 1.0;
03395
03396
03397 selectConstraintsOn = FALSE;
03398
03399
03400
03401 movingConstraintsOn = FALSE;
03402
03403
03404 rotConstraintsOn = FALSE;
03405
03406 }
03407
03408 else {
03409 if (rotConstraintsOn) {
03410 rotConsAxis = rotConsAxis.unit();
03411 }
03412 }
03413 if(opts.defined("rotConstraints")
03414 && opts.defined("movingConstraints")) {
03415 NAMD_die("Rotating and moving constraints are mutually exclusive!");
03416 }
03417
03418
03419
03420 if(opts.defined("selectConstraints") && !opts.defined("selectConstrX")
03421 && !opts.defined("selectConstrY") && !opts.defined("selectConstrZ")) {
03422 NAMD_die("selectConstraints was specified, but no Cartesian components were defined!");
03423 }
03424 if (!opts.defined("selectConstraints")) {
03425 constrXOn = FALSE;
03426 constrYOn = FALSE;
03427 constrZOn = FALSE;
03428 }
03429
03430
03431
03432
03433
03434 if (!opts.defined("SMD")) {
03435 SMDOn = FALSE;
03436 }
03437
03438 if (SMDOn) {
03439
03440 if (SMDDir.length2() == 0) {
03441 NAMD_die("SMD direction vector must be non-zero");
03442 }
03443 else {
03444 SMDDir = SMDDir.unit();
03445 }
03446
03447 if (SMDOutputFreq > 0 && SMDOutputFreq < stepsPerCycle
03448 || SMDOutputFreq % stepsPerCycle != 0) {
03449 NAMD_die("SMDOutputFreq must be a multiple of stepsPerCycle");
03450 }
03451 }
03452
03453
03454
03455 if (!sphericalBCOn)
03456 {
03457 sphericalBCr1 = 0.0;
03458 sphericalBCk1 = 0.0;
03459 sphericalBCexp1 = 0;
03460 sphericalBCr2 = 0.0;
03461 sphericalBCk2 = 0.0;
03462 sphericalBCexp2 = 0;
03463 }
03464 else if (!opts.defined("sphericalBCr2"))
03465 {
03466 sphericalBCr2 = -1.0;
03467 sphericalBCk2 = 0.0;
03468 sphericalBCexp2 = 0;
03469 }
03470
03471 if (!cylindricalBCOn)
03472 {
03473 cylindricalBCr1 = 0.0;
03474 cylindricalBCk1 = 0.0;
03475 cylindricalBCexp1 = 0;
03476 cylindricalBCr2 = 0.0;
03477 cylindricalBCk2 = 0.0;
03478 cylindricalBCexp2 = 0;
03479 cylindricalBCl1 = 0.0;
03480 cylindricalBCl2 = 0.0;
03481 }
03482 else if (!opts.defined("cylindricalBCr2"))
03483 {
03484 cylindricalBCr2 = -1.0;
03485 cylindricalBCk2 = 0.0;
03486 cylindricalBCexp2 = 0;
03487 cylindricalBCl2 = 0.0;
03488 }
03489
03490 if (!eFieldOn)
03491 {
03492 eField.x = 0.0;
03493 eField.y = 0.0;
03494 eField.z = 0.0;
03495 eFieldFreq = 0.0;
03496 eFieldPhase = 0.0;
03497 }
03498 else
03499 {
03500 if (!opts.defined("eFieldFreq")) eFieldFreq = 0.0;
03501 if (!opts.defined("eFieldPhase")) eFieldPhase = 0.0;
03502 }
03503
03504 if (!stirOn)
03505 {
03506 stirFilename[0] = STRINGNULL;
03507 stirStartingTheta = 0.0;
03508 stirVel = 0.0;
03509 stirK = 0.0;
03510 stirAxis.x = 0.0;
03511 stirAxis.y = 0.0;
03512 stirAxis.z = 0.0;
03513 stirPivot.x = 0.0;
03514 stirPivot.y = 0.0;
03515 stirPivot.z = 0.0;
03516 }
03517
03518 if (!opts.defined("langevin"))
03519 {
03520 langevinTemp = 0.0;
03521 }
03522
03523
03524 if (!opts.defined("loweAndersen"))
03525 {
03526 loweAndersenTemp = 0.0;
03527 }
03528
03529
03530 if (!opts.defined("tcouple"))
03531 {
03532 tCoupleTemp = 0.0;
03533 }
03534
03535 if (HydrogenBonds)
03536 {
03537 if (daCutoffDist > pairlistDist)
03538 NAMD_die("Hydrogen bond cutoff distance must be <= pairlist distance");
03539 }
03540
03541
03542
03543
03544
03545 if (pairInteractionOn) {
03546 if (outputEnergies != 1) {
03547 iout << iWARN << "Setting outputEnergies to 1 due to\n";
03548 iout << iWARN << "pairInteraction calculations\n" << endi;
03549 outputEnergies = 1;
03550 }
03551 }
03552 if (pairInteractionOn || pressureProfileOn) {
03553 if (nonbondedFrequency != 1) {
03554 iout << iWARN << "Setting nonbondedFreq to 1 due to\n";
03555 iout << iWARN << "pairInteraction or pressure profile calculations\n" << endi;
03556 }
03557 }
03558
03559
03560 if (!opts.defined("outputTiming"))
03561 {
03562 outputTiming = firstLdbStep;
03563 int ot2 = 10 * outputEnergies;
03564 if ( outputTiming < ot2 ) outputTiming = ot2;
03565 }
03566
03567 }
03568
03569 void SimParameters::print_config(ParseOptions &opts, ConfigList *config, char *&cwd) {
03570
03571 StringList *current;
03572
03573
03574
03575 iout << iINFO << "SIMULATION PARAMETERS:\n";
03576 iout << iINFO << "TIMESTEP " << dt << "\n" << endi;
03577 iout << iINFO << "NUMBER OF STEPS " << N << "\n";
03578 iout << iINFO << "STEPS PER CYCLE " << stepsPerCycle << "\n";
03579 iout << endi;
03580
03581 if ( lattice.a_p() || lattice.b_p() || lattice.c_p() ) {
03582 if ( lattice.a_p() )
03583 iout << iINFO << "PERIODIC CELL BASIS 1 " << lattice.a() << "\n";
03584 if ( lattice.b_p() )
03585 iout << iINFO << "PERIODIC CELL BASIS 2 " << lattice.b() << "\n";
03586 if ( lattice.c_p() )
03587 iout << iINFO << "PERIODIC CELL BASIS 3 " << lattice.c() << "\n";
03588 iout << iINFO << "PERIODIC CELL CENTER " << lattice.origin() << "\n";
03589 if (wrapWater) {
03590 iout << iINFO << "WRAPPING WATERS AROUND PERIODIC BOUNDARIES ON OUTPUT.\n";
03591 }
03592 if (wrapAll) {
03593 iout << iINFO << "WRAPPING ALL CLUSTERS AROUND PERIODIC BOUNDARIES ON OUTPUT.\n";
03594 }
03595 if (wrapNearest) {
03596 iout << iINFO << "WRAPPING TO IMAGE NEAREST TO PERIODIC CELL CENTER.\n";
03597 }
03598 iout << endi;
03599 }
03600
03601 if ( CkNumPes() > 512 ) ldbUnloadOne = TRUE;
03602 if ( ldbUnloadOne || CkNumPes() > 128 ) ldbUnloadZero = TRUE;
03603
03604 if (ldBalancer == LDBAL_NONE) {
03605 iout << iINFO << "LOAD BALANCER None\n" << endi;
03606 } else {
03607 if (ldBalancer == LDBAL_CENTRALIZED) {
03608 iout << iINFO << "LOAD BALANCER Centralized\n" << endi;
03609 } else if (ldBalancer == LDBAL_HYBRID) {
03610 iout << iINFO << "LOAD BALANCER Hybrid\n" << endi;
03611 }
03612
03613 if (ldbStrategy == LDBSTRAT_DEFAULT) {
03614 iout << iINFO << "LOAD BALANCING STRATEGY New Load Balancers -- DEFAULT\n";
03615 } else if (ldbStrategy == LDBSTRAT_REFINEONLY) {
03616 iout << iINFO << "LOAD BALANCING STRATEGY Refinement Only\n";
03617 } else if (ldbStrategy == LDBSTRAT_COMPREHENSIVE) {
03618 iout << iINFO << "LOAD BALANCING STRATEGY Comprehensive\n";
03619 } else if (ldbStrategy == LDBSTRAT_OLD) {
03620 iout << iINFO << "LOAD BALANCING STRATEGY Old Load Balancers\n";
03621 }
03622
03623 iout << iINFO << "LDB PERIOD " << ldbPeriod << " steps\n";
03624 iout << iINFO << "FIRST LDB TIMESTEP " << firstLdbStep << "\n";
03625 if (ldBalancer == LDBAL_HYBRID)
03626 iout << iINFO << "HYBRIDLB GROUP SIZE " << hybridGroupSize << "\n";
03627 iout << iINFO << "LAST LDB TIMESTEP " << lastLdbStep << "\n";
03628 if ( ldbRelativeGrainsize > 0. )
03629 iout << iINFO << "LDB RELATIVE GRAINSIZE " << ldbRelativeGrainsize << "\n";
03630 iout << iINFO << "LDB BACKGROUND SCALING " << ldbBackgroundScaling << "\n";
03631 iout << iINFO << "HOM BACKGROUND SCALING " << ldbHomeBackgroundScaling << "\n";
03632 if ( PMEOn ) {
03633 iout << iINFO << "PME BACKGROUND SCALING "
03634 << ldbPMEBackgroundScaling << "\n";
03635 if ( ldbUnloadPME )
03636 iout << iINFO << "REMOVING LOAD FROM PME NODES" << "\n";
03637 }
03638 if ( ldbUnloadZero ) iout << iINFO << "REMOVING LOAD FROM NODE 0\n";
03639 if ( ldbUnloadOne ) iout << iINFO << "REMOVING LOAD FROM NODE 1\n";
03640 iout << endi;
03641 }
03642
03643 if ( ldbUnloadOne || CkNumPes() > 256 ) noPatchesOnOne = TRUE;
03644 if ( ldbUnloadZero || noPatchesOnOne ||
03645 CkNumPes() > 64 || ( IMDon && CkNumPes() > 8 ) ) {
03646 noPatchesOnZero = TRUE;
03647 }
03648 #ifdef NAMD_CUDA
03649 noPatchesOnZero = FALSE;
03650 noPatchesOnOne = FALSE;
03651 #endif
03652 if ( noPatchesOnZero ) iout << iINFO << "REMOVING PATCHES FROM PROCESSOR 0\n";
03653 if ( noPatchesOnOne ) iout << iINFO << "REMOVING PATCHES FROM PROCESSOR 1\n";
03654 iout << endi;
03655
03656 #ifdef NAMD_CUDA
03657 maxSelfPart = maxPairPart = 1;
03658 #endif
03659
03660 if (ldBalancer == LDBAL_HYBRID) {
03661 iout << iINFO << "MAX SELF PARTITIONS " << maxSelfPart << "\n"
03662 << iINFO << "MAX PAIR PARTITIONS " << maxPairPart << "\n"
03663 << iINFO << "SELF PARTITION ATOMS " << numAtomsSelf << "\n"
03664 << iINFO << "SELF2 PARTITION ATOMS " << numAtomsSelf2 << "\n"
03665 << iINFO << "PAIR PARTITION ATOMS " << numAtomsPair << "\n"
03666 << iINFO << "PAIR2 PARTITION ATOMS " << numAtomsPair2 << "\n";
03667 }
03668 iout << iINFO << "MIN ATOMS PER PATCH " << minAtomsPerPatch << "\n"
03669 << endi;
03670
03671 if (initialTemp < 0)
03672 {
03673 current = config->find("velocities");
03674
03675 if (current == NULL)
03676 {
03677 current = config->find("binvelocities");
03678 }
03679
03680 iout << iINFO << "VELOCITY FILE " << current->data << "\n";
03681 }
03682 else
03683 {
03684 iout << iINFO << "INITIAL TEMPERATURE "
03685 << initialTemp << "\n";
03686 }
03687 iout << endi;
03688
03689 iout << iINFO << "CENTER OF MASS MOVING INITIALLY? ";
03690
03691 if (comMove)
03692 {
03693 iout << "YES\n";
03694 }
03695 else
03696 {
03697 iout << "NO\n";
03698 }
03699 iout << endi;
03700
03701 if ( zeroMomentum ) {
03702 iout << iINFO << "REMOVING CENTER OF MASS DRIFT DURING SIMULATION";
03703 if ( zeroMomentumAlt ) iout << " (ALT METHOD)";
03704 iout << "\n" << endi;
03705 }
03706
03707 iout << iINFO << "DIELECTRIC "
03708 << dielectric << "\n";
03709
03710 if ( nonbondedScaling != 1.0 )
03711 {
03712 iout << iINFO << "NONBONDED SCALING " << nonbondedScaling << "\n" << endi;
03713 }
03714 iout << iINFO << "EXCLUDE ";
03715
03716 switch (exclude)
03717 {
03718 case NONE:
03719 iout << "NONE\n";
03720 break;
03721 case ONETWO:
03722 iout << "ONETWO\n";
03723 break;
03724 case ONETHREE:
03725 iout << "ONETHREE\n";
03726 break;
03727 case ONEFOUR:
03728 iout << "ONE-FOUR\n";
03729 break;
03730 default:
03731 iout << "SCALED ONE-FOUR\n";
03732 break;
03733 }
03734 iout << endi;
03735
03736 if (exclude == SCALED14)
03737 {
03738 iout << iINFO << "1-4 ELECTROSTATICS SCALED BY " << scale14 << "\n";
03739 iout << iINFO << "MODIFIED 1-4 VDW PARAMETERS WILL BE USED\n" << endi;
03740 } else {
03741 iout << iWARN << "MODIFIED 1-4 VDW PARAMETERS WILL BE IGNORED\n" << endi;
03742 }
03743
03744 if (dcdFrequency > 0)
03745 {
03746 iout << iINFO << "DCD FILENAME "
03747 << dcdFilename << "\n";
03748 iout << iINFO << "DCD FREQUENCY "
03749 << dcdFrequency << "\n";
03750 iout << iINFO << "DCD FIRST STEP "
03751 << ( ((firstTimestep + dcdFrequency)/dcdFrequency)*dcdFrequency ) << "\n";
03752 if ( dcdUnitCell ) {
03753 iout << iINFO << "DCD FILE WILL CONTAIN UNIT CELL DATA\n";
03754 }
03755 }
03756 else
03757 {
03758 iout << iINFO << "NO DCD TRAJECTORY OUTPUT\n";
03759 }
03760 iout << endi;
03761
03762 if (xstFrequency > 0)
03763 {
03764 iout << iINFO << "XST FILENAME "
03765 << xstFilename << "\n";
03766 iout << iINFO << "XST FREQUENCY "
03767 << xstFrequency << "\n";
03768 }
03769 else
03770 {
03771 iout << iINFO << "NO EXTENDED SYSTEM TRAJECTORY OUTPUT\n";
03772 }
03773 iout << endi;
03774
03775 if (velDcdFrequency > 0)
03776 {
03777 iout << iINFO << "VELOCITY DCD FILENAME "
03778 << velDcdFilename << "\n";
03779 iout << iINFO << "VELOCITY DCD FREQUENCY "
03780 << velDcdFrequency << "\n";
03781 iout << iINFO << "VELOCITY DCD FIRST STEP "
03782 << ( ((firstTimestep + velDcdFrequency)/velDcdFrequency)*velDcdFrequency ) << "\n";
03783 }
03784 else
03785 {
03786 iout << iINFO << "NO VELOCITY DCD OUTPUT\n";
03787 }
03788 iout << endi;
03789
03790 if (forceDcdFrequency > 0)
03791 {
03792 iout << iINFO << "FORCE DCD FILENAME "
03793 << forceDcdFilename << "\n";
03794 iout << iINFO << "FORCE DCD FREQUENCY "
03795 << forceDcdFrequency << "\n";
03796 iout << iINFO << "FORCE DCD FIRST STEP "
03797 << ( ((firstTimestep + forceDcdFrequency)/forceDcdFrequency)*forceDcdFrequency ) << "\n";
03798 }
03799 else
03800 {
03801 iout << iINFO << "NO FORCE DCD OUTPUT\n";
03802 }
03803 iout << endi;
03804
03805 iout << iINFO << "OUTPUT FILENAME "
03806 << outputFilename << "\n" << endi;
03807 if (binaryOutput)
03808 {
03809 iout << iINFO << "BINARY OUTPUT FILES WILL BE USED\n" << endi;
03810 }
03811 #ifdef MEM_OPT_VERSION
03812 if(!binaryOutput){
03813 iout << iWARN <<"SINCE MEMORY OPTIMIZED VERSION IS USED, OUTPUT IN TEXT FORMAT IS DISABLED!\n" << endi;
03814 binaryOutput = TRUE;
03815 }
03816 #endif
03817
03818 if (! restartFrequency)
03819 {
03820 iout << iINFO << "NO RESTART FILE\n";
03821 }
03822 else
03823 {
03824 iout << iINFO << "RESTART FILENAME "
03825 << restartFilename << "\n";
03826 iout << iINFO << "RESTART FREQUENCY "
03827 << restartFrequency << "\n";
03828 if (restartSave) {
03829 iout << iINFO << "RESTART FILES WILL NOT BE OVERWRITTEN\n";
03830 }
03831
03832 if (binaryRestart)
03833 {
03834 iout << iINFO << "BINARY RESTART FILES WILL BE USED\n";
03835 }
03836 }
03837 iout << endi;
03838
03839 if (switchingActive)
03840 {
03841 iout << iINFO << "SWITCHING ACTIVE\n";
03842 if ( vdwForceSwitching ) {
03843 iout << iINFO << "VDW FORCE SWITCHING ACTIVE\n";
03844 }
03845 if ( martiniSwitching ) {
03846 iout << iINFO << "MARTINI RESIDUE-BASED COARSE-GRAIN SWITCHING ACTIVE\n";
03847 }
03848 iout << iINFO << "SWITCHING ON "
03849 << switchingDist << "\n";
03850 iout << iINFO << "SWITCHING OFF "
03851 << cutoff << "\n";
03852 }
03853 else
03854 {
03855 iout << iINFO << "CUTOFF "
03856 << cutoff << "\n";
03857 }
03858
03859 iout << iINFO << "PAIRLIST DISTANCE " << pairlistDist << "\n";
03860 iout << iINFO << "PAIRLIST SHRINK RATE " << pairlistShrink << "\n";
03861 iout << iINFO << "PAIRLIST GROW RATE " << pairlistGrow << "\n";
03862 iout << iINFO << "PAIRLIST TRIGGER " << pairlistTrigger << "\n";
03863 iout << iINFO << "PAIRLISTS PER CYCLE " << pairlistsPerCycle << "\n";
03864 if ( outputPairlists )
03865 iout << iINFO << "PAIRLIST OUTPUT STEPS " << outputPairlists << "\n";
03866 iout << endi;
03867
03868 if ( pairlistMinProcs > 1 )
03869 iout << iINFO << "REQUIRING " << pairlistMinProcs << " PROCESSORS FOR PAIRLISTS\n";
03870 usePairlists = ( CkNumPes() >= pairlistMinProcs );
03871
03872 #ifdef OPENATOM_VERSION
03873 if ( openatomOn )
03874 {
03875 iout << iINFO << "OPENATOM QM/MM CAR-PARINELLO ACTIVE\n";
03876 iout << iINFO << "OPENATOM CONFIG FILE: " << openatomConfig << "\n";
03877 iout << iINFO << "OPENATOM STRUCT FILE: " << openatomStruct << "\n";
03878 iout << iINFO << "OPENATOM PDB FILE: " << openatomPDB << "\n";
03879 }
03880 #endif // OPENATOM_VERSION
03881
03882
03883 if ( (alchOn) && (!usePairlists)) {
03884 NAMD_die("Sorry, Alchemical simulations require pairlists to be enabled\n");
03885 }
03886 #ifdef NAMD_CUDA
03887 if ( ! usePairlists ) {
03888 usePairlists = 1;
03889 iout << iINFO << "CUDA ACCELERATION REQUIRES PAIRLISTS\n";
03890 }
03891 #endif
03892
03893 iout << iINFO << "PAIRLISTS " << ( usePairlists ? "ENABLED" : "DISABLED" )
03894 << "\n" << endi;
03895
03896 iout << iINFO << "MARGIN " << margin << "\n";
03897
03898 if ( splitPatch == SPLIT_PATCH_HYDROGEN ) {
03899 iout << iINFO << "HYDROGEN GROUP CUTOFF " << hgroupCutoff << "\n";
03900 }
03901
03902 iout << iINFO << "PATCH DIMENSION "
03903 << patchDimension << "\n";
03904
03905 iout << endi;
03906
03907 if (outputEnergies != 1)
03908 {
03909 iout << iINFO << "ENERGY OUTPUT STEPS "
03910 << outputEnergies << "\n";
03911 iout << endi;
03912 }
03913
03914 if (mergeCrossterms) {
03915 iout << iINFO << "CROSSTERM ENERGY INCLUDED IN DIHEDRAL\n" << endi;
03916 }
03917
03918 if (outputMomenta != 0)
03919 {
03920 iout << iINFO << "MOMENTUM OUTPUT STEPS "
03921 << outputMomenta << "\n";
03922 iout << endi;
03923 }
03924
03925 if (outputTiming != 0)
03926 {
03927 iout << iINFO << "TIMING OUTPUT STEPS "
03928 << outputTiming << "\n";
03929 iout << endi;
03930 }
03931
03932 if (outputCudaTiming != 0)
03933 {
03934 iout << iINFO << "CUDA TIMING OUTPUT STEPS "
03935 << outputCudaTiming << "\n";
03936 iout << endi;
03937 }
03938
03939 if (outputPressure != 0)
03940 {
03941 iout << iINFO << "PRESSURE OUTPUT STEPS "
03942 << outputPressure << "\n";
03943 iout << endi;
03944 }
03945
03946 if (fixedAtomsOn)
03947 {
03948 iout << iINFO << "FIXED ATOMS ACTIVE\n";
03949 #ifdef NAMD_CUDA
03950 iout << iWARN << "FOR CUDA VERSION PRESSURE WILL BE INCORRECT WITHOUT FIXEDATOMSFORCES\n";
03951 if (berendsenPressureOn || langevinPistonOn) {
03952 fixedAtomsForces = TRUE;
03953 iout << iWARN << "ENABLING FIXEDATOMSFORCES DUE TO PRESSURE CONTROL\n";
03954 }
03955 #endif
03956 if ( fixedAtomsForces )
03957 iout << iINFO << "FORCES BETWEEN FIXED ATOMS ARE CALCULATED\n";
03958 iout << endi;
03959 }
03960
03961 if (constraintsOn)
03962 {
03963 iout << iINFO << "HARMONIC CONSTRAINTS ACTIVE\n";
03964
03965 iout << iINFO << "HARMONIC CONS EXP "
03966 << constraintExp << "\n";
03967
03968 if (constraintScaling != 1.0) {
03969 iout << iINFO << "HARMONIC CONS SCALING "
03970 << constraintScaling << "\n";
03971 }
03972
03973
03974
03975 if (selectConstraintsOn) {
03976 iout << iINFO << "SELECTED CARTESIAN COMPONENTS OF HARMONIC RESTRAINTS ACTIVE\n";
03977
03978 if (constrXOn)
03979 iout << iINFO << "RESTRAINING X-COMPONENTS OF CARTESIAN COORDINATES!\n";
03980
03981 if (constrYOn)
03982 iout << iINFO << "RESTRAINING Y-COMPONENTS OF CARTESIAN COORDINATES!\n";
03983
03984 if (constrZOn)
03985 iout << iINFO << "RESTRAINING Z-COMPONENTS OF CARTESIAN COORDINATES!\n";
03986 }
03987 iout << endi;
03988
03989
03990
03991
03992 if (movingConstraintsOn) {
03993 iout << iINFO << "MOVING HARMONIC CONSTRAINTS ACTIVE\n";
03994
03995 iout << iINFO << "MOVING CONSTRAINT VELOCITY "
03996 << movingConsVel << " ANGSTROM/TIMESTEP\n";
03997
03998 iout << iINFO << "ALL CONSTRAINED ATOMS WILL MOVE\n";
03999 }
04000
04001 iout << endi;
04002
04003
04004
04005 if (rotConstraintsOn) {
04006 iout << iINFO << "ROTATING HARMONIC CONSTRAINTS ACTIVE\n";
04007
04008 iout << iINFO << "AXIS OF ROTATION "
04009 << rotConsAxis << "\n";
04010
04011 iout << iINFO << "PIVOT OF ROTATION "
04012 << rotConsPivot << "\n";
04013
04014 iout << iINFO << "ROTATING CONSTRAINT VELOCITY "
04015 << rotConsVel << " DEGREES/TIMESTEP\n";
04016 }
04017 iout << endi;
04018
04019 }
04020
04021
04022 if (movDragOn) {
04023 iout << iINFO << "MOVING DRAG ACTIVE.\n";
04024
04025 iout << iINFO << "MOVING DRAG MAIN PDB FILE "
04026 << movDragFile << "\n";
04027
04028 iout << iINFO << "MOVING DRAG GLOBAL VELOCITY (A/step) "
04029 << movDragGlobVel << "\n";
04030
04031 iout << iINFO << "MOVING DRAG LINEAR VELOCITY FILE "
04032 << movDragVelFile << "\n";
04033
04034 iout << endi;
04035 }
04036
04037
04038 if (rotDragOn) {
04039 iout << iINFO << "ROTATING DRAG ACTIVE.\n";
04040
04041 iout << iINFO << "ROTATING DRAG MAIN PDB FILE "
04042 << rotDragFile << "\n";
04043
04044 iout << iINFO << "ROTATING DRAG AXIS FILE "
04045 << rotDragAxisFile << "\n";
04046
04047 iout << iINFO << "ROTATING DRAG PIVOT POINT FILE "
04048 << rotDragPivotFile << "\n";
04049
04050 iout << iINFO << "ROTATING DRAG GLOBAL ANGULAR VELOCITY (deg/step) "
04051 << rotDragGlobVel << "\n";
04052
04053 iout << iINFO << "ROTATING DRAG ANGULAR VELOCITY FILE "
04054 << rotDragVelFile << "\n";
04055
04056 iout << endi;
04057 }
04058
04059
04060
04061 if (consTorqueOn) {
04062 iout << iINFO << "\"CONSTANT\" TORQUE ACTIVE.\n";
04063
04064 iout << iINFO << "\"CONSTANT\" TORQUE MAIN PDB FILE "
04065 << consTorqueFile << "\n";
04066
04067 iout << iINFO << "\"CONSTANT\" TORQUE AXIS FILE "
04068 << consTorqueAxisFile << "\n";
04069
04070 iout << iINFO << "\"CONSTANT\" TORQUE PIVOT POINT FILE "
04071 << consTorquePivotFile << "\n";
04072
04073 iout << iINFO << "\"CONSTANT\" TORQUE GLOBAL VALUE (Kcal/(mol*A^2)) "
04074 << consTorqueGlobVal << "\n";
04075
04076 iout << iINFO << "\"CONSTANT\" TORQUE DACTORS FILE "
04077 << consTorqueValFile << "\n";
04078
04079 iout << endi;
04080 }
04081
04082 if (mgridforceOn) {
04083 iout << iINFO << "GRID FORCE ACTIVE\n";
04084 iout << iINFO << " Please include this reference in published work using\n";
04085 iout << iINFO << " the Gridforce module of NAMD: David Wells, Volha Abramkina,\n";
04086 iout << iINFO << " and Aleksei Aksimentiev, J. Chem. Phys. 127:125101-10 (2007).\n";
04087 print_mgrid_params();
04088 }
04089
04090
04091
04092 if (SMDOn) {
04093 iout << iINFO << "SMD ACTIVE\n";
04094
04095 iout << iINFO << "SMD VELOCITY "
04096 << SMDVel << " ANGSTROM/TIMESTEP\n";
04097
04098 iout << iINFO << "SMD DIRECTION "
04099 << SMDDir << "\n";
04100
04101 iout << iINFO << "SMD K "
04102 << SMDk << "\n";
04103
04104 iout << iINFO << "SMD K2 "
04105 << SMDk2 << "\n";
04106
04107 iout << iINFO << "SMD OUTPUT FREQUENCY "
04108 << SMDOutputFreq << " TIMESTEPS\n";
04109
04110 iout << iINFO << "SMD FILE " << SMDFile << "\n";
04111
04112 iout << endi;
04113 }
04114
04115
04116
04117 if (TMDOn) {
04118 iout << iINFO << "TMD ACTIVE BETWEEN STEPS " << TMDFirstStep
04119 << " and " << TMDLastStep << "\n";
04120 iout << iINFO << "TMD K " << TMDk << "\n";
04121 iout << iINFO << "TMD FILE " << TMDFile << "\n";
04122 iout << iINFO << "TMD OUTPUT FREQUENCY " << TMDOutputFreq << "\n";
04123 if (TMDInitialRMSD) {
04124 iout << iINFO << "TMD TARGET RMSD AT FIRST STEP " << TMDInitialRMSD << "\n";
04125 } else {
04126 iout << iINFO << "TMD TARGET RMSD AT FIRST STEP COMPUTED FROM INITIAL COORDINATES\n";
04127 }
04128 iout << iINFO << "TMD TARGET RMSD AT FINAL STEP " << TMDFinalRMSD << "\n";
04129 iout << endi;
04130 }
04131
04132 if (symmetryOn) {
04133 if (symmetryLastStep == -1){
04134 iout << iINFO << "SYMMETRY RESTRAINTS ACTIVE BETWEEN STEPS " << symmetryFirstStep << " and " << "INFINITY" << "\n";
04135 }
04136 else{
04137 iout << iINFO << "SYMMETRY RESTRAINTS ACTIVE BETWEEN STEPS " << symmetryFirstStep << " and " << symmetryLastStep << "\n";
04138 }
04139
04140
04141 current = config->find("symmetryFile");
04142 for ( ; current; current = current->next ) {
04143 iout << iINFO << "SYMMETRY FILE " << current->data << "\n";
04144 }
04145
04146 current = config->find("symmetryMatrixFile");
04147 for ( ; current; current = current->next ) {
04148 iout << iINFO << "SYMMETRY MATRIX FILE " << current->data << "\n";
04149 }
04150 iout << iINFO << "SYMMETRY FORCE CONSTANT " << symmetryk << "\n";
04151 if (symmetryScaleForces){
04152 iout << iINFO << "SYMMETRY SCALE FORCES ON\n";
04153 }
04154 iout << iINFO << "SYMMETRY FIRST FULL STEP " << symmetryFirstFullStep << "\n";
04155 if (symmetryLastFullStep == -1){
04156 iout << iINFO << "SYMMETRY LAST FULL STEP " << "INFINITY" << "\n";
04157
04158 }
04159 else {
04160 iout << iINFO << "SYMMETRY LAST FULL STEP " << symmetryLastFullStep << "\n";
04161
04162 }
04163
04164 iout << endi;
04165 }
04166
04167
04168
04169
04170 if (alchFepOn)
04171 {
04172 iout << iINFO << "ALCHEMICAL FEP ON\n";
04173 iout << iINFO << "FEP CURRENT LAMBDA VALUE "
04174 << alchLambda << "\n";
04175 iout << iINFO << "FEP COMPARISON LAMBDA VALUE "
04176 << alchLambda2 << "\n";
04177 iout << iINFO << "FEP VDW SHIFTING COEFFICIENT "
04178 << alchVdwShiftCoeff << "\n";
04179 iout << iINFO << "FEP ELEC. ACTIVE FOR ANNIHILATED "
04180 << "PARTICLES BETWEEN LAMBDA = 0 AND LAMBDA = "
04181 << (1 - alchElecLambdaStart) << "\n";
04182 iout << iINFO << "FEP ELEC. ACTIVE FOR EXNIHILATED "
04183 << "PARTICLES BETWEEN LAMBDA = "
04184 << alchElecLambdaStart << " AND LAMBDA = 1\n";
04185 iout << iINFO << "FEP VDW ACTIVE FOR ANNIHILATED "
04186 << "PARTICLES BETWEEN LAMBDA = "
04187 << (1 - alchVdwLambdaEnd) << " AND LAMBDA = 1\n";
04188 iout << iINFO << "FEP VDW ACTIVE FOR EXNIHILATED "
04189 << "PARTICLES BETWEEN LAMBDA = 0 AND LAMBDA = "
04190 << alchVdwLambdaEnd << "\n";
04191 if (alchFepWCADispOn)
04192 {
04193 iout << iINFO << "FEP WEEKS-CHANDLER-ANDERSEN DECOMPOSITION (DISPERSION) ON\n";
04194 }
04195 if (alchFepWCARepuOn)
04196 {
04197 iout << iINFO << "FEP WEEKS-CHANDLER-ANDERSEN DECOMPOSITION (REPULSION) ON\n";
04198 iout << iINFO << "FEP WEEKS-CHANDLER-ANDERSEN RCUT1 = " << alchFepWCArcut1 << " AND RCUT2 = " << alchFepWCArcut2 << "\n";
04199 }
04200 }
04201
04202
04203 if (alchThermIntOn)
04204 {
04205 iout << iINFO << "THERMODYNAMIC INTEGRATION (TI) ON\n";
04206 iout << iINFO << "TI LAMBDA VALUE "
04207 << alchLambda << "\n";
04208 iout << iINFO << "TI VDW SHIFTING COEFFICIENT "
04209 << alchVdwShiftCoeff << "\n";
04210 iout << iINFO << "TI ELEC. ACTIVE FOR ANNIHILATED "
04211 << "PARTICLES BETWEEN LAMBDA = 0 AND LAMBDA = "
04212 << (1 - alchElecLambdaStart) << "\n";
04213 iout << iINFO << "TI ELEC. ACTIVE FOR EXNIHILATED "
04214 << "PARTICLES BETWEEN LAMBDA = "
04215 << alchElecLambdaStart << " AND LAMBDA = 1\n";
04216 iout << iINFO << "TI VDW ACTIVE FOR ANNIHILATED "
04217 << "PARTICLES BETWEEN LAMBDA = "
04218 << (1 - alchVdwLambdaEnd) << " AND LAMBDA = 1\n";
04219 iout << iINFO << "TI VDW ACTIVE FOR EXNIHILATED "
04220 << "PARTICLES BETWEEN LAMBDA = 0 AND LAMBDA = "
04221 << alchVdwLambdaEnd << "\n";
04222 }
04223
04224
04225 if ( lesOn ) {
04226 iout << iINFO << "LOCALLY ENHANCED SAMPLING ACTIVE\n";
04227 iout << iINFO << "LOCAL ENHANCEMENT FACTOR IS "
04228 << lesFactor << "\n";
04229 if ( lesReduceTemp ) iout << iINFO
04230 << "SCALING ENHANCED ATOM TEMPERATURE BY 1/" << lesFactor << "\n";
04231 if ( lesReduceMass ) iout << iINFO
04232 << "SCALING ENHANCED ATOM MASS BY 1/" << lesFactor << "\n";
04233 }
04234
04235 if ( pairInteractionOn ) {
04236 iout << iINFO << "PAIR INTERACTION CALCULATIONS ACTIVE\n";
04237 iout << iINFO << "USING FLAG " << pairInteractionGroup1
04238 << " FOR GROUP 1\n";
04239 if (pairInteractionSelf) {
04240 iout << iINFO << "COMPUTING ONLY SELF INTERACTIONS FOR GROUP 1 ATOMS\n";
04241 } else {
04242 iout << iINFO << "USING FLAG " << pairInteractionGroup2
04243 << " FOR GROUP 2\n";
04244 }
04245 }
04246
04247 if (consForceOn) {
04248 iout << iINFO << "CONSTANT FORCE ACTIVE\n";
04249 if ( consForceScaling != 1.0 ) {
04250 iout << iINFO << "CONSTANT FORCE SCALING "
04251 << consForceScaling << "\n" << endi;
04252 }
04253 }
04254
04255
04256
04257
04258 if (extForcesOn) {
04259 iout << iINFO << "EXTERNAL COMMAND FORCES ACTIVE\n";
04260 iout << iINFO << "EXT FORCES COMMAND: " << extForcesCommand << "\n";
04261 iout << iINFO << "EXT COORD FILENAME: " << extCoordFilename << "\n";
04262 iout << iINFO << "EXT FORCE FILENAME: " << extForceFilename << "\n";
04263 iout << endi;
04264 }
04265
04266
04267
04268 if (GBISserOn) {
04269 GBISOn = 0;
04270 iout << iINFO<< "GBIS GENERALIZED BORN IMPLICIT SOLVENT ACTIVE (SERIAL)\n";
04271 }
04272 if (GBISOn) {
04273 iout << iINFO << "GBIS GENERALIZED BORN IMPLICIT SOLVENT ACTIVE\n";
04274 }
04275 if (GBISOn || GBISserOn) {
04276 iout << iINFO << "GBIS SOLVENT DIELECTRIC: " << solvent_dielectric<< "\n";
04277 iout << iINFO << "GBIS PROTEIN DIELECTRIC: " << dielectric<< "\n";
04278 iout <<iINFO<<"GBIS COULOMB RADIUS OFFSET: "<< coulomb_radius_offset<<" Ang\n";
04279 iout << iINFO << "GBIS ION CONCENTRATION: " << ion_concentration << " M\n";
04280 iout << iINFO << "GBIS DEBYE SCREENING LENGTH: " << 1.0/kappa << " Ang\n";
04281 iout << iINFO << "GBIS DELTA: " << gbis_delta << "\n";
04282 iout << iINFO << "GBIS BETA: " << gbis_beta << "\n";
04283 iout << iINFO << "GBIS GAMMA: " << gbis_gamma << "\n";
04284 iout << iINFO << "GBIS BORN RADIUS CUTOFF: " << alpha_cutoff << " Ang\n";
04285 iout << iINFO << "GBIS MAX BORN RADIUS: " << alpha_max << " Ang\n";
04286 iout << endi;
04287 }
04288
04289 if (LCPOOn) {
04290 iout << iINFO << "SASA SURFACE TENSION: " << surface_tension<< " kcal/mol/Ang^2\n";
04291 }
04292
04293 tclBCScript = 0;
04294 if (tclBCOn) {
04295 iout << iINFO << "TCL BOUNDARY FORCES ACTIVE\n";
04296 current = config->find("tclBCScript");
04297 if ( current ) {
04298 tclBCScript = current->data;
04299 iout << iINFO << "TCL BOUNDARY FORCES SCRIPT " << current->data << "\n";
04300 }
04301 iout << iINFO << "TCL BOUNDARY FORCES ARGS " << tclBCArgs << "\n";
04302 iout << endi;
04303 }
04304
04305
04306
04307 globalForcesOn = ( tclForcesOn || freeEnergyOn || miscForcesOn ||
04308 (IMDon && ! IMDignore) || SMDOn || TMDOn || colvarsOn || symmetryOn );
04309
04310
04311 if (tclForcesOn)
04312 {
04313 iout << iINFO << "TCL GLOBAL FORCES ACTIVE\n";
04314
04315 current = config->find("tclForcesScript");
04316
04317 for ( ; current; current = current->next ) {
04318
04319 if ( strstr(current->data,"\n") ) {
04320 iout << iINFO << "TCL GLOBAL FORCES SCRIPT INLINED IN CONFIG FILE\n";
04321 continue;
04322 }
04323
04324 iout << iINFO << "TCL GLOBAL FORCES SCRIPT " << current->data << "\n";
04325
04326 }
04327 iout << endi;
04328 }
04329
04330 if (miscForcesOn)
04331 {
04332 iout << iINFO << "MISC FORCES ACTIVE\n";
04333
04334 current = config->find("miscForcesScript");
04335
04336 for ( ; current; current = current->next ) {
04337
04338 if ( strstr(current->data,"\n") ) {
04339 iout << iINFO << "MISC FORCES SCRIPT INLINED IN CONFIG FILE\n";
04340 continue;
04341 }
04342
04343 iout << iINFO << "MISC FORCES SCRIPT " << current->data << "\n";
04344
04345 }
04346 iout << endi;
04347 }
04348
04349 if (freeEnergyOn)
04350 {
04351 iout << iINFO << "FREE ENERGY PERTURBATION ACTIVE\n";
04352
04353 current = config->find("freeEnergyConfig");
04354
04355 for ( ; current; current = current->next ) {
04356
04357 if ( strstr(current->data,"\n") ) {
04358 iout << iINFO << "FREE ENERGY PERTURBATION SCRIPT INLINED IN CONFIG FILE\n";
04359 continue;
04360 }
04361
04362 iout << iINFO << "FREE ENERGY PERTURBATION SCRIPT " << current->data << "\n";
04363
04364 }
04365 iout << endi;
04366 }
04367
04368 if (colvarsOn)
04369 {
04370 iout << iINFO << "COLLECTIVE VARIABLES CALCULATION REQUESTED\n";
04371
04372 current = config->find ("colvarsConfig");
04373 for ( ; current; current = current->next ) {
04374 if ( strstr(current->data,"\n") ) {
04375 iout << iINFO << "COLLECTIVE VARIABLES CONFIGURATION INLINED IN CONFIG FILE\n";
04376 continue;
04377 }
04378 iout << iINFO << "COLLECTIVE VARIABLES CONFIGURATION " << current->data << "\n";
04379 }
04380
04381 current = config->find ("colvarsInput");
04382 for ( ; current; current = current->next ) {
04383 if ( strstr(current->data,"\n") ) {
04384 iout << iINFO << "COLLECTIVE VARIABLES RESTART INFORMATION INLINED IN CONFIG FILE\n";
04385 continue;
04386 }
04387 iout << iINFO << "COLLECTIVE VARIABLES RESTART INFORMATION " << current->data << "\n";
04388 }
04389
04390 iout << endi;
04391 }
04392
04393 if (IMDon)
04394 {
04395 iout << iINFO << "INTERACTIVE MD ACTIVE\n";
04396 iout << iINFO << "INTERACTIVE MD PORT " << IMDport << "\n";
04397 iout << iINFO << "INTERACTIVE MD FREQ " << IMDfreq << "\n";
04398 if (IMDignore) {
04399 iout << iINFO << "INTERACTIVE MD WILL NOT INFLUENCE SIMULATION\n";
04400 } else {
04401 if (IMDwait) iout << iINFO << "WILL AWAIT INTERACTIVE MD CONNECTION\n";
04402 }
04403 iout << endi;
04404 }
04405
04406 if (globalOn && !dihedralOn)
04407 {
04408 iout << iINFO << "GLOBAL INTEGRATION TEST MODE ACTIVE\n";
04409 }
04410
04411
04412 if (dihedralOn)
04413 {
04414 iout << iINFO << "DIHEDRAL ANGLE DYNAMICS ACTIVE\n";
04415 if (!COLDOn)
04416 {
04417 iout << iINFO << "*** DIHEDRAL ANGLE DYNAMICS IS HIGHLY EXPERIMENTAL ***\n";
04418 iout << iINFO << "PLEASE CONSIDER USING THE COLD OPTION AS WELL\n";
04419 }
04420 }
04421
04422
04423
04424
04425
04426
04427
04428
04429 if (COLDOn)
04430 {
04431 iout << iINFO << "COLD (CONSTRAINED OVERDAMPED LANGEVIN DYNAMICS) ACTIVE\n";
04432
04433 iout << iINFO << "COLD TARGET TEMP "
04434 << COLDTemp << "\n";
04435
04436 iout << iINFO << "COLD COLLISION RATE "
04437 << COLDRate << "\n";
04438 }
04439
04440 if (cylindricalBCOn)
04441 {
04442 iout << iINFO << "CYLINDRICAL BOUNDARY CONDITIONS ACTIVE\n";
04443 iout << iINFO << "AXIS " << cylindricalBCAxis << "\n";
04444 iout << iINFO << "RADIUS #1 " << cylindricalBCr1 << "\n";
04445 iout << iINFO << "FORCE CONSTANT #1 " << cylindricalBCk1 << "\n";
04446 iout << iINFO << "EXPONENT #1 " << cylindricalBCexp1 << "\n";
04447 iout << iINFO << "LENGTH #1 " << cylindricalBCl1 << "\n";
04448 if (cylindricalBCr2 > 0.0)
04449 {
04450 iout << iINFO << "RADIUS #2 " << cylindricalBCr2 << "\n";
04451 iout << iINFO << "FORCE CONSTANT #2 " << cylindricalBCk2 << "\n";
04452 iout << iINFO << "EXPONENT #2 " << cylindricalBCexp2 << "\n";
04453 iout << iINFO << "LENGTH #2 " << cylindricalBCl2 << "\n";
04454 }
04455 iout << iINFO << "CYLINDER BOUNDARY CENTER(" << cylindricalCenter.x << ", "
04456 << cylindricalCenter.y << ", " << cylindricalCenter.z << ")\n";
04457 iout << endi;
04458 }
04459
04460 if (sphericalBCOn)
04461 {
04462 iout << iINFO << "SPHERICAL BOUNDARY CONDITIONS ACTIVE\n";
04463
04464 iout << iINFO << "RADIUS #1 "
04465 << sphericalBCr1 << "\n";
04466 iout << iINFO << "FORCE CONSTANT #1 "
04467 << sphericalBCk1 << "\n";
04468 iout << iINFO << "EXPONENT #1 "
04469 << sphericalBCexp1 << "\n";
04470
04471 if (sphericalBCr2 > 0)
04472 {
04473 iout << iINFO << "RADIUS #2 "
04474 << sphericalBCr2 << "\n";
04475 iout << iINFO << "FORCE CONSTANT #2 "
04476 << sphericalBCk2 << "\n";
04477 iout << iINFO << "EXPONENT #2 "
04478 << sphericalBCexp2 << "\n";
04479 }
04480
04481 iout << iINFO << "SPHERE BOUNDARY CENTER(" << sphericalCenter.x << ", "
04482 << sphericalCenter.y << ", " << sphericalCenter.z << ")\n";
04483 iout << endi;
04484 }
04485
04486 if (eFieldOn)
04487 {
04488 iout << iINFO << "ELECTRIC FIELD ACTIVE\n";
04489
04490 iout << iINFO << "E-FIELD VECTOR ("
04491 << eField.x << ", " << eField.y
04492 << ", " << eField.z << ")\n";
04493 iout << iINFO << "E-FIELD FREQUENCY IS (1/ps) " << eFieldFreq << "\n";
04494 iout << iINFO << "E-FIELD PHASE IS (deg) " << eFieldPhase << "\n";
04495
04496 iout << endi;
04497 }
04498
04499 if (stirOn)
04500 {
04501 iout << iINFO << "STIRRING TORQUES ACTIVE\n";
04502
04503 iout << iINFO << "STIR STARTING THETA (deg) "<< stirStartingTheta << "\n";
04504 iout << iINFO << "STIR ANGULAR VELOCITY (deg/ts) " << stirVel <<"\n";
04505 iout << iINFO << "STIR FORCE HARMONIC SPRING CONSTANT "<< stirK << "\n";
04506 iout << iINFO << "STIR AXIS OF ROTATION (DIRECTION) ("
04507 << stirAxis.x << ", " << stirAxis.y
04508 << ", " << stirAxis.z << ")\n";
04509 iout << iINFO << "STIR PIVOT POINT (COORDINATE) ("
04510 << stirPivot.x << ", " << stirPivot.y
04511 << ", " << stirPivot.z << ")\n";
04512 current = config->find("stirFilename");
04513
04514 iout << iINFO << "STIR ATOMS AND ORIGINAL POSITIONS FROM FILE " <<current ->data << '\n';
04515 current = config->find("stirredAtomsCol");
04516 iout << iINFO <<"STIR FILE COLUMN " << current ->data << '\n';
04517 iout << endi;
04518 }
04519
04520 if (drudeOn)
04521 {
04522 iout << iINFO << "DRUDE MODEL DUAL THERMOSTAT IS ACTIVE\n";
04523 iout << iINFO << "DRUDE BOND TEMPERATURE " << drudeTemp << "\n";
04524 if (drudeDamping > 0.0) {
04525 iout << iINFO << "DRUDE DAMPING COEFFICIENT IS "
04526 << drudeDamping << " INVERSE PS\n";
04527 }
04528 if (drudeBondConst > 0.0) {
04529 iout << iINFO << "DRUDE QUARTIC RESTRAINT IS ACTIVE FOR DRUDE BONDS\n";
04530 iout << iINFO << "DRUDE MAXIMUM BOND LENGTH BEFORE RESTRAINT IS "
04531 << drudeBondLen << "\n";
04532 iout << iINFO << "DRUDE BOND RESTRAINT CONSTANT IS "
04533 << drudeBondConst << "\n";
04534 }
04535 if (drudeNbtholeCut > 0.0) {
04536 iout << iINFO << "DRUDE NBTHOLE IS ACTIVE\n";
04537 iout << iINFO << "DRUDE NBTHOLE RADIUS IS "
04538 << drudeNbtholeCut << "\n";
04539 }
04540 }
04541
04542 if (langevinOn)
04543 {
04544 iout << iINFO << "LANGEVIN DYNAMICS ACTIVE\n";
04545 iout << iINFO << "LANGEVIN TEMPERATURE "
04546 << langevinTemp << "\n";
04547 if (langevinDamping > 0.0) {
04548 iout << iINFO << "LANGEVIN DAMPING COEFFICIENT IS "
04549 << langevinDamping << " INVERSE PS\n";
04550 if (langevinHydrogen)
04551 iout << iINFO << "LANGEVIN DYNAMICS APPLIED TO HYDROGENS\n";
04552 else
04553 iout << iINFO << "LANGEVIN DYNAMICS NOT APPLIED TO HYDROGENS\n";
04554 } else {
04555 iout << iINFO << "LANGEVIN DAMPING COEFFICIENTS DETERMINED FROM FILES\n";
04556 current = config->find("langevinFile");
04557 if ( current ) iout << iINFO << "LANGEVIN DAMPING FILE: " <<
04558 current->data << "\n";
04559 else iout << iINFO << "LANGEVIN DAMPING FILE IS COORDINATE PDB\n";
04560 current = config->find("langevinCol");
04561 if ( current ) iout << iINFO << "LANGEVIN DAMPING COLUMN: " <<
04562 current->data << "\n";
04563 else iout << iINFO << "LANGEVIN DAMPING COLUMN: DEFAULT (4TH, O)\n";
04564 }
04565 iout << endi;
04566 }
04567
04568
04569 if (loweAndersenOn)
04570 {
04571 iout << iINFO << "LOWE-ANDERSEN DYNAMICS ACTIVE\n";
04572 iout << iINFO << "LOWE-ANDERSEN TEMPERATURE "
04573 << loweAndersenTemp << " K\n";
04574 iout << iINFO << "LOWE-ANDERSEN RATE "
04575 << loweAndersenRate << " INVERSE PS\n";
04576 iout << iINFO << "LOWE-ANDERSEN CUTOFF "
04577 << loweAndersenCutoff << " ANGSTROMS\n";
04578 iout << endi;
04579 }
04580
04581
04582 if (tCoupleOn)
04583 {
04584 iout << iINFO << "TEMPERATURE COUPLING ACTIVE\n";
04585 iout << iINFO << "COUPLING TEMPERATURE "
04586 << tCoupleTemp << "\n";
04587 iout << endi;
04588 }
04589
04590 if (minimizeOn)
04591 {
04592 iout << iINFO << "OLD STYLE MINIMIZATION ACTIVE\n";
04593 iout << endi;
04594 }
04595
04596 if (minimizeCGOn)
04597 {
04598 iout << iINFO << "CONJUGATE GRADIENT MINIMIZATION ACTIVE\n";
04599 iout << iINFO << "LINE MINIMIZATION GOAL = " << minLineGoal << "\n";
04600 iout << iINFO << "BABY STEP SIZE = " << minBabyStep << "\n";
04601 iout << iINFO << "TINY STEP SIZE = " << minTinyStep << "\n";
04602 iout << endi;
04603 }
04604
04605 if (maximumMove)
04606 {
04607 iout << iINFO << "MAXIMUM MOVEMENT "
04608 << maximumMove << "\n";
04609 iout << endi;
04610 }
04611
04612 if (rescaleFreq > 0)
04613 {
04614 iout << iINFO << "VELOCITY RESCALE FREQ "
04615 << rescaleFreq << "\n";
04616 iout << iINFO << "VELOCITY RESCALE TEMP "
04617 << rescaleTemp << "\n";
04618 iout << endi;
04619 }
04620
04621 if (reassignFreq > 0)
04622 {
04623 iout << iINFO << "VELOCITY REASSIGNMENT FREQ "
04624 << reassignFreq << "\n";
04625 iout << iINFO << "VELOCITY REASSIGNMENT TEMP "
04626 << reassignTemp << "\n";
04627 if ( reassignIncr != 0. )
04628 iout << iINFO << "VELOCITY REASSIGNMENT INCR "
04629 << reassignIncr << "\n";
04630 if ( reassignHold != 0. )
04631 iout << iINFO << "VELOCITY REASSIGNMENT HOLD "
04632 << reassignHold << "\n";
04633 iout << endi;
04634 }
04635
04636 if (berendsenPressureOn && langevinPistonOn)
04637 {
04638 NAMD_die("Multiple pressure control algorithms selected!\n");
04639 }
04640
04641 if (excludeFromPressure) {
04642 iout << iINFO << "EXCLUDE FROM PRESSURE ACTIVE\n";
04643 }
04644 if (useConstantArea && useConstantRatio) {
04645 NAMD_die("useConstantArea and useConstantRatio are mutually exclusive.\n");
04646 }
04647 if (useConstantRatio && !useFlexibleCell) {
04648 NAMD_die("useConstantRatio requires useFlexibleCell.\n");
04649 }
04650 if (useConstantArea && surfaceTensionTarget) {
04651 NAMD_die("surfaceTensionTarget and useConstantArea are mutually exclusive.\n");
04652 }
04653 if (useConstantArea && !useFlexibleCell) {
04654 NAMD_die("useConstantArea requires useFlexibleCell.\n");
04655 }
04656
04657 if (berendsenPressureOn || langevinPistonOn) {
04658 if (rigidBonds != RIGID_NONE && useGroupPressure == FALSE) {
04659 useGroupPressure = TRUE;
04660 iout << iWARN << "Option useGroupPressure is being enabled "
04661 << "due to pressure control with rigidBonds.\n" << endi;
04662 }
04663 }
04664
04665 if (berendsenPressureOn)
04666 {
04667 if ( ! opts.defined("BerendsenPressureFreq") ) {
04668 berendsenPressureFreq = nonbondedFrequency;
04669 if ( fullElectFrequency )
04670 berendsenPressureFreq = fullElectFrequency;
04671 }
04672 if ( (berendsenPressureFreq % nonbondedFrequency) || ( fullElectFrequency
04673 && (berendsenPressureFreq % fullElectFrequency) ) )
04674 NAMD_die("berendsenPressureFreq must be a multiple of both fullElectFrequency and nonbondedFrequency\n");
04675 iout << iINFO << "BERENDSEN PRESSURE COUPLING ACTIVE\n";
04676 iout << iINFO << " TARGET PRESSURE IS "
04677 << berendsenPressureTarget << " BAR\n";
04678 iout << iINFO << " COMPRESSIBILITY ESTIMATE IS "
04679 << berendsenPressureCompressibility << " BAR^(-1)\n";
04680 iout << iINFO << " RELAXATION TIME IS "
04681 << berendsenPressureRelaxationTime << " FS\n";
04682 iout << iINFO << " APPLIED EVERY "
04683 << berendsenPressureFreq << " STEPS\n";
04684 iout << iINFO << " PRESSURE CONTROL IS "
04685 << (useGroupPressure?"GROUP":"ATOM") << "-BASED\n";
04686 iout << endi;
04687 berendsenPressureTarget /= PRESSUREFACTOR;
04688 berendsenPressureCompressibility *= PRESSUREFACTOR;
04689 }
04690
04691 if (langevinPistonOn)
04692 {
04693 iout << iINFO << "LANGEVIN PISTON PRESSURE CONTROL ACTIVE\n";
04694 iout << iINFO << " TARGET PRESSURE IS "
04695 << langevinPistonTarget << " BAR\n";
04696 iout << iINFO << " OSCILLATION PERIOD IS "
04697 << langevinPistonPeriod << " FS\n";
04698 iout << iINFO << " DECAY TIME IS "
04699 << langevinPistonDecay << " FS\n";
04700 iout << iINFO << " PISTON TEMPERATURE IS "
04701 << langevinPistonTemp << " K\n";
04702 iout << iINFO << " PRESSURE CONTROL IS "
04703 << (useGroupPressure?"GROUP":"ATOM") << "-BASED\n";
04704 iout << iINFO << " INITIAL STRAIN RATE IS "
04705 << strainRate << "\n";
04706 iout << endi;
04707 langevinPistonTarget /= PRESSUREFACTOR;
04708 }
04709
04710 if (berendsenPressureOn || langevinPistonOn) {
04711 iout << iINFO << " CELL FLUCTUATION IS "
04712 << (useFlexibleCell?"AN":"") << "ISOTROPIC\n";
04713 if (useConstantRatio)
04714 iout << iINFO << " SHAPE OF CELL IS CONSTRAINED IN X-Y PLANE\n";
04715 if (useConstantArea)
04716 iout << iINFO << " CONSTANT AREA PRESSURE CONTROL ACTIVE\n";
04717 }
04718
04719 if (surfaceTensionTarget != 0)
04720 {
04721 iout << iINFO << "SURFACE TENSION CONTROL ACTIVE\n";
04722 iout << iINFO << " TARGET SURFACE TENSION IS "
04723 << surfaceTensionTarget << " DYN/CM\n";
04724 iout << endi;
04725
04726
04727 surfaceTensionTarget *= 100.0 / PRESSUREFACTOR;
04728 }
04729
04730 if (pressureProfileOn) {
04731 if ((berendsenPressureOn || langevinPistonOn) && !dcdUnitCell) {
04732 #if 1
04733 iout << iWARN << "Turning on dcdUnitCell so that trajectory files contain unit cell data.\n" << endi;
04734 dcdUnitCell = 1;
04735 #else
04736 NAMD_die("Sorry, pressure profile not implemented for constant pressure.");
04737 #endif
04738 }
04739
04740 if (pressureProfileEwaldOn)
04741 pressureProfileOn = 0;
04742
04743 if (pressureProfileSlabs < 1)
04744 NAMD_die("pressureProfileSlabs must be positive.");
04745 iout << iINFO << "PRESSURE PROFILE CALCULATIONS ACTIVE\n";
04746 iout << iINFO << " NUMBER OF SLABS: " << pressureProfileSlabs << "\n";
04747 iout << iINFO << " SLAB THICKNESS: " << cellBasisVector3.z / pressureProfileSlabs
04748 << "\n";
04749 iout << iINFO << " TIMESTEPS BETWEEN DATA OUTPUT: "
04750 << pressureProfileFreq << "\n";
04751 iout << iINFO << " NUMBER OF ATOM TYPES: " << pressureProfileAtomTypes << "\n";
04752 iout << endi;
04753 } else {
04754 pressureProfileEwaldOn = 0;
04755 pressureProfileAtomTypes = 1;
04756 }
04757
04758 if (accelMDOn) {
04759 iout << iINFO << "ACCELERATED MD ACTIVE\n";
04760
04761 if ( accelMDdual) {
04762 accelMDdihe = FALSE;
04763 iout << iINFO << "APPLYING DUAL BOOST\n";
04764 }
04765 else if ( accelMDdihe ) {
04766 iout << iINFO << "BOOSTING DIHEDRAL POTENTIAL\n";
04767 } else {
04768 iout << iINFO << "BOOSTING TOTAL POTENTIAL\n";
04769 }
04770
04771 iout << iINFO << "accelMDE: " << accelMDE << " KCAL/MOL, accelMDalpha: " << accelMDalpha << " KCAL/MOL\n";
04772 if (accelMDdual) {
04773 iout << iINFO << "accelMDTE: " << accelMDTE << " KCAL/MOL, accelMDTalpha: " << accelMDTalpha << " KCAL/MOL\n";
04774 }
04775 if ( accelMDLastStep > 0) {
04776 iout << iINFO << "accelMD WILL BE DONE FROM STEP " << accelMDFirstStep << " TO STEP " << accelMDLastStep << "\n";
04777 } else {
04778 iout << iINFO << "accelMD WILL BE DONE FROM STEP " << accelMDFirstStep << " TO THE END OF THE SIMULATION \n";
04779 }
04780 iout << iINFO << "accelMD OUTPUT FREQUENCY " << accelMDOutFreq << "\n";
04781 iout << endi;
04782 }
04783
04784 if (adaptTempOn) {
04785 iout << iINFO << "ADAPTIVE TEMPERING ACTIVE:\n";
04786 iout << iINFO << " OUTPUT FREQUENCY: " << adaptTempOutFreq << "\n";
04787 iout << iINFO << " TEMPERATURE UPDATE FREQUENCY: " << adaptTempFreq << "\n";
04788 if ( adaptTempLastStep > 0 )
04789 iout << iINFO << " ADAPTIVE TEMPERING WILL BE DONE FROM STEP " << adaptTempFirstStep << " TO " << adaptTempLastStep << "\n";
04790 else
04791 iout << iINFO << " ADAPTIVE TEMPERING WILL BE DONE FROM STEP " << adaptTempFirstStep << "\n";
04792 if ( adaptTempLangevin )
04793 iout << iINFO << " ADAPTIVE TEMPERING COUPLED TO LANGEVIN THERMOSTAT\n";
04794 if ( adaptTempRescale )
04795 iout << iINFO << " ADAPTIVE TEMPERING COUPLED TO VELOCITY RESCALING\n";
04796 if (adaptTempRestartFreq > 0) {
04797 iout << iINFO << " WRITING RESTART INFORMATION TO " << adaptTempRestartFile << " EVERY " << adaptTempRestartFreq << " STEPS\n";
04798 }
04799
04800 }
04801
04802 if (FMAOn)
04803 {
04804 iout << iINFO << "FMA ACTIVE\n";
04805 iout << iINFO << "FMA THETA "
04806 << fmaTheta << "\n";
04807 iout << endi;
04808 }
04809
04810 FFTWWisdomString = 0;
04811 if (PMEOn)
04812 {
04813 iout << iINFO << "PARTICLE MESH EWALD (PME) ACTIVE\n";
04814 iout << iINFO << "PME TOLERANCE "
04815 << PMETolerance << "\n";
04816 iout << iINFO << "PME EWALD COEFFICIENT "
04817 << PMEEwaldCoefficient << "\n";
04818 iout << iINFO << "PME INTERPOLATION ORDER "
04819 << PMEInterpOrder << "\n";
04820 iout << iINFO << "PME GRID DIMENSIONS "
04821 << PMEGridSizeX << " "
04822 << PMEGridSizeY << " "
04823 << PMEGridSizeZ << "\n";
04824 iout << iINFO << "PME MAXIMUM GRID SPACING "
04825 << PMEGridSpacing << "\n";
04826 if ( PMEBarrier ) {
04827 iout << iINFO << "PME BARRIER ENABLED\n";
04828 }
04829 iout << endi;
04830 if ( useDPME ) iout << iINFO << "USING OLD DPME CODE\n";
04831 #ifdef NAMD_FFTW
04832 else if ( FFTWUseWisdom ) {
04833 #ifdef NAMD_FFTW_3
04834 iout << iINFO << "Attempting to read FFTW data from system" <<"\n" <<endi;
04835 fftwf_import_system_wisdom();
04836 #endif
04837 if (! opts.defined("FFTWWisdomFile")) {
04838 strcpy(FFTWWisdomFile,"FFTW_NAMD_");
04839 strcat(FFTWWisdomFile,NAMD_VERSION);
04840 strcat(FFTWWisdomFile,"_");
04841 strcat(FFTWWisdomFile,NAMD_PLATFORM);
04842 #ifdef NAMD_FFTW_3
04843 strcat(FFTWWisdomFile,"_FFTW3");
04844 #endif
04845 strcat(FFTWWisdomFile,".txt");
04846 }
04847
04848 iout << iINFO << "Attempting to read FFTW data from "
04849 << FFTWWisdomFile << "\n" << endi;
04850 FILE *wisdom_file = fopen(FFTWWisdomFile,"r");
04851 if ( wisdom_file ) {
04852 #ifdef NAMD_FFTW_3
04853 fftwf_import_wisdom_from_file(wisdom_file);
04854 #else
04855 fftw_import_wisdom_from_file(wisdom_file);
04856 #endif
04857 fclose(wisdom_file);
04858 }
04859 int nrp = 1;
04860
04861
04862 int minslices = PMEMinSlices;
04863 int dimx = PMEGridSizeX;
04864 int nrpx = ( dimx + minslices - 1 ) / minslices;
04865 if ( nrpx > nrp ) nrp = nrpx;
04866 int dimy = PMEGridSizeY;
04867 int nrpy = ( dimy + minslices - 1 ) / minslices;
04868 if ( nrpy > nrp ) nrp = nrpy;
04869
04870
04871 int nrpp = CkNumPes();
04872
04873 if ( nrpp < nrp ) nrp = nrpp;
04874
04875
04876 int nrps = PMEProcessors;
04877 if ( nrps > CkNumPes() ) nrps = CkNumPes();
04878 if ( nrps > 0 ) nrp = nrps;
04879
04880
04881 int bx = ( dimx + nrp - 1 ) / nrp;
04882 int nrpbx = ( dimx + bx - 1 ) / bx;
04883 int by = ( dimy + nrp - 1 ) / nrp;
04884 int nrpby = ( dimy + by - 1 ) / by;
04885 nrp = ( nrpby > nrpbx ? nrpby : nrpbx );
04886 if ( bx != ( dimx + nrp - 1 ) / nrp )
04887 NAMD_bug("Error in selecting number of PME processors.");
04888 if ( by != ( dimy + nrp - 1 ) / nrp )
04889 NAMD_bug("Error in selecting number of PME processors.");
04890
04891
04892
04893
04894 int block2 = (PMEGridSizeY + nrp - 1) / nrp;
04895 int block2_min = PMEGridSizeY % block2;
04896 if ( ! block2_min ) block2_min = block2;
04897 int dim3 = 2 * (PMEGridSizeZ/2 + 1);
04898
04899 int n[3]; n[0] = PMEGridSizeX; n[1] = PMEGridSizeY; n[2] = PMEGridSizeZ;
04900 fftw_complex *work = new fftw_complex[n[0]];
04901 float *grid1 = (float *) fftwf_malloc(sizeof(float) *n[1]*dim3);
04902 float *grid2 = (float *) fftwf_malloc(sizeof(float) *n[0]*block2*dim3*2);
04903 iout << iINFO << "Optimizing 6 FFT steps. 1..." << endi;
04904 #ifdef NAMD_FFTW_3
04905 int fftwFlags = FFTWPatient ? FFTW_PATIENT : FFTWEstimate ? FFTW_ESTIMATE : FFTW_MEASURE ;
04906 int planLineSizes[1];
04907 planLineSizes[0]=n[2];
04908 int nx= n[0];
04909 int ny=block2;
04910 int sizeLines=nx*ny;
04911 int zdim = dim3;
04912 int nz=zdim;
04913 int xStride=block2 * dim3 / 2;
04914 fftwf_destroy_plan(
04915 fftwf_plan_many_dft_r2c(1, planLineSizes, sizeLines,
04916 (float *) grid2, NULL, 1,
04917 dim3,
04918 (fftwf_complex *) grid2,
04919 NULL, 1,
04920 dim3/2,
04921 fftwFlags));
04922
04923 iout << " 2..." << endi;
04924 fftwf_destroy_plan(
04925 fftwf_plan_many_dft_c2r(1, planLineSizes, sizeLines,
04926 (fftwf_complex *) grid2,
04927 NULL, 1,
04928 dim3/2,
04929 (float *) grid2, NULL, 1,
04930 dim3,
04931 fftwFlags));
04932 iout << " 3..." << endi;
04933 sizeLines=nz;
04934 planLineSizes[0]=block2;
04935 fftwf_destroy_plan(fftwf_plan_many_dft(1, planLineSizes, sizeLines,
04936 (fftwf_complex *) grid2, NULL,
04937 sizeLines, 1,
04938 (fftwf_complex *) grid2, NULL,
04939 sizeLines, 1,
04940 FFTW_FORWARD,
04941 fftwFlags));
04942 iout << " 4..." << endi;
04943 fftwf_destroy_plan(fftwf_plan_many_dft(1, planLineSizes, sizeLines,
04944 (fftwf_complex *) grid2, NULL,
04945 sizeLines, 1,
04946 (fftwf_complex *) grid2, NULL,
04947 sizeLines, 1,
04948 FFTW_FORWARD,
04949 fftwFlags));
04950 iout << " 5..." << endi;
04951 sizeLines=ny*nz;
04952 planLineSizes[0]=n[0];
04953 fftwf_destroy_plan(fftwf_plan_many_dft(1, planLineSizes, sizeLines,
04954 (fftwf_complex *) grid2, NULL,
04955 sizeLines, 1,
04956 (fftwf_complex *) grid2, NULL,
04957 sizeLines, 1,
04958 FFTW_FORWARD,
04959 fftwFlags));
04960 iout << " 6..." << endi;
04961 fftwf_destroy_plan(fftwf_plan_many_dft(1, planLineSizes, sizeLines,
04962 (fftwf_complex *) grid2, NULL,
04963 sizeLines, 1,
04964 (fftwf_complex *) grid2, NULL,
04965 sizeLines, 1,
04966 FFTW_BACKWARD,
04967 fftwFlags));
04968
04969 #else
04970 rfftwnd_destroy_plan( rfftwnd_create_plan_specific(
04971 2, n+1, FFTW_REAL_TO_COMPLEX,
04972 ( FFTWEstimate ? FFTW_ESTIMATE : FFTW_MEASURE )
04973 | FFTW_IN_PLACE | FFTW_USE_WISDOM, grid1, 1, 0, 0) );
04974 iout << " 2..." << endi;
04975 fftw_destroy_plan( fftw_create_plan_specific(n[0], FFTW_REAL_TO_COMPLEX,
04976 ( FFTWEstimate ? FFTW_ESTIMATE : FFTW_MEASURE )
04977 | FFTW_IN_PLACE | FFTW_USE_WISDOM, (fftw_complex *) grid2,
04978 block2*dim3/2, work, 1) );
04979 iout << " 3..." << endi;
04980 fftw_destroy_plan( fftw_create_plan_specific(n[0], FFTW_REAL_TO_COMPLEX,
04981 ( FFTWEstimate ? FFTW_ESTIMATE : FFTW_MEASURE )
04982 | FFTW_IN_PLACE | FFTW_USE_WISDOM, (fftw_complex *) grid2,
04983 block2_min*dim3/2, work, 1) );
04984 iout << " 4..." << endi;
04985 fftw_destroy_plan( fftw_create_plan_specific(n[0], FFTW_COMPLEX_TO_REAL,
04986 ( FFTWEstimate ? FFTW_ESTIMATE : FFTW_MEASURE )
04987 | FFTW_IN_PLACE | FFTW_USE_WISDOM, (fftw_complex *) grid2,
04988 block2*dim3/2, work, 1) );
04989 iout << " 5..." << endi;
04990 fftw_destroy_plan( fftw_create_plan_specific(n[0], FFTW_COMPLEX_TO_REAL,
04991 ( FFTWEstimate ? FFTW_ESTIMATE : FFTW_MEASURE )
04992 | FFTW_IN_PLACE | FFTW_USE_WISDOM, (fftw_complex *) grid2,
04993 block2_min*dim3/2, work, 1) );
04994 iout << " 6..." << endi;
04995 rfftwnd_destroy_plan( rfftwnd_create_plan_specific(
04996 2, n+1, FFTW_COMPLEX_TO_REAL,
04997 ( FFTWEstimate ? FFTW_ESTIMATE : FFTW_MEASURE )
04998 | FFTW_IN_PLACE | FFTW_USE_WISDOM, grid1, 1, 0, 0) );
04999 #endif
05000 iout << " Done.\n" << endi;
05001 delete [] work;
05002 fftwf_free(grid1);
05003 fftwf_free(grid2);
05004
05005 iout << iINFO << "Writing FFTW data to "
05006 << FFTWWisdomFile << "\n" << endi;
05007 wisdom_file = fopen(FFTWWisdomFile,"w");
05008 if ( wisdom_file ) {
05009 #ifdef NAMD_FFTW_3
05010 fftwf_export_wisdom_to_file(wisdom_file);
05011 #else
05012 fftw_export_wisdom_to_file(wisdom_file);
05013 #endif
05014 fclose(wisdom_file);
05015 }
05016
05017 #ifdef NAMD_FFTW_3
05018 FFTWWisdomString = fftwf_export_wisdom_to_string();
05019 #else
05020 FFTWWisdomString = fftw_export_wisdom_to_string();
05021 #endif
05022 }
05023 #endif
05024 iout << endi;
05025 }
05026 if (fullDirectOn)
05027 {
05028 iout << iINFO << "DIRECT FULL ELECTROSTATIC CALCULATIONS ACTIVE\n";
05029 iout << endi;
05030 }
05031 if (MSMOn)
05032 {
05033 iout << iINFO
05034 << "MULTILEVEL SUMMATION METHOD (MSM) FOR ELECTROSTATICS ACTIVE\n";
05035 #if !defined(CHARM_HAS_MSA)
05036
05037 #endif
05038 if (MsmSerialOn) {
05039 iout << iINFO
05040 << "PERFORMING SERIAL MSM CALCULATION FOR LONG-RANGE PART\n";
05041 iout << endi;
05042 }
05043 }
05044
05045 if ( FMAOn || PMEOn || MSMOn || fullDirectOn || GBISOn)
05046 {
05047 iout << iINFO << "FULL ELECTROSTATIC EVALUATION FREQUENCY "
05048 << fullElectFrequency << "\n";
05049 iout << endi;
05050
05051 if ( ( outputEnergies % fullElectFrequency ) &&
05052 ( fullElectFrequency % outputEnergies ) )
05053 NAMD_die("Either outputEnergies must be a multiple of fullElectFrequency or vice versa.\n");
05054 }
05055
05056 if (MTSAlgorithm == NAIVE)
05057 {
05058 iout << iINFO << "USING NAIVE (CONSTANT FORCE) MTS SCHEME.\n" << endi;
05059 }
05060 if (MTSAlgorithm == VERLETI )
05061 {
05062 iout << iINFO << "USING VERLET I (r-RESPA) MTS SCHEME.\n" << endi;
05063 }
05064
05065 if (longSplitting == SHARP)
05066 iout << iINFO << "SHARP SPLITTING OF LONG RANGE ELECTROSTATICS\n";
05067 else if (longSplitting == XPLOR)
05068 iout << iINFO << "XPLOR SPLITTING OF LONG RANGE ELECTROSTATICS\n";
05069 else if (longSplitting == C1)
05070 iout << iINFO << "C1 SPLITTING OF LONG RANGE ELECTROSTATICS\n";
05071 else if (longSplitting == C2)
05072 iout << iINFO << "C2 SPLITTING OF LONG RANGE ELECTROSTATICS\n";
05073
05074 if (splitPatch == SPLIT_PATCH_POSITION)
05075 iout << iINFO << "PLACING ATOMS IN PATCHES BY POSITION\n";
05076 else if (splitPatch == SPLIT_PATCH_HYDROGEN)
05077 iout << iINFO << "PLACING ATOMS IN PATCHES BY HYDROGEN GROUPS\n";
05078
05079 iout << endi;
05080
05081 if (mollyOn)
05082 {
05083 iout << iINFO << "SLOW FORCE MOLLIFICATION : \n";
05084 iout << iINFO << " ERROR TOLERANCE : " << mollyTol << "\n";
05085 iout << iINFO << " MAX ITERATIONS : " << mollyIter << "\n";
05086 iout << endi;
05087 }
05088
05089 if (rigidBonds != RIGID_NONE)
05090 {
05091 iout << iINFO << "RIGID BONDS TO HYDROGEN : ";
05092 if (rigidBonds == RIGID_ALL) iout << "ALL\n";
05093 if (rigidBonds == RIGID_WATER) iout << "WATER\n";
05094 iout << iINFO << " ERROR TOLERANCE : " << rigidTol << "\n";
05095 iout << iINFO << " MAX ITERATIONS : " << rigidIter << "\n";
05096 if (useSettle) iout << iINFO << "RIGID WATER USING SETTLE ALGORITHM\n";
05097 iout << endi;
05098 }
05099
05100
05101 if (nonbondedFrequency != 1)
05102 {
05103 iout << iINFO << "NONBONDED FORCES EVALUATED EVERY " << nonbondedFrequency << " STEPS\n";
05104 }
05105
05106 iout << iINFO << "RANDOM NUMBER SEED "
05107 << randomSeed << "\n";
05108
05109 iout << endi;
05110
05111 iout << iINFO << "USE HYDROGEN BONDS? ";
05112 if (HydrogenBonds)
05113 {
05114 iout << "YES\n" << endi;
05115 iout << iINFO << "USE ANTECEDENT ATOMS? ";
05116 iout << (useAntecedent ? "YES" : "NO");
05117 iout << "\nHB DIST CUT, ON, OFF ";
05118 iout << daCutoffDist << " , " << daOnDist << " , " << daOffDist;
05119 iout << "\nHB ANGLE CUT, ON, OFF ";
05120 iout << dhaCutoffAngle << " , " << dhaOnAngle << " , ";
05121 iout << dhaOffAngle;
05122 iout << "\nHB ATT, REP exponents ";
05123 iout << distAttExp << " , " << distRepExp;
05124 iout << "\nHB AA, HA exponents ";
05125 iout << aaAngleExp << " , " << haAngleExp;
05126 iout << "\n" << endi;
05127 }
05128 else
05129 {
05130 iout << "NO\n" << endi;
05131 }
05132
05133
05134
05135 if (amberOn)
05136 { iout << iINFO << "Using AMBER format force field!\n";
05137 current = config->find("parmfile");
05138 iout << iINFO << "AMBER PARM FILE " << current->data << '\n';
05139 if (opts.defined("coordinates"))
05140 { current = config->find("coordinates");
05141 iout << iINFO << "COORDINATE PDB " << current->data << '\n';
05142 }
05143 else
05144 { current = config->find("ambercoor");
05145 iout << iINFO << "AMBER COORDINATE FILE " << current->data << '\n';
05146 }
05147 if (readExclusions)
05148 iout << iINFO << "Exclusions will be read from PARM file!\n";
05149 else
05150 iout << iINFO << "Exclusions in PARM file will be ignored!\n";
05151 iout << iINFO << "SCNB (VDW SCALING) " << vdwscale14 << "\n" << endi;
05152 }
05153 else if(gromacsOn)
05154 {
05155 iout << iINFO << "Using GROMACS format force field!\n";
05156
05157 current = config->find("grotopfile");
05158
05159 if (current == NULL)
05160 NAMD_die("no GROMACS topology file defined!?");
05161 iout << iINFO << "GROMACS TOPO FILE " << current->data << '\n';
05162
05163
05164 current = config->find("grocoorfile");
05165 if (current == NULL) {
05166 current = config->find("coordinates");
05167 if (current == NULL) {
05168 NAMD_die("no coordinate file defined!?");
05169 }
05170 }
05171 iout << iINFO << "GROMACS COOR FILE " << current->data << '\n'
05172 << endi;
05173
05174 }
05175 else {
05176 if ( !usePluginIO ) {
05177 if ( current = config->find("coordinates") )
05178 iout << iINFO << "COORDINATE PDB " << current->data << '\n' << endi;
05179 }
05180
05181 current = config->find("structure");
05182
05183 iout << iINFO << "STRUCTURE FILE "
05184 << current->data << "\n" << endi;
05185
05186 if (cosAngles)
05187 {
05188 iout << iINFO << "COSANGLES ON. SOME ANGLES WILL BE COSINE-BASED\n" << endi;
05189 }
05190
05191
05192 if (paraTypeXplorOn)
05193 {
05194 iout << iINFO << "PARAMETER file: XPLOR format! (default) \n" << endi;
05195 }
05196 else if (paraTypeCharmmOn)
05197 {
05198 iout << iINFO << "PARAMETER file: CHARMM format! \n" << endi;
05199 }
05200
05201
05202 current = config->find("parameters");
05203
05204 while (current != NULL)
05205 {
05206 iout << iINFO << "PARAMETERS "
05207 << current->data << "\n" << endi;
05208 current = current->next;
05209 }
05210 }
05211
05212 iout << iINFO << "USING " <<
05213 ( vdwGeometricSigma ? "GEOMETRIC" : "ARITHMETIC" ) <<
05214 " MEAN TO COMBINE L-J SIGMA PARAMETERS\n" << endi;
05215
05216 if (opts.defined("bincoordinates"))
05217 {
05218 current = config->find("bincoordinates");
05219
05220 iout << iINFO << "BINARY COORDINATES "
05221 << current->data << "\n";
05222 }
05223
05224
05225 if (firstTimestep)
05226 {
05227 iout << iINFO << "FIRST TIMESTEP "
05228 << firstTimestep << "\n" << endi;
05229 }
05230 }
05231
05232
05233
05234
05235
05236
05237
05238
05239
05240 void SimParameters::parse_mgrid_params(ConfigList *config)
05241 {
05242 StringList *current;
05243
05244 mgridforcelist.clear();
05245 char *key = new char[81];
05246 char *valstr = new char[256];
05247
05248 if (gridforceOn) {
05249 mgridforceOn = TRUE;
05250 char *default_key = MGRIDFORCEPARAMS_DEFAULTKEY;
05251 MGridforceParams* mgfp = NULL;
05252 mgfp = mgridforcelist.find_key(default_key);
05253 if (mgfp != NULL) {
05254 iout << iINFO << "MGRIDFORCEPOTFILE key "
05255 << key << " redefined for file " << valstr << "\n" << endi;
05256 } else {
05257 mgfp = mgridforcelist.add(default_key);
05258 }
05259 mgfp->gridforceVolts = gridforceVolts;
05260 mgfp->gridforceScale = gridforceScale;
05261
05262 parse_mgrid_string_param(config,"gridforcefile",&(mgfp->gridforceFile));
05263 parse_mgrid_string_param(config,"gridforcecol",&(mgfp->gridforceCol));
05264 parse_mgrid_string_param(config,"gridforcechargecol",&(mgfp->gridforceQcol));
05265 parse_mgrid_string_param(config,"gridforcepotfile",&(mgfp->gridforceVfile));
05266
05267 mgfp->gridforceCont[0] = gridforceContA1;
05268 mgfp->gridforceCont[1] = gridforceContA2;
05269 mgfp->gridforceCont[2] = gridforceContA3;
05270 mgfp->gridforceVOffset = gridforceVOffset;
05271
05272 mgfp->gridforceLite = gridforceLite;
05273 }
05274
05275
05276 current = config->find("mgridforcepotfile");
05277 while (current != NULL) {
05278 int curlen = strlen(current->data);
05279
05280
05281 sscanf(current->data,"%80s%255s",key,valstr);
05282
05283 MGridforceParams* mgfp = NULL;
05284 mgfp = mgridforcelist.find_key(key);
05285 if ( mgfp != NULL) {
05286 iout << iINFO << "MGRIDFORCEPOTFILE key "
05287 << key << " redefined for file " << valstr << "\n" << endi;
05288 } else {
05289 mgfp = mgridforcelist.add(key);
05290 }
05291 int fnamelen = strlen(valstr);
05292 mgfp->gridforceVfile = new char[fnamelen+1];
05293 strncpy(mgfp->gridforceVfile,valstr,fnamelen+1);
05294 mgfp->gridforceScale.x =
05295 mgfp->gridforceScale.y =
05296 mgfp->gridforceScale.z = 1.;
05297 mgfp->gridforceVOffset.x =
05298 mgfp->gridforceVOffset.y =
05299 mgfp->gridforceVOffset.z = 0.;
05300
05301 current = current->next;
05302 }
05303
05304 current = config->find("mgridforcefile");
05305 while (current != NULL) {
05306 int curlen = strlen(current->data);
05307
05308
05309 sscanf(current->data,"%80s%255s",key,valstr);
05310
05311 MGridforceParams* mgfp = NULL;
05312 mgfp = mgridforcelist.find_key(key);
05313 if ( mgfp == NULL) {
05314 iout << iINFO << "MGRIDFORCEFILE no key "
05315 << key << " defined for file " << valstr << "\n" << endi;
05316 } else {
05317 int fnamelen = strlen(valstr);
05318 if (mgfp->gridforceFile != NULL) {
05319 delete [] mgfp->gridforceFile;
05320 }
05321 mgfp->gridforceFile = new char[fnamelen+1];
05322 strncpy(mgfp->gridforceFile,valstr,fnamelen+1);
05323 }
05324
05325 current = current->next;
05326 }
05327
05328 current = config->find("mgridforcevolts");
05329 while (current != NULL) {
05330
05331
05332 int curlen = strlen(current->data);
05333 sscanf(current->data,"%80s%255s",key,valstr);
05334
05335 MGridforceParams* mgfp = NULL;
05336 mgfp = mgridforcelist.find_key(key);
05337 if ( mgfp == NULL) {
05338 iout << iINFO << "MGRIDFORCEVOLTS no key "
05339 << key << " defined for file " << valstr << "\n" << endi;
05340 } else {
05341 int boolval = MGridforceParamsList::atoBool(valstr);
05342 if (boolval == -1) {
05343 iout << iINFO << "MGRIDFORCEVOLTS key "
05344 << key << " boolval " << valstr << " badly defined" << endi;
05345 } else {
05346 mgfp->gridforceVolts = (boolval == 1);
05347 }
05348 }
05349
05350 current = current->next;
05351 }
05352
05353 current = config->find("mgridforcescale");
05354 while (current != NULL) {
05355
05356
05357 int curlen = strlen(current->data);
05358 int nread;
05359 sscanf(current->data,"%80s%n",key,&nread);
05360 char *val = current->data + nread + 1;
05361
05362 MGridforceParams* mgfp = NULL;
05363 mgfp = mgridforcelist.find_key(key);
05364 if ( mgfp == NULL) {
05365 iout << iINFO << "MGRIDFORCESCALE no key "
05366 << key << " defined for vector " << val << "\n" << endi;
05367 } else {
05368 mgfp->gridforceScale.set(val);
05369 }
05370
05371 current = current->next;
05372 }
05373
05374 current = config->find("mgridforcevoff");
05375 while (current != NULL) {
05376
05377
05378 int curlen = strlen(current->data);
05379 int nread;
05380 sscanf(current->data,"%80s%n",key,&nread);
05381 char *val = current->data + nread + 1;
05382
05383 MGridforceParams* mgfp = NULL;
05384 mgfp = mgridforcelist.find_key(key);
05385 if ( mgfp == NULL) {
05386 iout << iINFO << "MGRIDFORCEVOFF no key "
05387 << key << " defined for vector " << val << "\n" << endi;
05388 } else {
05389 mgfp->gridforceVOffset.set(val);
05390 }
05391
05392 current = current->next;
05393 }
05394
05395 current = config->find("mgridforcecol");
05396 while (current != NULL) {
05397
05398
05399 int curlen = strlen(current->data);
05400 sscanf(current->data,"%80s%255s",key,valstr);
05401
05402 MGridforceParams* mgfp = NULL;
05403 mgfp = mgridforcelist.find_key(key);
05404 if ( mgfp == NULL) {
05405 iout << iINFO << "MGRIDFORCECOL no key "
05406 << key << " defined for file " << valstr << "\n" << endi;
05407 } else {
05408 int collen = strlen(valstr);
05409 if (mgfp->gridforceCol != NULL) {
05410 delete [] mgfp->gridforceCol;
05411 }
05412 mgfp->gridforceCol = new char[collen+1];
05413 strncpy(mgfp->gridforceCol,valstr,collen+1);
05414 }
05415
05416 current = current->next;
05417 }
05418
05419 current = config->find("mgridforcechargecol");
05420 while (current != NULL) {
05421
05422
05423 int curlen = strlen(current->data);
05424 sscanf(current->data,"%80s%255s",key,valstr);
05425
05426 MGridforceParams* mgfp = NULL;
05427 mgfp = mgridforcelist.find_key(key);
05428 if ( mgfp == NULL) {
05429 iout << iINFO << "MGRIDFORCECHARGECOL no key "
05430 << key << " defined for file " << valstr << "\n" << endi;
05431 } else {
05432 int collen = strlen(valstr);
05433 if (mgfp->gridforceQcol != NULL) {
05434 delete [] mgfp->gridforceQcol;
05435 }
05436 mgfp->gridforceQcol = new char[collen+1];
05437 strncpy(mgfp->gridforceQcol,valstr,collen+1);
05438 }
05439
05440 current = current->next;
05441 }
05442
05443 current = config->find("mgridforcecont1");
05444 while (current != NULL) {
05445
05446
05447 int curlen = strlen(current->data);
05448 sscanf(current->data,"%80s%255s",key,valstr);
05449
05450 MGridforceParams* mgfp = NULL;
05451 mgfp = mgridforcelist.find_key(key);
05452 if ( mgfp == NULL) {
05453 iout << iINFO << "MGRIDFORCECONT1 no key "
05454 << key << " defined for file " << valstr << "\n" << endi;
05455 } else {
05456 int boolval = MGridforceParamsList::atoBool(valstr);
05457 if (boolval == -1) {
05458 iout << iINFO << "MGRIDFORCECONT1 key "
05459 << key << " boolval " << valstr << " badly defined" << endi;
05460 } else {
05461 mgfp->gridforceCont[0] = (boolval == 1);
05462 }
05463 }
05464
05465 current = current->next;
05466 }
05467
05468 current = config->find("mgridforcecont2");
05469 while (current != NULL) {
05470
05471
05472 int curlen = strlen(current->data);
05473 sscanf(current->data,"%80s%255s",key,valstr);
05474
05475 MGridforceParams* mgfp = NULL;
05476 mgfp = mgridforcelist.find_key(key);
05477 if ( mgfp == NULL) {
05478 iout << iINFO << "MGRIDFORCECONT2 no key "
05479 << key << " defined for file " << valstr << "\n" << endi;
05480 } else {
05481 int boolval = MGridforceParamsList::atoBool(valstr);
05482 if (boolval == -1) {
05483 iout << iINFO << "MGRIDFORCECONT2 key "
05484 << key << " boolval " << valstr << " badly defined" << endi;
05485 } else {
05486 mgfp->gridforceCont[1] = (boolval == 1);
05487 }
05488 }
05489
05490 current = current->next;
05491 }
05492 current = config->find("mgridforcecont3");
05493 while (current != NULL) {
05494
05495
05496 int curlen = strlen(current->data);
05497 sscanf(current->data,"%80s%255s",key,valstr);
05498
05499 MGridforceParams* mgfp = NULL;
05500 mgfp = mgridforcelist.find_key(key);
05501 if ( mgfp == NULL) {
05502 iout << iINFO << "MGRIDFORCECONT3 no key "
05503 << key << " defined for file " << valstr << "\n" << endi;
05504 NAMD_die("MGRIDFORCE error");
05505 } else {
05506 int boolval = MGridforceParamsList::atoBool(valstr);
05507 if (boolval == -1) {
05508 iout << iINFO << "MGRIDFORCECONT3 key "
05509 << key << " boolval " << valstr << " badly defined" << endi;
05510 } else {
05511 mgfp->gridforceCont[2] = (boolval == 1);
05512 }
05513 }
05514
05515 current = current->next;
05516 }
05517
05518 current = config->find("mgridforcelite");
05519 while (current != NULL) {
05520
05521
05522 int curlen = strlen(current->data);
05523 sscanf(current->data,"%80s%255s",key,valstr);
05524
05525 MGridforceParams* mgfp = NULL;
05526 mgfp = mgridforcelist.find_key(key);
05527 if ( mgfp == NULL) {
05528 iout << iINFO << "MGRIDFORCELITE no key "
05529 << key << " defined for file " << valstr << "\n" << endi;
05530 } else {
05531 int boolval = MGridforceParamsList::atoBool(valstr);
05532 if (boolval == -1) {
05533 iout << iINFO << "MGRIDFORCELITE key "
05534 << key << " boolval " << valstr << " badly defined" << endi;
05535 } else {
05536 mgfp->gridforceLite = (boolval == 1);
05537 }
05538 }
05539
05540 current = current->next;
05541 }
05542
05543 delete [] valstr;
05544 delete [] key;
05545
05546
05547
05548 MGridforceParams* params = mgridforcelist.get_first();
05549
05550 while (params != NULL) {
05551 if (params->gridforceFile == NULL) {
05552 char errmsg[255];
05553 sprintf(errmsg,"Value undefined for gridforceFile for key %s\n",
05554 params->gridforceKey);
05555 NAMD_die(errmsg);
05556 }
05557 if (params->gridforceCol == NULL) {
05558 char errmsg[255];
05559 sprintf(errmsg,"Value undefined for gridforceCol for key %s\n",
05560 params->gridforceKey);
05561 NAMD_die(errmsg);
05562 }
05563 if (params->gridforceLite) {
05564 for (int i = 0; i < 3; i++) {
05565 if (params->gridforceCont[i]) {
05566 char errmsg[255];
05567 sprintf(errmsg, "Cannot use gridforceCont and gridforceLite together for key %s\n",
05568 params->gridforceKey);
05569 NAMD_die(errmsg);
05570 }
05571 }
05572 }
05573 params = params->next;
05574 }
05575
05576 }
05577
05578 void SimParameters::parse_mgrid_string_param(ConfigList *cl,
05579 const char *fieldname,
05580 char **dest)
05581 {
05582 StringList *vallist = cl->find(fieldname);
05583 char *val = NULL;
05584
05585 if (vallist != NULL) {
05586 val = vallist->data;
05587 } else {
05588 return;
05589 }
05590
05591 int len = 0;
05592 if (val == NULL) {
05593 *dest = NULL;
05594 } else {
05595 len = strlen(val);
05596 if (len == 0) {
05597 *dest = NULL;
05598 } else {
05599 *dest = new char[len+1];
05600 strncpy(*dest,val,len+1);
05601 }
05602 }
05603 }
05604
05605
05606
05607 void SimParameters::create_output_directories(const char *dirname){
05608
05609
05610
05611
05612 int baselen = strlen(outputFilename);
05613 char *filename = new char[baselen+32];
05614 memset(filename, 0, baselen+32);
05615 strcpy(filename, outputFilename);
05616 struct stat st;
05617 if(stat(filename, &st)!=0) {
05618 int ret = MKDIR(filename);
05619 if(ret!=0) {
05620 char errmsg[512];
05621 sprintf(errmsg, "Error in creating top-level directory %s!", filename);
05622 NAMD_die(errmsg);
05623 }
05624 }
05625
05626
05627 strcat(filename, PATHSEPSTR);
05628 strcat(filename, dirname);
05629
05630 if(stat(filename, &st)!=0) {
05631 int ret = MKDIR(filename);
05632 if(ret!=0) {
05633 char errmsg[512];
05634 sprintf(errmsg, "Error in creating middle-level directory %s!", filename);
05635 NAMD_die(errmsg);
05636 }
05637 }
05638
05639
05640 char tmpstr[256];
05641 for(int i=0; i<numoutputprocs; i++) {
05642 memset(tmpstr, 0, 256);
05643 sprintf(tmpstr, "%s%s%d", filename, PATHSEPSTR, i);
05644 if(stat(tmpstr, &st)!=0) {
05645 int ret = MKDIR(tmpstr);
05646 if(ret!=0) {
05647 char errmsg[512];
05648 sprintf(errmsg, "Error in creating last-level directory %s!", tmpstr);
05649 NAMD_die(errmsg);
05650 }
05651 }
05652 }
05653 }
05654
05655
05656
05657
05658
05659
05660
05661 #define BoolToString(b) ((b) ? "TRUE" : "FALSE")
05662
05663 void SimParameters::print_mgrid_params()
05664 {
05665 const MGridforceParams* params = mgridforcelist.get_first();
05666
05667 while (params != NULL) {
05668 iout << iINFO << "MGRIDFORCE key " << params->gridforceKey << "\n" << endi;
05669 iout << iINFO << " Potfile " << params->gridforceVfile
05670 << "\n" << endi;
05671 iout << iINFO << " Scale " << params->gridforceScale
05672 << "\n" << endi;
05673 iout << iINFO << " File " << params->gridforceFile
05674 << "\n" << endi;
05675 iout << iINFO << " Col " << params->gridforceCol
05676 << "\n" << endi;
05677
05678 char *qcol_msg = "Use atom charge";
05679 if (params->gridforceQcol != NULL) {
05680 qcol_msg = params->gridforceQcol;
05681 }
05682 iout << iINFO << " ChargeCol " << qcol_msg
05683 << "\n" << endi;
05684 iout << iINFO << " VOffset " << params->gridforceVOffset
05685 << "\n" << endi;
05686 iout << iINFO << " Continuous K1 "
05687 << BoolToString(params->gridforceCont[0])
05688 << "\n" << endi;
05689 iout << iINFO << " Continuous K2 "
05690 << BoolToString(params->gridforceCont[1])
05691 << "\n" << endi;
05692 iout << iINFO << " Continuous K3 "
05693 << BoolToString(params->gridforceCont[2])
05694 << "\n" << endi;
05695 iout << iINFO << " Volts "
05696 << BoolToString(params->gridforceVolts)
05697 << "\n" << endi;
05698 iout << iINFO << " Gridforce-Lite "
05699 << BoolToString(params->gridforceLite)
05700 << "\n" << endi;
05701 params = params->next;
05702 }
05703 }
05704
05705
05706
05707
05708
05709
05710
05711
05712
05713
05714
05715
05716
05717
05718 void SimParameters::send_SimParameters(MOStream *msg)
05719
05720 {
05721
05722
05723
05724
05725
05726
05727 msg->put(sizeof(SimParameters),(char*)this);
05728 if ( FFTWWisdomString ) {
05729 int fftwlen = strlen(FFTWWisdomString) + 1;
05730 msg->put(fftwlen);
05731 msg->put(fftwlen,FFTWWisdomString);
05732 }
05733 if ( tclBCScript ) {
05734 int tcllen = strlen(tclBCScript) + 1;
05735 msg->put(tcllen);
05736 msg->put(tcllen,tclBCScript);
05737 }
05738
05739 #ifdef MEM_OPT_VERSION
05740 int filelen = strlen(binAtomFile)+1;
05741 msg->put(filelen);
05742 msg->put(filelen, binAtomFile);
05743
05744 filelen = strlen(binCoorFile)+1;
05745 msg->put(filelen);
05746 msg->put(filelen, binCoorFile);
05747
05748 if(binVelFile) {
05749 filelen = strlen(binVelFile)+1;
05750 msg->put(filelen);
05751 msg->put(filelen, binVelFile);
05752 }
05753 #endif
05754
05755 mgridforcelist.pack_data(msg);
05756
05757 msg->end();
05758 }
05759
05760
05761
05762
05763
05764
05765
05766
05767
05768
05769
05770 void SimParameters::receive_SimParameters(MIStream *msg)
05771
05772 {
05773 msg->get(sizeof(SimParameters),(char*)this);
05774 if ( FFTWWisdomString ) {
05775 int fftwlen;
05776 msg->get(fftwlen);
05777 FFTWWisdomString = new char[fftwlen];
05778 msg->get(fftwlen,FFTWWisdomString);
05779 #ifdef NAMD_FFTW
05780 #ifdef NAMD_FFTW_3
05781 fftwf_import_wisdom_from_string(FFTWWisdomString);
05782 #else
05783 fftw_import_wisdom_from_string(FFTWWisdomString);
05784 #endif
05785 #endif
05786 }
05787 if ( tclBCScript ) {
05788 int tcllen;
05789 msg->get(tcllen);
05790 tclBCScript = new char[tcllen];
05791 msg->get(tcllen,tclBCScript);
05792 }
05793
05794 #ifdef MEM_OPT_VERSION
05795 int filelen;
05796 msg->get(filelen);
05797 binAtomFile = new char[filelen];
05798 msg->get(filelen, binAtomFile);
05799
05800 msg->get(filelen);
05801 binCoorFile = new char[filelen];
05802 msg->get(filelen, binCoorFile);
05803
05804 if(binVelFile) {
05805 msg->get(filelen);
05806 binVelFile = new char[filelen];
05807 msg->get(filelen, binVelFile);
05808 }
05809 #endif
05810
05811
05812
05813
05814 mgridforcelist.clear();
05815 mgridforcelist.unpack_data(msg);
05816
05817 delete msg;
05818 }
05819
05820
