APBS  1.5
vfetk.h
Go to the documentation of this file.
1 
62 #ifndef _VFETK_H_
63 #define _VFETK_H_
64 
65 #include "apbscfg.h"
66 
67 #include "maloc/maloc.h"
68 #include "mc/mc.h"
69 
70 #include "generic/vhal.h"
71 #include "generic/vatom.h"
72 // #include "generic/valist.h"
73 #include "generic/vpbe.h"
74 #include "generic/vunit.h"
75 #include "generic/vgreen.h"
76 #include "generic/vcap.h"
77 #include "generic/pbeparm.h"
78 #include "generic/femparm.h"
79 #include "fem/vcsm.h"
80 
87  VLT_SLU=0,
88  VLT_MG=1,
89  VLT_CG=2,
91 };
92 
98 
99 
108 };
109 
115 
125 };
126 
132 
142 };
143 
149 
159 };
160 
166 
176 struct sVfetk {
177 
178  Vmem *vmem;
179  Gem *gm;
182  AM *am;
183  Aprx *aprx;
184  PDE *pde;
188  int lmax;
189  double ltol;
191  int nmax;
192  double ntol;
195  int pjac;
200  int level;
202 };
203 
207 typedef struct sVfetk Vfetk;
208 
216  double nvec[VAPBS_DIM];
217  double vx[4][VAPBS_DIM];
218  double xq[VAPBS_DIM];
219  double U[MAXV];
220  double dU[MAXV][VAPBS_DIM];
221  double W;
222  double dW[VAPBS_DIM];
223  double d2W;
224  int sType;
225  int fType;
226  double diel;
227  double ionacc;
228  double A;
229  double F;
230  double B;
231  double DB;
232  double jumpDiel;
235  int initGreen;
237  SS *simp;
239  VV *verts[4];
240  int nverts;
241  double ionConc[MAXION];
242  double ionQ[MAXION];
243  double ionRadii[MAXION];
244  double zkappa2;
245  double zks2;
246  double ionstr;
247  int nion;
248  double Fu_v;
249  double DFu_wv;
250  double delta;
251  double u_D;
252  double u_T;
253 };
254 
260 
261 #if !defined(VINLINE_VFETK)
262 
268  VEXTERNC Gem* Vfetk_getGem(
269  Vfetk *thee
270  );
271 
277  VEXTERNC AM* Vfetk_getAM(
278  Vfetk *thee
279  );
280 
286  VEXTERNC Vpbe* Vfetk_getVpbe(
287  Vfetk *thee
288  );
289 
295  VEXTERNC Vcsm* Vfetk_getVcsm(
296  Vfetk *thee
297  );
298 
305  VEXTERNC int Vfetk_getAtomColor(
306  Vfetk *thee,
307  int iatom
308  );
309 
310 #else /* if defined(VINLINE_VFETK) */
311 # define Vfetk_getGem(thee) ((thee)->gm)
312 # define Vfetk_getAM(thee) ((thee)->am)
313 # define Vfetk_getVpbe(thee) ((thee)->pbe)
314 # define Vfetk_getVcsm(thee) ((thee)->csm)
315 # define Vfetk_getAtomColor(thee, iatom) (Vatom_getPartID(Valist_getAtom(Vpbe_getValist(thee->pbe), iatom)))
316 #endif /* if !defined(VINLINE_VFETK) */
317 
318 /* ///////////////////////////////////////////////////////////////////////////
319 // Class Vfetk: Non-Inlineable methods (vfetk.c)
321 
331 VEXTERNC Vfetk* Vfetk_ctor(
332  Vpbe *pbe, /**< Vpbe (PBE manager object) */
333  Vhal_PBEType type
334  );
335 
345 VEXTERNC int Vfetk_ctor2(
346  Vfetk *thee,
347  Vpbe *pbe,
348  Vhal_PBEType type
349  );
350 
356 VEXTERNC void Vfetk_dtor(
357  Vfetk **thee
358  );
359 
365 VEXTERNC void Vfetk_dtor2(
366  Vfetk *thee
367  );
368 
378 VEXTERNC double* Vfetk_getSolution(
379  Vfetk *thee,
380  int *length
381  );
382 
388 VEXTERNC void Vfetk_setParameters(
389  Vfetk *thee,
390  PBEparm *pbeparm,
391  FEMparm *feparm
392  );
393 
412 VEXTERNC double Vfetk_energy(
413  Vfetk *thee,
414  int color,
418  int nonlin
420  );
421 
451 VEXTERNC double Vfetk_dqmEnergy(
452  Vfetk *thee,
453  int color
457  );
458 
476 VEXTERNC double Vfetk_qfEnergy(
477  Vfetk *thee,
478  int color
480  );
481 
489 VEXTERNC unsigned long int Vfetk_memChk(
490  Vfetk *thee
491  );
492 
508 VEXTERNC void Vfetk_setAtomColors(
509  Vfetk *thee
510  );
511 
520 VEXTERNC void Bmat_printHB(
521  Bmat *thee,
522  char *fname
523  );
524 
530 VEXTERNC Vrc_Codes Vfetk_genCube(
531  Vfetk *thee,
532  double center[3],
533  double length[3],
534  Vfetk_MeshLoad meshType
535  );
536 
542 VEXTERNC Vrc_Codes Vfetk_loadMesh(
543  Vfetk *thee,
544  double center[3],
545  double length[3],
546  Vfetk_MeshLoad meshType,
547  Vio *sock
548  );
549 
556 VEXTERNC PDE* Vfetk_PDE_ctor(
557  Vfetk *fetk
558  );
559 
566 VEXTERNC int Vfetk_PDE_ctor2(
567  PDE *thee,
568  Vfetk *fetk
569  );
570 
577 VEXTERNC void Vfetk_PDE_dtor(
578  PDE **thee
579  );
580 
587 VEXTERNC void Vfetk_PDE_dtor2(
588  PDE *thee
589  );
590 
596 VEXTERNC void Vfetk_PDE_initAssemble(
597  PDE *thee,
598  int ip[],
599  double rp[]
600  );
601 
608 VEXTERNC void Vfetk_PDE_initElement(
609  PDE *thee,
610  int elementType,
611  int chart,
614  double tvx[][VAPBS_DIM],
615  void *data
616  );
617 
623 VEXTERNC void Vfetk_PDE_initFace(
624  PDE *thee,
625  int faceType,
627  int chart,
629  double tnvec[]
630  );
631 
639 VEXTERNC void Vfetk_PDE_initPoint(
640  PDE *thee,
641  int pointType,
642  int chart,
644  double txq[],
645  double tU[],
646  double tdU[][VAPBS_DIM]
647  );
648 
666 VEXTERNC void Vfetk_PDE_Fu(
667  PDE *thee,
668  int key,
670  double F[]
671  );
672 
683 VEXTERNC double Vfetk_PDE_Fu_v(
684  PDE *thee,
685  int key,
687  double V[],
688  double dV[][VAPBS_DIM]
689  );
690 
702 VEXTERNC double Vfetk_PDE_DFu_wv(
703  PDE *thee,
704  int key,
706  double W[],
707  double dW[][VAPBS_DIM],
708  double V[],
709  double dV[][VAPBS_DIM]
710  );
711 
718 VEXTERNC void Vfetk_PDE_delta(
719  PDE *thee,
720  int type,
721  int chart,
722  double txq[],
723  void *user,
724  double F[]
725  );
726 
734 VEXTERNC void Vfetk_PDE_u_D(
735  PDE *thee,
736  int type,
737  int chart,
738  double txq[],
739  double F[]
740  );
741 
749 VEXTERNC void Vfetk_PDE_u_T(
750  PDE *thee,
751  int type,
752  int chart,
753  double txq[],
754  double F[]
755  );
756 
762 VEXTERNC void Vfetk_PDE_bisectEdge(
763  int dim,
764  int dimII,
765  int edgeType,
766  int chart[],
768  double vx[][VAPBS_DIM]
769  );
770 
776 VEXTERNC void Vfetk_PDE_mapBoundary(
777  int dim,
778  int dimII,
779  int vertexType,
780  int chart,
781  double vx[VAPBS_DIM]
782  );
783 
792 VEXTERNC int Vfetk_PDE_markSimplex(
793  int dim,
794  int dimII,
795  int simplexType,
796  int faceType[VAPBS_NVS],
797  int vertexType[VAPBS_NVS],
798  int chart[],
799  double vx[][VAPBS_DIM],
800  void *simplex
801  );
802 
808 VEXTERNC void Vfetk_PDE_oneChart(
809  int dim,
810  int dimII,
811  int objType,
812  int chart[],
813  double vx[][VAPBS_DIM],
814  int dimV
815  );
816 
826 VEXTERNC double Vfetk_PDE_Ju(
827  PDE *thee,
828  int key
829  );
830 
838 VEXTERNC void Vfetk_externalUpdateFunction(
839  SS **simps,
841  int num
842  );
843 
844 
907 VEXTERNC int Vfetk_PDE_simplexBasisInit(
908  int key,
910  int dim,
911  int comp,
913  int *ndof,
914  int dof[]
915  );
916 
924 VEXTERNC void Vfetk_PDE_simplexBasisForm(
925  int key,
927  int dim,
928  int comp ,
929  int pdkey,
938  double xq[],
939  double basis[]
941  );
942 
948 VEXTERNC void Vfetk_readMesh(
949  Vfetk *thee,
950  int skey,
951  Vio *sock
952  );
953 
959 VEXTERNC void Vfetk_dumpLocalVar();
960 
968 VEXTERNC int Vfetk_fillArray(
969  Vfetk *thee,
970  Bvec *vec,
971  Vdata_Type type
972  );
973 
988 VEXTERNC int Vfetk_write(
989  Vfetk *thee,
990  const char *iodev,
992  const char *iofmt,
994  const char *thost,
995  const char *fname,
996  Bvec *vec,
997  Vdata_Format format
998  );
999 
1005 VEXTERNC Vrc_Codes Vfetk_loadGem(
1006  Vfetk *thee,
1007  Gem *gm
1008  );
1009 
1010 
1011 #endif /* ifndef _VFETK_H_ */
double ionstr
Definition: vfetk.h:246
Definition: vfetk.h:89
double xq[VAPBS_DIM]
Definition: vfetk.h:218
VEXTERNC void Vfetk_PDE_initPoint(PDE *thee, int pointType, int chart, double txq[], double tU[], double tdU[][VAPBS_DIM])
Do once-per-point initialization.
double jumpDiel
Definition: vfetk.h:232
Definition: vfetk.h:122
double DFu_wv
Definition: vfetk.h:249
Contains public data members for Vgreen class/module.
Definition: vgreen.h:82
VEXTERNC void Vfetk_PDE_delta(PDE *thee, int type, int chart, double txq[], void *user, double F[])
Evaluate a (discretized) delta function source term at the given point.
Definition: vfetk.c:1780
VEXTERNC int Vfetk_PDE_markSimplex(int dim, int dimII, int simplexType, int faceType[VAPBS_NVS], int vertexType[VAPBS_NVS], int chart[], double vx[][VAPBS_DIM], void *simplex)
User-defined error estimator – in our case, a geometry-based refinement method; forcing simplex ref...
Contains declarations for the Vcsm class.
Contains declarations for class Vgreen.
FEMparm * feparm
Definition: vfetk.h:198
enum eVdata_Format Vdata_Format
Declaration of the Vdata_Format type as the Vdata_Format enum.
Definition: vhal.h:323
PBEparm * pbeparm
Definition: vfetk.h:197
int lmax
Definition: vfetk.h:188
double u_T
Definition: vfetk.h:252
VEXTERNC void Vfetk_setAtomColors(Vfetk *thee)
Transfer color (partition ID) information frmo a partitioned mesh to the atoms.
Definition: vfetk.c:849
Vfetk LocalVar subclass.
Definition: vfetk.h:215
VEXTERNC void Bmat_printHB(Bmat *thee, char *fname)
Writes a Bmat to disk in Harwell-Boeing sparse matrix format.
Definition: vfetk.c:1026
Gem * gm
Definition: vfetk.h:179
VEXTERNC Vcsm * Vfetk_getVcsm(Vfetk *thee)
Get a pointer to the Vcsm (charge-simplex map) object.
Definition: vfetk.c:510
Vfetk_LsolvType lkey
Definition: vfetk.h:187
double vx[4][VAPBS_DIM]
Definition: vfetk.h:217
VEXTERNC void Vfetk_PDE_simplexBasisForm(int key, int dim, int comp, int pdkey, double xq[], double basis[])
Evaluate the bases for the trial or test space, for a particular component of the system,...
Definition: vfetk.c:2203
Contains public data members for Vpbe class/module.
Definition: vpbe.h:84
eVfetk_PrecType
Preconditioner type.
Definition: vfetk.h:155
double ionConc[MAXION]
Definition: vfetk.h:241
double A
Definition: vfetk.h:228
VEXTERNC void Vfetk_externalUpdateFunction(SS **simps, int num)
External hook to simplex subdivision routines in Gem. Called each time a simplex is subdivided (we us...
Definition: vfetk.c:2078
VEXTERNC void Vfetk_dtor2(Vfetk *thee)
FORTRAN stub object destructor.
Definition: vfetk.c:633
enum eVfetk_NsolvType Vfetk_NsolvType
Declare FEMparm_NsolvType type.
Definition: vfetk.h:131
VEXTERNC int Vfetk_PDE_simplexBasisInit(int key, int dim, int comp, int *ndof, int dof[])
Initialize the bases for the trial or the test space, for a particular component of the system,...
Definition: vfetk.c:2140
VEXTERNC void Vfetk_PDE_mapBoundary(int dim, int dimII, int vertexType, int chart, double vx[VAPBS_DIM])
Map a boundary point to some pre-defined shape.
VEXTERNC void Vfetk_PDE_dtor(PDE **thee)
Destroys FEtk PDE object.
Definition: vfetk.c:1231
VEXTERNC Vpbe * Vfetk_getVpbe(Vfetk *thee)
Get a pointer to the Vpbe (PBE manager) object.
Definition: vfetk.c:503
Contains declarations for class APOLparm.
double U[MAXV]
Definition: vfetk.h:219
VEXTERNC void Vfetk_PDE_initAssemble(PDE *thee, int ip[], double rp[])
Do once-per-assembly initialization.
Definition: vfetk.c:1508
int nmax
Definition: vfetk.h:191
VEXTERNC PDE * Vfetk_PDE_ctor(Vfetk *fetk)
Constructs the FEtk PDE object.
Definition: vfetk.c:1176
double u_D
Definition: vfetk.h:251
double d2W
Definition: vfetk.h:223
VEXTERNC double Vfetk_PDE_DFu_wv(PDE *thee, int key, double W[], double dW[][VAPBS_DIM], double V[], double dV[][VAPBS_DIM])
This is the linearization of the weak form of the PBE; e.g., for use in a Newton iteration....
Contains declarations for class Vpbe.
Contains declarations for class Vcap.
eVfetk_GuessType
Initial guess type.
Definition: vfetk.h:138
Vmem * vmem
Definition: vfetk.h:178
VEXTERNC void Vfetk_dumpLocalVar()
Debugging routine to print out local variables used by PDE object.
Definition: vfetk.c:2255
VEXTERNC AM * Vfetk_getAM(Vfetk *thee)
Get a pointer to the AM (algebra manager) object.
Definition: vfetk.c:497
VEXTERNC int Vfetk_getAtomColor(Vfetk *thee, int iatom)
Get the partition information for a particular atom.
Definition: vfetk.c:517
double nvec[VAPBS_DIM]
Definition: vfetk.h:216
VEXTERNC void Vfetk_PDE_u_D(PDE *thee, int type, int chart, double txq[], double F[])
Evaluate the Dirichlet boundary condition at the given point.
Definition: vfetk.c:1867
VEXTERNC Vrc_Codes Vfetk_loadMesh(Vfetk *thee, double center[3], double length[3], Vfetk_MeshLoad meshType, Vio *sock)
Loads a mesh into the Vfetk (and associated) object(s).
Definition: vfetk.c:980
Vfetk_PrecType lprec
Definition: vfetk.h:194
Definition: vfetk.h:123
VEXTERNC int Vfetk_ctor2(Vfetk *thee, Vpbe *pbe, Vhal_PBEType type)
FORTRAN stub constructor for Vfetk object.
Definition: vfetk.c:545
VEXTERNC Vrc_Codes Vfetk_loadGem(Vfetk *thee, Gem *gm)
Load a Gem geometry manager object into Vfetk.
Contains a collection of useful constants and conversion factors.
VEXTERNC void Vfetk_dtor(Vfetk **thee)
Object destructor.
Definition: vfetk.c:625
VEXTERNC void Vfetk_setParameters(Vfetk *thee, PBEparm *pbeparm, FEMparm *feparm)
Set the parameter objects.
Definition: vfetk.c:615
VEXTERNC double Vfetk_PDE_Fu_v(PDE *thee, int key, double V[], double dV[][VAPBS_DIM])
This is the weak form of the PBE; i.e. the strong form integrated with a test function to give: wher...
Definition: vfetk.c:1703
eVfetk_MeshLoad
Mesh loading operation.
Definition: vfetk.h:104
Vfetk_NsolvType nkey
Definition: vfetk.h:190
int level
Definition: vfetk.h:200
Contains declarations for class Vatom.
Vfetk_GuessType gues
Definition: vfetk.h:193
eVfetk_LsolvType
Linear solver type.
Definition: vfetk.h:86
enum eVhal_PBEType Vhal_PBEType
Declaration of the Vhal_PBEType type as the Vhal_PBEType enum.
Definition: vhal.h:151
enum eVfetk_LsolvType Vfetk_LsolvType
Declare FEMparm_LsolvType type.
Definition: vfetk.h:97
eVfetk_NsolvType
Non-linear solver type.
Definition: vfetk.h:121
double ionRadii[MAXION]
Definition: vfetk.h:243
Charge-simplex map class.
Definition: vcsm.h:89
VEXTERNC void Vfetk_PDE_initElement(PDE *thee, int elementType, int chart, double tvx[][VAPBS_DIM], void *data)
Do once-per-element initialization.
Vgreen * green
Definition: vfetk.h:234
VEXTERNC double * Vfetk_getSolution(Vfetk *thee, int *length)
Create an array containing the solution (electrostatic potential in units of ) at the finest mesh lev...
Definition: vfetk.c:641
VEXTERNC double Vfetk_PDE_Ju(PDE *thee, int key)
Energy functional. This returns the energy (less delta function terms) in the form: for a 1:1 electr...
Definition: vfetk.c:2003
double F
Definition: vfetk.h:229
VV * verts[4]
Definition: vfetk.h:239
double ionQ[MAXION]
Definition: vfetk.h:242
Parameter structure for FEM-specific variables from input files.
Definition: femparm.h:133
Vhal_PBEType type
Definition: vfetk.h:199
double B
Definition: vfetk.h:230
VEXTERNC double Vfetk_qfEnergy(Vfetk *thee, int color)
Get the "fixed charge" contribution to the electrostatic energy.
Definition: vfetk.c:732
PDE * pde
Definition: vfetk.h:184
double zks2
Definition: vfetk.h:245
Definition: vfetk.h:124
VEXTERNC double Vfetk_dqmEnergy(Vfetk *thee, int color)
Get the "mobile charge" and "polarization" contributions to the electrostatic energy.
Definition: vfetk.c:842
Vcsm * csm
Definition: vfetk.h:186
Aprx * aprx
Definition: vfetk.h:183
double delta
Definition: vfetk.h:250
Contains generic macro definitions for APBS.
VEXTERNC void Vfetk_readMesh(Vfetk *thee, int skey, Vio *sock)
Read in mesh and initialize associated internal structures.
#define VAPBS_DIM
Our dimension.
Definition: vhal.h:402
VEXTERNC void Vfetk_PDE_Fu(PDE *thee, int key, double F[])
Evaluate strong form of PBE. For interior points, this is: where is the (possibly nonlinear) mobile...
Definition: vfetk.c:1695
Contains declarations for class PBEparm.
VEXTERNC void Vfetk_PDE_bisectEdge(int dim, int dimII, int edgeType, int chart[], double vx[][VAPBS_DIM])
Define the way manifold edges are bisected.
double zkappa2
Definition: vfetk.h:244
Definition: vfetk.h:90
Parameter structure for PBE variables from input files.
Definition: pbeparm.h:117
Definition: vfetk.h:158
double dW[VAPBS_DIM]
Definition: vfetk.h:222
VEXTERNC void Vfetk_PDE_dtor2(PDE *thee)
FORTRAN stub: destroys FEtk PDE object.
Definition: vfetk.c:1246
Definition: vfetk.h:87
AM * am
Definition: vfetk.h:182
VEXTERNC Gem * Vfetk_getGem(Vfetk *thee)
Get a pointer to the Gem (grid manager) object.
Definition: vfetk.c:490
double DB
Definition: vfetk.h:231
VEXTERNC void Vfetk_PDE_initFace(PDE *thee, int faceType, int chart, double tnvec[])
Do once-per-face initialization.
Definition: vfetk.c:1559
enum eVdata_Type Vdata_Type
Declaration of the Vdata_Type type as the Vdata_Type enum.
Definition: vhal.h:302
double Fu_v
Definition: vfetk.h:248
VEXTERNC void Vfetk_PDE_u_T(PDE *thee, int type, int chart, double txq[], double F[])
Evaluate the "true solution" at the given point for comparison with the numerical solution.
Definition: vfetk.c:1886
enum eVfetk_PrecType Vfetk_PrecType
Declare FEMparm_GuessType type.
Definition: vfetk.h:165
VEXTERNC void Vfetk_PDE_oneChart(int dim, int dimII, int objType, int chart[], double vx[][VAPBS_DIM], int dimV)
Unify the chart for different coordinate systems – a no-op for us.
VEXTERNC Vrc_Codes Vfetk_genCube(Vfetk *thee, double center[3], double length[3], Vfetk_MeshLoad meshType)
Construct a rectangular mesh (in the current Vfetk object)
Definition: vfetk.c:885
Contains public data members for Vfetk class/module.
Definition: vfetk.h:176
Definition: vfetk.h:88
double ltol
Definition: vfetk.h:189
VEXTERNC int Vfetk_PDE_ctor2(PDE *thee, Vfetk *fetk)
Intializes the FEtk PDE object.
Definition: vfetk.c:1187
Vpbe * pbe
Definition: vfetk.h:185
int initGreen
Definition: vfetk.h:235
#define MAXION
The maximum number of ion species that can be involved in a single PBE calculation.
Definition: vhal.h:377
VPUBLIC double Vfetk_energy(Vfetk *thee, int color, int nonlin)
Return the total electrostatic energy.
Definition: vfetk.c:693
enum eVfetk_MeshLoad Vfetk_MeshLoad
Declare FEMparm_GuessType type.
Definition: vfetk.h:114
#define VAPBS_NVS
Number of vertices per simplex (hard-coded to 3D)
Definition: vhal.h:397
double dU[MAXV][VAPBS_DIM]
Definition: vfetk.h:220
VEXTERNC int Vfetk_fillArray(Vfetk *thee, Bvec *vec, Vdata_Type type)
Fill an array with the specified data.
Definition: vfetk.c:2299
int pjac
Definition: vfetk.h:195
VEXTERNC int Vfetk_write(Vfetk *thee, const char *iodev, const char *iofmt, const char *thost, const char *fname, Bvec *vec, Vdata_Format format)
Write out data.
Definition: vfetk.c:2464
VEXTERNC unsigned long int Vfetk_memChk(Vfetk *thee)
Return the memory used by this structure (and its contents) in bytes.
Definition: vfetk.c:867
Vfetk * fetk
Definition: vfetk.h:233
double ntol
Definition: vfetk.h:192
enum eVfetk_GuessType Vfetk_GuessType
Declare FEMparm_GuessType type.
Definition: vfetk.h:148
double W
Definition: vfetk.h:221