APBS 3.0.0
Loading...
Searching...
No Matches
main.c
Go to the documentation of this file.
1
68#include <time.h>
69
70#include "routines.h"
71
72VEMBED(rcsid="$Id$")
73
74
80int main(
81 int argc,
82 char **argv
83 )
84{
85 // PCE: Adding below variables temporarily
86 clock_t ts, te;
87 // End PCE
88
89 NOsh *nosh = VNULL;
90
91 MGparm *mgparm = VNULL;
92 FEMparm *feparm = VNULL;
93#ifdef ENABLE_BEM
94 BEMparm *bemparm = VNULL;
95#endif
96 GEOFLOWparm *geoflowparm = VNULL;
97 PBEparm *pbeparm = VNULL;
98 APOLparm *apolparm = VNULL;
99 Vparam *param = VNULL;
100#if defined(ENABLE_PBAM) || defined(ENABLE_PBSAM)
101 PBAMparm *pbamparm = VNULL;
102#endif
103
104#ifdef ENABLE_PBSAM
105 PBSAMparm *pbsamparm = VNULL;
106#endif
107
108 Vmem *mem = VNULL;
109 Vcom *com = VNULL;
110 Vio *sock = VNULL;
111#ifdef HAVE_MC_H
112 Vfetk *fetk[NOSH_MAXCALC];
113 Gem *gm[NOSH_MAXMOL];
114 int isolve;
115#else
116 void *fetk[NOSH_MAXCALC];
117 void *gm[NOSH_MAXMOL];
118#endif
119 Vpmg *pmg[NOSH_MAXCALC];
120 Vpmgp *pmgp[NOSH_MAXCALC];
121 Vpbe *pbe[NOSH_MAXCALC];
122 Valist *alist[NOSH_MAXMOL];
123 Vgrid *dielXMap[NOSH_MAXMOL],
124 *dielYMap[NOSH_MAXMOL],
125 *dielZMap[NOSH_MAXMOL],
126 *kappaMap[NOSH_MAXMOL],
127 *potMap[NOSH_MAXMOL],
128 *chargeMap[NOSH_MAXMOL];
129 char *input_path = VNULL,
130 *output_path = VNULL;
131 int i,
132 rank, // proc id
133 size, // total num of procs
134 k;
135 size_t bytesTotal,
136 highWater;
137 Voutput_Format outputformat;
138
139 int rc = 0;
140
141 /* The energy double arrays below store energies from various calculations. */
142 double qfEnergy[NOSH_MAXCALC],
143 qmEnergy[NOSH_MAXCALC];
144 double dielEnergy[NOSH_MAXCALC],
145 totEnergy[NOSH_MAXCALC];
146 double *atomEnergy[NOSH_MAXCALC];
147 AtomForce *atomForce[NOSH_MAXCALC]; /* Stores forces from various calculations. */
148 int nenergy[NOSH_MAXCALC], /* Stores either a flag (0,1) displaying whether
149 * energies were calculated, or, if PCE_COMPS
150 * was used, the number of atom energies stored
151 * for the given calculation. */
152 nforce[NOSH_MAXCALC]; /* Stores an integer which either says no
153 * calculation was performed (0) or gives the
154 * number of entries in the force array for each
155 * calculation. */
156
157 /* The real partition centers */
158 double realCenter[3];
159
160 /* Instructions: */
161 char header[] = {"\n\n\
162----------------------------------------------------------------------\n\
163 APBS -- Adaptive Poisson-Boltzmann Solver\n\
164 Version " PACKAGE_STRING "\n\
165 \n\
166 Nathan A. Baker (nathan.baker@pnnl.gov)\n\
167 Pacific Northwest National Laboratory\n\
168 \n\
169 Additional contributing authors listed in the code documentation.\n\
170 \n\
171 Copyright (c) 2010-2020 Battelle Memorial Institute. Developed at the Pacific\n\
172 Northwest National Laboratory, operated by Battelle Memorial Institute, Pacific\n\
173 Northwest Division for the U.S. Department of Energy.\n\
174 \n\
175 Portions Copyright (c) 2002-2010, Washington University in St. Louis.\n\
176 Portions Copyright (c) 2002-2020, Nathan A. Baker.\n\
177 Portions Copyright (c) 1999-2002, The Regents of the University of California.\n\
178 Portions Copyright (c) 1995, Michael Holst.\n\
179 All rights reserved.\n\
180 \n\
181 Redistribution and use in source and binary forms, with or without\n\
182 modification, are permitted provided that the following conditions are met:\n\
183 \n\
184 * Redistributions of source code must retain the above copyright notice, this\n\
185 list of conditions and the following disclaimer.\n\
186 \n\
187 * Redistributions in binary form must reproduce the above copyright notice,\n\
188 this list of conditions and the following disclaimer in the documentation\n\
189 and/or other materials provided with the distribution.\n\
190 \n\
191 * Neither the name of the developer nor the names of its contributors may be\n\
192 used to endorse or promote products derived from this software without\n\
193 specific prior written permission.\n\
194 \n\
195 THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS \"AS IS\" AND\n\
196 ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED\n\
197 WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE\n\
198 DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR\n\
199 ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES\n\
200 (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;\n\
201 LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND\n\
202 ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT\n\
203 (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS\n\
204 SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.\n\
205----------------------------------------------------------------------\n\
206 APBS uses FETK (the Finite Element ToolKit) to solve the\n\
207 Poisson-Boltzmann equation numerically. FETK is a portable collection\n\
208 of finite element modeling class libraries developed by the Michael Holst\n\
209 research group and written in an object-oriented form of C. FEtk is\n\
210 designed to solve general coupled systems of nonlinear partial differential\n\
211 equations using adaptive finite element methods, inexact Newton methods,\n\
212 and algebraic multilevel methods. More information about FEtk may be found\n\
213 at <http://www.FEtk.ORG>.\n\
214----------------------------------------------------------------------\n\
215 APBS also uses Aqua to solve the Poisson-Boltzmann equation numerically. \n\
216 Aqua is a modified form of the Holst group PMG library <http://www.FEtk.ORG>\n\
217 which has been modified by Patrice Koehl\n\
218 <http://koehllab.genomecenter.ucdavis.edu/> for improved efficiency and\n\
219 memory usage when solving the Poisson-Boltzmann equation.\n\
220----------------------------------------------------------------------\n\
221 Please cite your use of APBS as:\n\n\
222 Baker NA, Sept D, Joseph S, Holst MJ, McCammon JA. Electrostatics of\n\
223 nanosystems: application to microtubules and the ribosome. Proc.\n\
224 Natl. Acad. Sci. USA 98, 10037-10041 2001.\n\
225 \n\n"};
226 char *usage =
227{"\n\n\
228----------------------------------------------------------------------\n\
229 This driver program calculates electrostatic potentials, energies,\n\
230 and forces using both multigrid and finite element methods.\n\
231 It is invoked as:\n\n\
232 apbs [options] apbs.in\n\n\
233 where apbs.in is a formatted input file and [options] are:\n\n\
234--output-file=<name> Enables output logging to the path\n\
235 listed in <name>. Uses flat-file\n\
236 format is --output-format is not used.\n\
237--output-format=<type> Specifies format for logging. Options\n\
238 for type are either \"xml\" or \"flat\".\n\
239--help Display this help information.\n\
240--version Display the current APBS version.\n\
241----------------------------------------------------------------------\n\n"};
242
243 /* ************** CHECK PARALLEL STATUS *************** */
244 VASSERT(Vcom_init(&argc, &argv));
245 com = Vcom_ctor(1);
246 rank = Vcom_rank(com);
247 size = Vcom_size(com);
248 startVio();
249 Vnm_setIoTag(rank, size);
250 Vnm_tprint( 0, "Hello world from PE %d\n", rank);
251
252 /* A bit of array/pointer initialization */
253 mem = Vmem_ctor("MAIN");
254 for (i=0; i<NOSH_MAXCALC; i++) {
255 pmg[i] = VNULL;
256 pmgp[i] = VNULL;
257 fetk[i] = VNULL;
258 pbe[i] = VNULL;
259 qfEnergy[i] = 0;
260 qmEnergy[i] = 0;
261 dielEnergy[i] = 0;
262 totEnergy[i] = 0;
263 atomForce[i] = VNULL;
264 nenergy[i] = 0;
265 nforce[i] = 0;
266 }
267 for (i=0; i<NOSH_MAXMOL; i++) {
268 alist[i] = VNULL;
269 dielXMap[i] = VNULL;
270 dielYMap[i] = VNULL;
271 dielZMap[i] = VNULL;
272 kappaMap[i] = VNULL;
273 potMap[i] = VNULL;
274 chargeMap[i] = VNULL;
275 }
276
277 /* ********* CHECK INVOCATION AND OPTIONS ************* */
278 Vnm_tstart(APBS_TIMER_WALL_CLOCK, "APBS WALL CLOCK");
279 Vnm_tprint( 1, "%s", header);
280
281#ifdef APBS_FAST
282 Vnm_tprint(, 2"WARNING: APBS was compiled with the --enable-fast option.\n"
283 "WARNING: This mode is experimental and subject to change in future releases.\n"
284 "WARNING: The fast mode enables: Gauess-Seidel Smoothing and \n"
285 "WARNING: Conjugate Gradient Multigrid methods.\n\n");
286#endif
287
288 Vnm_tprint( 1, "This executable compiled on %s at %s\n\n", __DATE__, __TIME__);
289
290#if defined(WITH_TINKER)
291 Vnm_tprint( 2, "This executable was compiled with TINKER support and is not intended for stand-alone execution.\n");
292 Vnm_tprint( 2, "Please compile another version without TINKER support.\n");
293 exit(2);
294#endif
295
296 /* Process program arguments */
297 i=0;
298 outputformat = OUTPUT_NULL;
299 while (i<argc){
300 if (strncmp(argv[i], "--", 2) == 0) {
301
302 /* Long Options */
303 if (Vstring_strcasecmp("--version", argv[i]) == 0){
304 Vnm_tprint(2, "%s\n", PACKAGE_STRING);
305 VJMPERR1(0);
306 } else if (Vstring_strcasecmp("--help", argv[i]) == 0){
307 Vnm_tprint(2, "%s\n", usage);
308 VJMPERR1(0);
309 } else if (strncmp(argv[i], "--output-format", 15) == 0) {
310 if (strstr(argv[i], "xml") != NULL) {
311 Vnm_tprint(2, "XML output format is now deprecated, please use --output-format=flat instead!\n\n");
312 VJMPERR1(0);
313 }
314 else if (strstr(argv[i], "flat") != NULL) {
315 outputformat = OUTPUT_FLAT;
316 } else {
317 Vnm_tprint(2, "Invalid output-format type!\n");
318 VJMPERR1(0);
319 }
320 } else if (strncmp(argv[i], "--output-file=", 14) == 0){
321 output_path = strstr(argv[i], "=");
322 ++output_path;
323 if (outputformat == OUTPUT_NULL) outputformat = OUTPUT_FLAT;
324 } else {
325 Vnm_tprint(2, "UNRECOGNIZED COMMAND LINE OPTION %s!\n", argv[i]);
326 Vnm_tprint(2, "%s\n", usage);
327 VJMPERR1(0);
328 }
329 } else {
330
331 /* Set the path to the input file */
332 if ((input_path == VNULL) && (i != 0))
333 input_path = argv[i];
334 else if (i != 0) {
335 Vnm_tprint(2, "ERROR -- CALLED WITH TOO MANY ARGUMENTS!\n", \
336 argc);
337 Vnm_tprint(2, "%s\n", usage);
338 VJMPERR1(0);
339 }
340 }
341 i++;
342 }
343
344 /* If we set an output format but no path, error. */
345 if ((outputformat != 0) && (output_path == NULL)) {
346 Vnm_tprint(2, "The --output-path variable must be set when using --output-format!\n");
347 VJMPERR1(0);
348 }
349
350 /* If we failed to specify an input file, error. */
351 if (input_path == NULL) {
352 Vnm_tprint(2, "ERROR -- APBS input file not specified!\n", argc);
353 Vnm_tprint(2, "%s\n", usage);
354 VJMPERR1(0);
355 }
356
357 /* Append rank info if a parallel run */
358 if ((size > 1) && (output_path != NULL))
359 printf(output_path, "%s_%d", output_path, rank);
360
361 /* *************** PARSE INPUT FILE ******************* */
362 nosh = NOsh_ctor(rank, size);
363 Vnm_tprint( 1, "Parsing input file %s...\n", input_path);
364 Vnm_tprint( 1, "rank %d size %d...\n", rank, size);
365 sock = Vio_ctor("FILE", "ASC", VNULL, input_path, "r");
366 if (sock == VNULL) {
367 Vnm_tprint(2, "Error while opening input file %s!\n", input_path);
368 VJMPERR1(0);
369 }
370 if (!NOsh_parseInput(nosh, sock)) {
371 Vnm_tprint( 2, "Error while parsing input file.\n");
372 VJMPERR1(0);
373 }
374 else
375 Vnm_tprint( 1, "Parsed input file.\n");
376 Vio_dtor(&sock);
377
378 /* *************** LOAD PARAMETERS AND MOLECULES ******************* */
379 param = loadParameter(nosh);
380 if (loadMolecules(nosh, param, alist) != 1) {
381 Vnm_tprint(2, "Error reading molecules!\n");
382 VJMPERR1(0);
383 }
384
385 /* *************** SETUP CALCULATIONS *************** */
386 if (NOsh_setupElecCalc(nosh, alist) != 1) {
387 Vnm_tprint(2, "Error setting up ELEC calculations\n");
388 VJMPERR1(0);
389 }
390
391 if ((rc = NOsh_setupApolCalc(nosh, alist)) == ACD_ERROR) {
392 Vnm_tprint(2, "Error setting up APOL calculations\n");
393 VJMPERR1(0);
394 }
395
396 /* ******************* CHECK APOL********************** */
397 /* if((nosh->gotparm == 0) && (rc == ACD_YES)){
398 Vnm_print(1,"\nError you must provide a parameter file if you\n" \
399 " are performing an APOLAR calculation\n");
400 VJMPERR1(0);
401 } */
402
403#if defined(DEBUG_MAC_OSX_OCL)
404#include "mach_chud.h"
405#include <stdint.h>
406 uint64_t mbeg;
407 machm_(&mbeg);
408
409 if(clFinish != NULL)
410 {
411 int ret = initOpenCL();
412 printf("OpenCL runtime present - initialized = %i\n",ret);
413 }
414 else
415 {
416 setkOpenCLAvailable_(0);
417 printf("OpenCL is not present!\n");
418 }
419#endif
420
421#if defined(DEBUG_MAC_OSX_STANDARD)
422#include "mach_chud.h"
423#include <stdint.h>
424 uint64_t mbeg;
425 machm_(&mbeg);
426#endif
427
428 /* *************** LOAD MAPS ******************* */
429 if (loadDielMaps(nosh, dielXMap, dielYMap, dielZMap) != 1) {
430 Vnm_tprint(2, "Error reading dielectric maps!\n");
431 VJMPERR1(0);
432 }
433 if (loadKappaMaps(nosh, kappaMap) != 1) {
434 Vnm_tprint(2, "Error reading kappa maps!\n");
435 VJMPERR1(0);
436 }
437 if (loadPotMaps(nosh, potMap) != 1) {
438 Vnm_tprint(2, "Error reading potential maps!\n");
439 VJMPERR1(0);
440 }
441 if (loadChargeMaps(nosh, chargeMap) != 1) {
442 Vnm_tprint(2, "Error reading charge maps!\n");
443 VJMPERR1(0);
444 }
445
446 /* *************** DO THE CALCULATIONS ******************* */
447 Vnm_tprint( 1, "Preparing to run %d PBE calculations.\n",
448 nosh->ncalc);
449 for (i=0; i<nosh->ncalc; i++) {
450 Vnm_tprint( 1, "----------------------------------------\n");
451
452 switch (nosh->calc[i]->calctype) {
453 /* Multigrid */
454 case NCT_MG:
455 /* What is this? This seems like a very awkward way to find
456 the right ELEC statement... */
457 for (k=0; k<nosh->nelec; k++) {
458 if (nosh->elec2calc[k] >= i) {
459 break;
460 }
461 }
462 if (Vstring_strcasecmp(nosh->elecname[k], "") == 0) {
463 Vnm_tprint( 1, "CALCULATION #%d: MULTIGRID\n", i+1);
464 } else {
465 Vnm_tprint( 1, "CALCULATION #%d (%s): MULTIGRID\n",
466 i+1, nosh->elecname[k]);
467 }
468 /* Useful local variables */
469 mgparm = nosh->calc[i]->mgparm;
470 pbeparm = nosh->calc[i]->pbeparm;
471
472 /* Set up problem */
473 Vnm_tprint( 1, " Setting up problem...\n");
474
475 if (!initMG(i, nosh, mgparm, pbeparm, realCenter, pbe,
476 alist, dielXMap, dielYMap, dielZMap, kappaMap,
477 chargeMap, pmgp, pmg, potMap)) {
478 Vnm_tprint( 2, "Error setting up MG calculation!\n");
479 VJMPERR1(0);
480 }
481
482 /* Print problem parameters */
483 printMGPARM(mgparm, realCenter);
484 printPBEPARM(pbeparm);
485
486 /* Solve PDE */
487 if (solveMG(nosh, pmg[i], mgparm->type) != 1) {
488 Vnm_tprint(2, "Error solving PDE!\n");
489 VJMPERR1(0);
490 }
491
492 /* Set partition information for observables and I/O */
493 if (setPartMG(nosh, mgparm, pmg[i]) != 1) {
494 Vnm_tprint(2, "Error setting partition info!\n");
495 VJMPERR1(0);
496 }
497
498 /* Write out energies */
499 energyMG(nosh, i, pmg[i],
500 &(nenergy[i]), &(totEnergy[i]), &(qfEnergy[i]),
501 &(qmEnergy[i]), &(dielEnergy[i]));
502
503 /* Write out forces */
504 forceMG(mem, nosh, pbeparm, mgparm, pmg[i], &(nforce[i]),
505 &(atomForce[i]), alist);
506
507 /* Write out data folks might want */
508 writedataMG(rank, nosh, pbeparm, pmg[i]);
509
510 /* Write matrix */
511 writematMG(rank, nosh, pbeparm, pmg[i]);
512
513 /* If needed, cache atom energies */
514 nenergy[i] = 0;
515 if ((pbeparm->calcenergy == PCE_COMPS) && (outputformat != OUTPUT_NULL)){
516 storeAtomEnergy(pmg[i], i, &(atomEnergy[i]), &(nenergy[i]));
517 }
518
519 fflush(stdout);
520 fflush(stderr);
521
522 break;
523
524 /* ***** Do FEM calculation ***** */
525 case NCT_FEM:
526#ifdef HAVE_MC_H
527 for (k=0; k<nosh->nelec; k++) {
528 if (nosh->elec2calc[k] >= i) break;
529 }
530 if (Vstring_strcasecmp(nosh->elecname[i+1], "") == 0) {
531 Vnm_tprint( 1, "CALCULATION #%d: FINITE ELEMENT\n", i+1);
532 } else {
533 Vnm_tprint( 1, "CALCULATION #%d (%s): FINITE ELEMENT\n", i+1, nosh->elecname[k+1]);
534 }
535
536 /* Useful local variables */
537 feparm = nosh->calc[i]->femparm;
538 pbeparm = nosh->calc[i]->pbeparm;
539
540 /* Warn the user about some things */
541 Vnm_tprint(2, "#################### WARNING ###################\n");
542 Vnm_tprint(2, "## FE support is currently very experimental! ##\n");
543 Vnm_tprint(2, "#################### WARNING ###################\n");
544
545 /* Set up problem */
546 Vnm_tprint( 1, " Setting up problem...\n");
547 /* Attempt to initialize and do an initial refinement of the mesh data. The mesh data
548 * will be stored in the Vfetk object fetk, which contains the appropriate geometry
549 * manager (Gem) object and Vcsm object describing the mesh structure. The mesh will
550 * either be loaded from an external source or generated from scratch. */
551 if (initFE(i, nosh, feparm, pbeparm, pbe, alist, fetk) != VRC_SUCCESS) {
552 Vnm_tprint( 2, "Error setting up FE calculation!\n");
553 VJMPERR1(0);
554 }
555
556 /* Print problem parameters */
557 printFEPARM(i, nosh, feparm, fetk);
558 printPBEPARM(pbeparm);
559
560 /* Refine mesh - this continues to run the AM_markRefine procedure already run
561 * in initFE() to arrive at some initial refinement, but does checks of the
562 * simplices so that it refines until the error or size tolerances are reached.
563 * Once this is done, we have a mesh that has been refined to the point where
564 * we can attempt to solve - further refinement may be needed in the loop
565 * below. */
566 if (!preRefineFE(i, feparm, fetk)) {
567 Vnm_tprint( 2, "Error pre-refining mesh!\n");
568 VJMPERR1(0);
569 }
570
571 /* Solve-estimate-refine */
572 Vnm_tprint(2, "\n\nWARNING! DO NOT EXPECT PERFORMANCE OUT OF THE APBS/FEtk\n");
573 Vnm_tprint(2, "INTERFACE AT THIS TIME. THE FINITE ELEMENT SOLVER IS\n");
574 Vnm_tprint(2, "CURRENTLY NOT OPTIMIZED FOR THE PB EQUATION. IF YOU WANT\n");
575 Vnm_tprint(2, "PERFORMANCE, PLEASE USE THE MULTIGRID-BASED METHODS, E.G.\n");
576 Vnm_tprint(2, "MG-AUTO, MG-PARA, and MG-MANUAL (SEE DOCS.)\n\n");
577 Vnm_tprint(1, " Beginning solve-estimate-refine cycle:\n");
578
579 for (isolve=0; isolve<feparm->maxsolve; isolve++) {
580 Vnm_tprint(1, " Solve #%d...\n", isolve);
581
582 /* Attempt to solve the mesh by using one of MC's solver types. */
583 if (!solveFE(i, pbeparm, feparm, fetk)) {
584 Vnm_tprint(2, "ERROR SOLVING EQUATION!\n");
585 VJMPERR1(0);
586 }
587
588 /* Calculate the total electrostatic energy. */
589 if (!energyFE(nosh, i, fetk, &(nenergy[i]),
590 &(totEnergy[i]), &(qfEnergy[i]),
591 &(qmEnergy[i]), &(dielEnergy[i]))) {
592 Vnm_tprint(2, "ERROR SOLVING EQUATION!\n");
593 VJMPERR1(0);
594 }
595
596 /* We're not going to refine if we've hit the max number
597 * of solves */
598 if (isolve < (feparm->maxsolve)-1) {
599 /* Do a final error estimation and mesh refinement. */
600 if (!postRefineFE(i, feparm, fetk)) {
601 break;
602 }
603 }
604 bytesTotal = Vmem_bytesTotal();
605 highWater = Vmem_highWaterTotal();
606 Vnm_tprint(1, " Current memory use: %g MB\n",
607 ((double)bytesTotal/(1024.)/(1024.)));
608 Vnm_tprint(1, " High-water memory use: %g MB\n",
609 ((double)highWater/(1024.)/(1024.)));
610 }
611
612 Vnm_tprint(1, " Writing FEM data to files.\n");
613
614 /* Save data. */
615 if (!writedataFE(rank, nosh, pbeparm, fetk[i])) {
616 Vnm_tprint(2, " Error while writing FEM data!\n");
617 }
618#else /* ifdef HAVE_MC_H */
619 Vnm_print(2, "Error! APBS not compiled with FEtk!\n");
620 exit(2);
621#endif /* ifdef HAVE_MC_H */
622 break;
623
624 /* Do an apolar calculation */
625 case NCT_APOL:
626 /* Copied from NCT_MG. See the note above (top of loop) for
627 information about this loop.
628 */
629 for (k=0; k<nosh->napol; k++) {
630 if (nosh->apol2calc[k] >= i) {
631 break;
632 }
633 }
634
635 if (Vstring_strcasecmp(nosh->apolname[k], "") == 0) {
636 Vnm_tprint( 1, "CALCULATION #%d: APOLAR\n", i+1);
637 } else {
638 Vnm_tprint( 1, "CALCULATION #%d (%s): APOLAR\n",
639 i+1, nosh->apolname[k]);
640 }
641
642 apolparm = nosh->calc[i]->apolparm;
643 // poor man's execution timer.
644 ts = clock();
645 rc = initAPOL(nosh, mem, param, apolparm, &(nforce[i]), &(atomForce[i]),
646 alist[(apolparm->molid)-1]);
647 Vnm_print(0, "initAPOL: Time elapsed: %f\n", ((double)clock() - ts) / CLOCKS_PER_SEC);
648 if(rc == 0) {
649 Vnm_tprint(2, "Error calculating apolar solvation quantities!\n");
650 VJMPERR1(0);
651 }
652 break;
653
654 /* Boundary Element (tabi) */
655 case NCT_BEM:
656#ifdef ENABLE_BEM
657 /* What is this? This seems like a very awkward way to find
658 the right ELEC statement... */
659 for (k=0; k<nosh->nelec; k++) {
660 if (nosh->elec2calc[k] >= i) {
661 break;
662 }
663 }
664 if (Vstring_strcasecmp(nosh->elecname[k], "") == 0) {
665 Vnm_tprint( 1, "CALCULATION #%d: BOUNDARY ELEMENT\n", i+1);
666 } else {
667 Vnm_tprint( 1, "CALCULATION #%d (%s): BOUNDARY ELEMENT\n",
668 i+1, nosh->elecname[k]);
669 }
670 /* Useful local variables */
671 bemparm = nosh->calc[i]->bemparm;
672 pbeparm = nosh->calc[i]->pbeparm;
673
674 /* Set up problem */
675 Vnm_tprint( 1, " Setting up problem...\n");
676
677 if (!initBEM(i,nosh, bemparm, pbeparm, pbe)) {
678 Vnm_tprint( 2, "Error setting up BEM calculation!\n");
679 VJMPERR1(0);
680 }
681
682 /* Print problem parameters */
683 printBEMPARM(bemparm);
684 printPBEPARM(pbeparm);
685
686 /* Solve PDE */
687 if (solveBEM(alist, nosh, pbeparm, bemparm, bemparm->type) != 1) {
688 Vnm_tprint(2, "Error solving PDE!\n");
689 VJMPERR1(0);
690 }
691
692 /* Write out energies */
693 energyBEM(nosh, i,
694 &(nenergy[i]), &(totEnergy[i]), &(qfEnergy[i]),
695 &(qmEnergy[i]), &(dielEnergy[i]));
696
697 /* Write out forces */
698 forceBEM(nosh, pbeparm, bemparm, &(nforce[i]),
699 &(atomForce[i]), alist);
700
701 /* Write out data folks might want */
702 writedataBEM(rank, nosh, pbeparm);
703
704 /* Write matrix */
705 writematBEM(rank, nosh, pbeparm);
706
707 /* If needed, cache atom energies */
708 nenergy[i] = 0;
709 if ((pbeparm->calcenergy == PCE_COMPS) && (outputformat != OUTPUT_NULL)){
710 storeAtomEnergy(pmg[i], i, &(atomEnergy[i]), &(nenergy[i]));
711 }
712
713 fflush(stdout);
714 fflush(stderr);
715#else /* ifdef ENABLE_BEM */
716 Vnm_print(2, "Error! APBS not compiled with BEM!\n");
717 exit(2);
718#endif
719 break;
720
721 /* geometric flow */
722 case NCT_GEOFLOW:
723#ifdef ENABLE_GEOFLOW
724 /* What is this? This seems like a very awkward way to find
725 the right ELEC statement... */
726 for (k=0; k<nosh->nelec; k++) {
727 if (nosh->elec2calc[k] >= i) {
728 break;
729 }
730 }
731 if (Vstring_strcasecmp(nosh->elecname[k], "") == 0) {
732 Vnm_tprint( 1, "CALCULATION #%d: GEOMETRIC FLOW\n", i+1);
733 } else {
734 Vnm_tprint( 1, "CALCULATION #%d (%s): GEOMETRIC FLOW\n",
735 i+1, nosh->elecname[k]);
736 }
737 /* Useful local variables */
738 geoflowparm = nosh->calc[i]->geoflowparm;
739 apolparm = nosh->calc[i]->apolparm;
740 pbeparm = nosh->calc[i]->pbeparm;
741
742 /* Set up problem */
743 Vnm_tprint( 1, " Setting up problem...\n");
744
745
746 /* Solve PDE */
747 if (solveGeometricFlow(alist, nosh, pbeparm, apolparm, geoflowparm) != 1) {
748 Vnm_tprint(2, "Error solving GEOFLOW!\n");
749 VJMPERR1(0);
750 }
751
752 fflush(stdout);
753 fflush(stderr);
754 break;
755#else /* ifdef ENABLE_GEOFLOW*/
756 Vnm_print(2, "Error! APBS not compiled with GEOFLOW!\n");
757 exit(2);
758#endif
759
760 /* Poisson-boltzmann analytical method */
761 case NCT_PBAM:
762#ifdef ENABLE_PBAM
763 /* What is this? This seems like a very awkward way to find
764 the right ELEC statement... */
765 //Vnm_tprint( 1, "Made it to start\n");
766 for (k=0; k<nosh->nelec; k++) {
767 if (nosh->elec2calc[k] >= i) {
768 break;
769 }
770 }
771 if (Vstring_strcasecmp(nosh->elecname[k], "") == 0) {
772 Vnm_tprint( 1, "CALCULATION #%d: PBAM\n", i+1);
773 } else {
774 Vnm_tprint( 1, "CALCULATION #%d (%s): PBAM\n",
775 i+1, nosh->elecname[k]);
776 }
777 /* Useful local variables */
778 pbamparm = nosh->calc[i]->pbamparm;
779 pbeparm = nosh->calc[i]->pbeparm;
780
781 /* Set up problem */
782 Vnm_tprint( 1, " Setting up problem...\n");
783
784 /* Solve LPBE with PBAM method */
785 if (solvePBAM(alist, nosh, pbeparm, pbamparm) != 1) {
786 Vnm_tprint(2, "Error solving PBAM!\n");
787 VJMPERR1(0);
788 }
789
790 fflush(stdout);
791 fflush(stderr);
792 break;
793#else /* ifdef ENABLE_PBAM*/
794 Vnm_print(2, "Error! APBS not compiled with PBAM!\n");
795 exit(2);
796#endif
797
798 case NCT_PBSAM:
799#ifdef ENABLE_PBSAM
800 Vnm_tprint( 1, "Made it to start\n");
801 for (k=0; k<nosh->nelec; k++) {
802 if (nosh->elec2calc[k] >= i) {
803 break;
804 }
805 }
806 if (Vstring_strcasecmp(nosh->elecname[k], "") == 0) {
807 Vnm_tprint( 1, "CALCULATION #%d: PBSAM\n", i+1);
808 } else {
809 Vnm_tprint( 1, "CALCULATION #%d (%s): PBSAM\n",
810 i+1, nosh->elecname[k]);
811 }
812 /* Useful local variables */
813 pbamparm = nosh->calc[i]->pbamparm;
814 pbsamparm = nosh->calc[i]->pbsamparm;
815 pbeparm = nosh->calc[i]->pbeparm;
816
817 /* Set up problem */
818 Vnm_tprint( 1, " Setting up problem...\n");
819
820
821 /* Solve LPBE with PBSAM method */
822 if (solvePBSAM(alist, nosh, pbeparm, pbamparm, pbsamparm) != 1) {
823 Vnm_tprint(2, "Error solving PBSAM!\n");
824 VJMPERR1(0);
825 }
826
827 fflush(stdout);
828 fflush(stderr);
829 break;
830#else /* ifdef ENABLE_PBSAM*/
831 Vnm_print(2, "Error! APBS not compiled with PBSAM!\n");
832 exit(2);
833#endif
834
835 default:
836 Vnm_tprint(2, " Unknown calculation type (%d)!\n", nosh->calc[i]->calctype);
837 exit(2);
838 break;
839 }
840
841 }
842
843 //Clear out the parameter file memory
844 if(param != VNULL) Vparam_dtor(&param);
845
846 /* *************** HANDLE PRINT STATEMENTS ******************* */
847 if (nosh->nprint > 0) {
848 Vnm_tprint( 1, "----------------------------------------\n");
849 Vnm_tprint( 1, "PRINT STATEMENTS\n");
850 }
851 for (i=0; i<nosh->nprint; i++) {
852 /* Print energy */
853 if (nosh->printwhat[i] == NPT_ENERGY) {
854 printEnergy(com, nosh, totEnergy, i);
855 /* Print force */
856 } else if (nosh->printwhat[i] == NPT_FORCE) {
857 printForce(com, nosh, nforce, atomForce, i);
858 } else if (nosh->printwhat[i] == NPT_ELECENERGY) {
859 printElecEnergy(com, nosh, totEnergy, i);
860 } else if (nosh->printwhat[i] == NPT_ELECFORCE) {
861 printElecForce(com, nosh, nforce, atomForce, i);
862 } else if (nosh->printwhat[i] == NPT_APOLENERGY) {
863 printApolEnergy(nosh, i);
864 } else if (nosh->printwhat[i] == NPT_APOLFORCE) {
865 printApolForce(com, nosh, nforce, atomForce, i);
866 } else {
867 Vnm_tprint( 2, "Undefined PRINT keyword!\n");
868 break;
869 }
870 }
871 Vnm_tprint( 1, "----------------------------------------\n");
872
873 /* *************** HANDLE LOGGING *********************** */
874
875 if (outputformat == OUTPUT_FLAT) {
876 Vnm_tprint(2," Writing data to flat file %s...\n\n", output_path);
877 writedataFlat(nosh, com, output_path, totEnergy, qfEnergy, qmEnergy,
878 dielEnergy, nenergy, atomEnergy, nforce, atomForce);
879 }
880
881 /* Destroy energy arrays if they still exist */
882
883 for (i=0; i<nosh->ncalc; i++) {
884 if (nenergy[i] > 0) Vmem_free(mem, nenergy[i], sizeof(double),
885 (void **)&(atomEnergy[i]));
886 }
887
888 /* *************** GARBAGE COLLECTION ******************* */
889
890 Vnm_tprint( 1, "CLEANING UP AND SHUTTING DOWN...\n");
891 /* Clean up APBS structures */
892 killForce(mem, nosh, nforce, atomForce);
893 killEnergy();
894 killMG(nosh, pbe, pmgp, pmg);
895#ifdef HAVE_MC_H
896 killFE(nosh, pbe, fetk, gm);
897#endif
898 killChargeMaps(nosh, chargeMap);
899 killKappaMaps(nosh, kappaMap);
900 killDielMaps(nosh, dielXMap, dielYMap, dielZMap);
901 killMolecules(nosh, alist);
902 NOsh_dtor(&nosh);
903
904 /* Memory statistics */
905 bytesTotal = Vmem_bytesTotal();
906 highWater = Vmem_highWaterTotal();
907 Vnm_tprint( 1, "Final memory usage: %4.3f MB total, %4.3f MB high water\n",
908 (double)(bytesTotal)/(1024.*1024.),
909 (double)(highWater)/(1024.*1024.));
910
911 /* Clean up MALOC structures */
912 Vcom_dtor(&com);
913 Vmem_dtor(&mem);
914
915 /* And now it's time to so "so long"... */
916 Vnm_tprint(1, "\n\n");
917 Vnm_tprint( 1, "Thanks for using APBS!\n\n");
918
919#if defined(DEBUG_MAC_OSX_OCL)
920 mets_(&mbeg, "Main Program CL");
921#endif
922#if defined(DEBUG_MAC_OSX_STANDARD)
923 mets_(&mbeg, "Main Program Standard");
924#endif
925
926 /* This should be last */
927 Vnm_tstop(APBS_TIMER_WALL_CLOCK, "APBS WALL CLOCK");
928 Vnm_flush(1);
929 Vnm_flush(2);
930 Vcom_finalize();
931
932 fflush(NULL);
933
934 return 0;
935
936 VERROR1:
937 Vcom_finalize();
938 Vcom_dtor(&com);
939 Vmem_dtor(&mem);
940 return APBSRC;
941}
@ ACD_ERROR
Definition apolparm.h:114
VPUBLIC void printPBEPARM(PBEparm *pbeparm)
Print out generic PBE params loaded from input.
Definition routines.c:1002
VPUBLIC int loadPotMaps(NOsh *nosh, Vgrid *map[NOSH_MAXMOL])
Load the potential maps given in NOsh into grid objects.
Definition routines.c:793
VPUBLIC void killForce(Vmem *mem, NOsh *nosh, int nforce[NOSH_MAXCALC], AtomForce *atomForce[NOSH_MAXCALC])
Free memory from MG force calculation.
Definition routines.c:1782
VPUBLIC void storeAtomEnergy(Vpmg *pmg, int icalc, double **atomEnergy, int *nenergy)
Store energy in arrays for future use.
Definition routines.c:1870
VPUBLIC int solveMG(NOsh *nosh, Vpmg *pmg, MGparm_CalcType type)
Solve the PBE with MG.
Definition routines.c:1487
VPUBLIC int postRefineFE(int icalc, FEMparm *feparm, Vfetk *fetk[NOSH_MAXCALC])
Estimate error, mark mesh, and refine mesh after solve.
Definition routines.c:4239
VPUBLIC Vrc_Codes initFE(int icalc, NOsh *nosh, FEMparm *feparm, PBEparm *pbeparm, Vpbe *pbe[NOSH_MAXCALC], Valist *alist[NOSH_MAXMOL], Vfetk *fetk[NOSH_MAXCALC])
Initialize FE solver objects.
Definition routines.c:3711
VPUBLIC int preRefineFE(int icalc, FEMparm *feparm, Vfetk *fetk[NOSH_MAXCALC])
Pre-refine mesh before solve.
Definition routines.c:4046
VPUBLIC void startVio()
Wrapper to start MALOC Vio layer.
Definition routines.c:58
VPUBLIC void killKappaMaps(NOsh *nosh, Vgrid *map[NOSH_MAXMOL])
Destroy the loaded kappa maps.
Definition routines.c:776
VPUBLIC int printElecEnergy(Vcom *com, NOsh *nosh, double totEnergy[NOSH_MAXCALC], int iprint)
Combine and pretty-print energy data.
Definition routines.c:2853
int main(int argc, char **argv)
The main APBS function.
Definition main.c:80
VPUBLIC int printForce(Vcom *com, NOsh *nosh, int nforce[NOSH_MAXCALC], AtomForce *atomForce[NOSH_MAXCALC], int iprint)
Combine and pretty-print force data (deprecated...see printElecForce)
Definition routines.c:2980
VPUBLIC int setPartMG(NOsh *nosh, MGparm *mgparm, Vpmg *pmg)
Set MG partitions for calculating observables and performing I/O.
Definition routines.c:1523
VPUBLIC void killDielMaps(NOsh *nosh, Vgrid *dielXMap[NOSH_MAXMOL], Vgrid *dielYMap[NOSH_MAXMOL], Vgrid *dielZMap[NOSH_MAXMOL])
Destroy the loaded dielectric.
Definition routines.c:639
VPUBLIC void killMG(NOsh *nosh, Vpbe *pbe[NOSH_MAXCALC], Vpmgp *pmgp[NOSH_MAXCALC], Vpmg *pmg[NOSH_MAXCALC])
Kill structures initialized during an MG calculation.
Definition routines.c:1461
VPUBLIC int writedataFlat(NOsh *nosh, Vcom *com, const char *fname, double totEnergy[NOSH_MAXCALC], double qfEnergy[NOSH_MAXCALC], double qmEnergy[NOSH_MAXCALC], double dielEnergy[NOSH_MAXCALC], int nenergy[NOSH_MAXCALC], double *atomEnergy[NOSH_MAXCALC], int nforce[NOSH_MAXCALC], AtomForce *atomForce[NOSH_MAXCALC])
Write out information to a flat file.
Definition routines.c:1887
VPUBLIC int solveFE(int icalc, PBEparm *pbeparm, FEMparm *feparm, Vfetk *fetk[NOSH_MAXCALC])
Solve-estimate-refine.
Definition routines.c:4116
VPUBLIC void killEnergy()
Kill arrays allocated for energies.
Definition routines.c:1774
VPUBLIC int loadChargeMaps(NOsh *nosh, Vgrid *map[NOSH_MAXMOL])
Load the charge maps given in NOsh into grid objects.
Definition routines.c:884
VPUBLIC int printElecForce(Vcom *com, NOsh *nosh, int nforce[NOSH_MAXCALC], AtomForce *atomForce[NOSH_MAXCALC], int iprint)
Combine and pretty-print force data.
Definition routines.c:3228
VPUBLIC void killMolecules(NOsh *nosh, Valist *alist[NOSH_MAXMOL])
Destroy the loaded molecules.
Definition routines.c:233
VPUBLIC void printFEPARM(int icalc, NOsh *nosh, FEMparm *feparm, Vfetk *fetk[NOSH_MAXCALC])
Print out FE-specific params loaded from input.
Definition routines.c:3903
VPUBLIC int writedataMG(int rank, NOsh *nosh, PBEparm *pbeparm, Vpmg *pmg)
Write out observables from MG calculation to file.
Definition routines.c:2383
VPUBLIC int printEnergy(Vcom *com, NOsh *nosh, double totEnergy[NOSH_MAXCALC], int iprint)
Combine and pretty-print energy data (deprecated...see printElecEnergy)
Definition routines.c:2785
VPUBLIC int initMG(int icalc, NOsh *nosh, MGparm *mgparm, PBEparm *pbeparm, double realCenter[3], Vpbe *pbe[NOSH_MAXCALC], Valist *alist[NOSH_MAXMOL], Vgrid *dielXMap[NOSH_MAXMOL], Vgrid *dielYMap[NOSH_MAXMOL], Vgrid *dielZMap[NOSH_MAXMOL], Vgrid *kappaMap[NOSH_MAXMOL], Vgrid *chargeMap[NOSH_MAXMOL], Vpmgp *pmgp[NOSH_MAXCALC], Vpmg *pmg[NOSH_MAXCALC], Vgrid *potMap[NOSH_MAXMOL])
Initialize an MG calculation.
Definition routines.c:1212
VPUBLIC int printApolForce(Vcom *com, NOsh *nosh, int nforce[NOSH_MAXCALC], AtomForce *atomForce[NOSH_MAXCALC], int iprint)
Combine and pretty-print force data.
Definition routines.c:3471
VPUBLIC int writematMG(int rank, NOsh *nosh, PBEparm *pbeparm, Vpmg *pmg)
Write out operator matrix from MG calculation to file.
Definition routines.c:1799
VPUBLIC int writedataFE(int rank, NOsh *nosh, PBEparm *pbeparm, Vfetk *fetk)
Write FEM data to files.
Definition routines.c:4300
VPUBLIC int initAPOL(NOsh *nosh, Vmem *mem, Vparam *param, APOLparm *apolparm, int *nforce, AtomForce **atomForce, Valist *alist)
Upperlevel routine to the non-polar energy and force routines.
Definition routines.c:4475
VPUBLIC void printMGPARM(MGparm *mgparm, double realCenter[3])
Print out MG-specific params loaded from input.
Definition routines.c:1179
VPUBLIC void killFE(NOsh *nosh, Vpbe *pbe[NOSH_MAXCALC], Vfetk *fetk[NOSH_MAXCALC], Gem *gm[NOSH_MAXMOL])
Kill structures initialized during an FE calculation.
Definition routines.c:3683
VPUBLIC int energyMG(NOsh *nosh, int icalc, Vpmg *pmg, int *nenergy, double *totEnergy, double *qfEnergy, double *qmEnergy, double *dielEnergy)
Calculate electrostatic energies from MG solution.
Definition routines.c:1569
VPUBLIC int forceMG(Vmem *mem, NOsh *nosh, PBEparm *pbeparm, MGparm *mgparm, Vpmg *pmg, int *nforce, AtomForce **atomForce, Valist *alist[NOSH_MAXMOL])
Calculate forces from MG solution.
Definition routines.c:1636
VPUBLIC int energyFE(NOsh *nosh, int icalc, Vfetk *fetk[NOSH_MAXCALC], int *nenergy, double *totEnergy, double *qfEnergy, double *qmEnergy, double *dielEnergy)
Calculate electrostatic energies from FE solution.
Definition routines.c:4179
#define APBSRC
Return code for APBS during failure.
Definition routines.h:93
VPUBLIC int printApolEnergy(NOsh *nosh, int iprint)
Combine and pretty-print energy data.
Definition routines.c:2918
VPUBLIC int loadKappaMaps(NOsh *nosh, Vgrid *map[NOSH_MAXMOL])
Load the kappa maps given in NOsh into grid objects.
Definition routines.c:664
VPUBLIC int loadDielMaps(NOsh *nosh, Vgrid *dielXMap[NOSH_MAXMOL], Vgrid *dielYMap[NOSH_MAXMOL], Vgrid *dielZMap[NOSH_MAXMOL])
Load the dielectric maps given in NOsh into grid objects.
Definition routines.c:250
VPUBLIC Vparam * loadParameter(NOsh *nosh)
Loads and returns parameter object.
Definition routines.c:60
VPUBLIC void killChargeMaps(NOsh *nosh, Vgrid *map[NOSH_MAXMOL])
Destroy the loaded charge maps.
Definition routines.c:984
VPUBLIC int loadMolecules(NOsh *nosh, Vparam *param, Valist *alist[NOSH_MAXMOL])
Load the molecules given in NOsh into atom lists.
Definition routines.c:95
#define NOSH_MAXCALC
Maximum number of calculations in a run.
Definition nosh.h:87
VPUBLIC void NOsh_dtor(NOsh **thee)
Object destructor.
Definition nosh.c:354
#define NOSH_MAXMOL
Maximum number of molecules in a run.
Definition nosh.h:83
VPUBLIC int NOsh_setupElecCalc(NOsh *thee, Valist *alist[NOSH_MAXMOL])
Setup the series of electrostatics calculations.
Definition nosh.c:1374
VPUBLIC int NOsh_parseInput(NOsh *thee, Vio *sock)
Parse an input file from a socket.
Definition nosh.c:513
VPUBLIC NOsh * NOsh_ctor(int rank, int size)
Construct NOsh.
Definition nosh.c:308
VPUBLIC int NOsh_setupApolCalc(NOsh *thee, Valist *alist[NOSH_MAXMOL])
Setup the series of non-polar calculations.
Definition nosh.c:1469
@ NCT_BEM
Definition nosh.h:121
@ NCT_APOL
Definition nosh.h:120
@ NCT_PBSAM
Definition nosh.h:124
@ NCT_MG
Definition nosh.h:118
@ NCT_PBAM
Definition nosh.h:123
@ NCT_FEM
Definition nosh.h:119
@ NCT_GEOFLOW
Definition nosh.h:122
@ NPT_ELECFORCE
Definition nosh.h:156
@ NPT_FORCE
Definition nosh.h:154
@ NPT_APOLENERGY
Definition nosh.h:157
@ NPT_ELECENERGY
Definition nosh.h:155
@ NPT_ENERGY
Definition nosh.h:153
@ NPT_APOLFORCE
Definition nosh.h:158
@ PCE_COMPS
Definition pbeparm.h:84
#define VEMBED(rctag)
Allows embedding of RCS ID tags in object files.
Definition vhal.h:556
#define APBS_TIMER_WALL_CLOCK
APBS total execution timer ID.
Definition vhal.h:329
enum eVoutput_Format Voutput_Format
Declaration of the Voutput_Format type as the VOutput_Format enum.
Definition vhal.h:200
@ VRC_SUCCESS
Definition vhal.h:70
@ OUTPUT_NULL
Definition vhal.h:192
@ OUTPUT_FLAT
Definition vhal.h:193
VPUBLIC void Vparam_dtor(Vparam **thee)
Destroy object.
Definition vparam.c:213
VPUBLIC int Vstring_strcasecmp(const char *s1, const char *s2)
Case-insensitive string comparison (BSD standard)
Definition vstring.c:66
Header file for front end auxiliary routines.
Structure to hold atomic forces.
Definition routines.h:99
Reads and assigns charge/radii parameters.
Definition vparam.h:135
Parameter structure for APOL-specific variables from input files.
Definition apolparm.h:129
int molid
Definition apolparm.h:136
Parameter structure for BEM-specific variables from input files.
Definition bemparm.h:96
BEMparm_CalcType type
Definition bemparm.h:98
Parameter structure for FEM-specific variables from input files.
Definition femparm.h:133
int maxsolve
Definition femparm.h:165
Parameter structure for GEOFLOW-specific variables from input files.
Definition geoflowparm.h:98
Parameter structure for MG-specific variables from input files.
Definition mgparm.h:114
MGparm_CalcType type
Definition mgparm.h:116
FEMparm * femparm
Definition nosh.h:174
BEMparm * bemparm
Definition nosh.h:175
GEOFLOWparm * geoflowparm
Definition nosh.h:176
NOsh_CalcType calctype
Definition nosh.h:181
MGparm * mgparm
Definition nosh.h:173
APOLparm * apolparm
Definition nosh.h:180
PBSAMparm * pbsamparm
Definition nosh.h:178
PBAMparm * pbamparm
Definition nosh.h:177
PBEparm * pbeparm
Definition nosh.h:179
Class for parsing fixed format input files.
Definition nosh.h:195
char apolname[NOSH_MAXCALC][VMAX_ARGLEN]
Definition nosh.h:269
int nelec
Definition nosh.h:205
char elecname[NOSH_MAXCALC][VMAX_ARGLEN]
Definition nosh.h:267
int ncalc
Definition nosh.h:200
NOsh_calc * calc[NOSH_MAXCALC]
Definition nosh.h:197
int napol
Definition nosh.h:211
int apol2calc[NOSH_MAXCALC]
Definition nosh.h:229
NOsh_PrintType printwhat[NOSH_MAXPRINT]
Definition nosh.h:260
int elec2calc[NOSH_MAXCALC]
Definition nosh.h:221
int nprint
Definition nosh.h:259
Parameter structure for PBAM-specific variables from input files.
Definition pbamparm.h:105
Parameter structure for PBE variables from input files.
Definition pbeparm.h:117
PBEparm_calcEnergy calcenergy
Definition pbeparm.h:165
Parameter structure for PBSAM-specific variables from input files.
Definition pbsamparm.h:105
Container class for list of atom objects.
Definition valist.h:78
Contains public data members for Vfetk class/module.
Definition vfetk.h:176
Electrostatic potential oracle for Cartesian mesh data.
Definition vgrid.h:81
Contains public data members for Vpbe class/module.
Definition vpbe.h:84
Contains public data members for Vpmg class/module.
Definition vpmg.h:116
Contains public data members for Vpmgp class/module.
Definition vpmgp.h:80