APBS
1.4.1
|
FEtk master class (interface between FEtk and APBS) More...
Files | |
file | vfetk.c |
Class Vfetk methods. | |
file | vfetk.h |
Contains declarations for class Vfetk. | |
Data Structures | |
struct | sVfetk |
Contains public data members for Vfetk class/module. More... | |
struct | sVfetk_LocalVar |
Vfetk LocalVar subclass. More... | |
Macros | |
#define | VRINGMAX 1000 |
Maximum number of simplices in a simplex ring. | |
#define | VATOMMAX 1000000 |
Maximum number of atoms associated with a vertex. | |
Typedefs | |
typedef enum eVfetk_LsolvType | Vfetk_LsolvType |
Declare FEMparm_LsolvType type. | |
typedef enum eVfetk_MeshLoad | Vfetk_MeshLoad |
Declare FEMparm_GuessType type. | |
typedef enum eVfetk_NsolvType | Vfetk_NsolvType |
Declare FEMparm_NsolvType type. | |
typedef enum eVfetk_GuessType | Vfetk_GuessType |
Declare FEMparm_GuessType type. | |
typedef enum eVfetk_PrecType | Vfetk_PrecType |
Declare FEMparm_GuessType type. | |
typedef struct sVfetk_LocalVar | Vfetk_LocalVar |
Declaration of the Vfetk_LocalVar subclass as the Vfetk_LocalVar structure. | |
typedef struct sVfetk | Vfetk |
Declaration of the Vfetk class as the Vfetk structure. | |
Enumerations | |
enum | eVfetk_LsolvType { VLT_SLU =0, VLT_MG =1, VLT_CG =2, VLT_BCG =3 } |
Linear solver type. More... | |
enum | eVfetk_MeshLoad { VML_DIRICUBE, VML_NEUMCUBE, VML_EXTERNAL } |
Mesh loading operation. More... | |
enum | eVfetk_NsolvType { VNT_NEW =0, VNT_INC =1, VNT_ARC =2 } |
Non-linear solver type. More... | |
enum | eVfetk_GuessType { VGT_ZERO =0, VGT_DIRI =1, VGT_PREV =2 } |
Initial guess type. More... | |
enum | eVfetk_PrecType { VPT_IDEN =0, VPT_DIAG =1, VPT_MG =2 } |
Preconditioner type. More... | |
Functions | |
VPUBLIC double | Vfetk_energy (Vfetk *thee, int color, int nonlin) |
Return the total electrostatic energy. More... | |
VEXTERNC Gem * | Vfetk_getGem (Vfetk *thee) |
Get a pointer to the Gem (grid manager) object. More... | |
VEXTERNC AM * | Vfetk_getAM (Vfetk *thee) |
Get a pointer to the AM (algebra manager) object. More... | |
VEXTERNC Vpbe * | Vfetk_getVpbe (Vfetk *thee) |
Get a pointer to the Vpbe (PBE manager) object. More... | |
VEXTERNC Vcsm * | Vfetk_getVcsm (Vfetk *thee) |
Get a pointer to the Vcsm (charge-simplex map) object. More... | |
VEXTERNC int | Vfetk_getAtomColor (Vfetk *thee, int iatom) |
Get the partition information for a particular atom. More... | |
VEXTERNC Vfetk * | Vfetk_ctor (Vpbe *pbe, Vhal_PBEType type) |
Constructor for Vfetk object. More... | |
VEXTERNC int | Vfetk_ctor2 (Vfetk *thee, Vpbe *pbe, Vhal_PBEType type) |
FORTRAN stub constructor for Vfetk object. More... | |
VEXTERNC void | Vfetk_dtor (Vfetk **thee) |
Object destructor. More... | |
VEXTERNC void | Vfetk_dtor2 (Vfetk *thee) |
FORTRAN stub object destructor. More... | |
VEXTERNC double * | Vfetk_getSolution (Vfetk *thee, int *length) |
Create an array containing the solution (electrostatic potential in units of ![]() | |
VEXTERNC void | Vfetk_setParameters (Vfetk *thee, PBEparm *pbeparm, FEMparm *feparm) |
Set the parameter objects. More... | |
VEXTERNC double | Vfetk_dqmEnergy (Vfetk *thee, int color) |
Get the "mobile charge" and "polarization" contributions to the electrostatic energy. More... | |
VEXTERNC double | Vfetk_qfEnergy (Vfetk *thee, int color) |
Get the "fixed charge" contribution to the electrostatic energy. More... | |
VEXTERNC unsigned long int | Vfetk_memChk (Vfetk *thee) |
Return the memory used by this structure (and its contents) in bytes. More... | |
VEXTERNC void | Vfetk_setAtomColors (Vfetk *thee) |
Transfer color (partition ID) information frmo a partitioned mesh to the atoms. More... | |
VEXTERNC void | Bmat_printHB (Bmat *thee, char *fname) |
Writes a Bmat to disk in Harwell-Boeing sparse matrix format. More... | |
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) More... | |
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). More... | |
VEXTERNC PDE * | Vfetk_PDE_ctor (Vfetk *fetk) |
Constructs the FEtk PDE object. More... | |
VEXTERNC int | Vfetk_PDE_ctor2 (PDE *thee, Vfetk *fetk) |
Intializes the FEtk PDE object. More... | |
VEXTERNC void | Vfetk_PDE_dtor (PDE **thee) |
Destroys FEtk PDE object. More... | |
VEXTERNC void | Vfetk_PDE_dtor2 (PDE *thee) |
FORTRAN stub: destroys FEtk PDE object. More... | |
VEXTERNC void | Vfetk_PDE_initAssemble (PDE *thee, int ip[], double rp[]) |
Do once-per-assembly initialization. More... | |
VEXTERNC void | Vfetk_PDE_initElement (PDE *thee, int elementType, int chart, double tvx[][VAPBS_DIM], void *data) |
Do once-per-element initialization. More... | |
VEXTERNC void | Vfetk_PDE_initFace (PDE *thee, int faceType, int chart, double tnvec[]) |
Do once-per-face initialization. More... | |
VEXTERNC void | Vfetk_PDE_initPoint (PDE *thee, int pointType, int chart, double txq[], double tU[], double tdU[][VAPBS_DIM]) |
Do once-per-point initialization. More... | |
VEXTERNC void | Vfetk_PDE_Fu (PDE *thee, int key, double F[]) |
Evaluate strong form of PBE. For interior points, this is:
where
where | |
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:
where | |
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. This is the functional linearization of the strong form integrated with a test function to give:
where | |
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. More... | |
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. More... | |
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. More... | |
VEXTERNC void | Vfetk_PDE_bisectEdge (int dim, int dimII, int edgeType, int chart[], double vx[][VAPBS_DIM]) |
Define the way manifold edges are bisected. More... | |
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. More... | |
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 refinement at the dielectric boundary and (for non-regularized PBE) the charges. More... | |
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. More... | |
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 electrolyte where | |
VEXTERNC void | Vfetk_externalUpdateFunction (SS **simps, int num) |
External hook to simplex subdivision routines in Gem. Called each time a simplex is subdivided (we use it to update the charge-simplex map) More... | |
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, at all quadrature points on the master simplex element. More... | |
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, at all quadrature points on the master simplex element. More... | |
VEXTERNC void | Vfetk_readMesh (Vfetk *thee, int skey, Vio *sock) |
Read in mesh and initialize associated internal structures. More... | |
VEXTERNC void | Vfetk_dumpLocalVar () |
Debugging routine to print out local variables used by PDE object. More... | |
VEXTERNC int | Vfetk_fillArray (Vfetk *thee, Bvec *vec, Vdata_Type type) |
Fill an array with the specified data. More... | |
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. More... | |
VEXTERNC Vrc_Codes | Vfetk_loadGem (Vfetk *thee, Gem *gm) |
Load a Gem geometry manager object into Vfetk. More... | |
FEtk master class (interface between FEtk and APBS)
enum eVfetk_GuessType |
enum eVfetk_LsolvType |
enum eVfetk_MeshLoad |
enum eVfetk_NsolvType |
enum eVfetk_PrecType |
VEXTERNC void Bmat_printHB | ( | Bmat * | thee, |
char * | fname | ||
) |
VEXTERNC Vfetk* Vfetk_ctor | ( | Vpbe * | pbe, |
Vhal_PBEType | type | ||
) |
Constructor for Vfetk object.
pbe | Vpbe (PBE manager object) |
type | Version of PBE to solve |
VEXTERNC int Vfetk_ctor2 | ( | Vfetk * | thee, |
Vpbe * | pbe, | ||
Vhal_PBEType | type | ||
) |
FORTRAN stub constructor for Vfetk object.
thee | Vfetk object memory |
pbe | PBE manager object |
type | Version of PBE to solve |
VEXTERNC double Vfetk_dqmEnergy | ( | Vfetk * | thee, |
int | color | ||
) |
Get the "mobile charge" and "polarization" contributions to the electrostatic energy.
Using the solution at the finest mesh level, get the electrostatic energy due to the interaction of the mobile charges with the potential and polarization of the dielectric medium:
for the NPBE and
for the LPBE. Here denotes the counterion species,
is the bulk ionic strength,
is the modified Debye-Huckel parameter,
is the concentration of species
,
is the charge of species
,
is the dielectric function, and
is the dimensionless electrostatic potential. The energy is scaled to units of
.
thee | Vfetk object |
color | Partition restriction for energy evaluation, only used if non-negative |
thee | The Vfetk object |
color | Partition restriction for energy calculation; if non-negative, energy calculation is restricted to the specified partition (indexed by simplex and atom colors |
VEXTERNC void Vfetk_dtor | ( | Vfetk ** | thee | ) |
VEXTERNC void Vfetk_dtor2 | ( | Vfetk * | thee | ) |
VEXTERNC void Vfetk_dumpLocalVar | ( | ) |
VEXTERNC double Vfetk_energy | ( | Vfetk * | thee, |
int | color, | ||
int | nonlin | ||
) |
Return the total electrostatic energy.
Using the solution at the finest mesh level, get the electrostatic energy using the free energy functional for the Poisson-Boltzmann equation without removing any self-interaction terms (i.e., removing the reference state of isolated charges present in an infinite dielectric continuum with the same relative permittivity as the interior of the protein) and return the result in units of . The argument color allows the user to control the partition on which this energy is calculated; if (color == -1) no restrictions are used. The solution is obtained from the finest level of the passed AM object, but atomic data from the Vfetk object is used to calculate the energy.
< Total energy
<
<
thee | THe Vfetk object |
color | Partition restriction for energy calculation; if non-negative, energy calculation is restricted to the specified partition (indexed by simplex and atom colors |
nonlin | If 1, the NPBE energy functional is used; otherwise, the LPBE energy functional is used. If -2, SMPBE is used. |
VEXTERNC void Vfetk_externalUpdateFunction | ( | SS ** | simps, |
int | num | ||
) |
External hook to simplex subdivision routines in Gem. Called each time a simplex is subdivided (we use it to update the charge-simplex map)
simps | List of parent (simps[0]) and children (remainder) simplices |
num | Number of simplices in list |
VEXTERNC int Vfetk_fillArray | ( | Vfetk * | thee, |
Bvec * | vec, | ||
Vdata_Type | type | ||
) |
Fill an array with the specified data.
thee | The Vfetk object with the data |
vec | The vector to hold the data |
type | THe type of data to write |
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)
Generates a new cube mesh within the provided Vfetk object based on the specified mesh type. Creates a new copy of the mesh based on the global variables at the top of the file and the mesh type, then recenters the mesh based on the center and length variables provided to the function.
thee | Vfetk object |
center | Center for mesh, which the new mesh will adjust to |
length | Mesh lengths, which the new mesh will adjust to |
meshType | Mesh boundary conditions |
VEXTERNC AM* Vfetk_getAM | ( | Vfetk * | thee | ) |
VEXTERNC int Vfetk_getAtomColor | ( | Vfetk * | thee, |
int | iatom | ||
) |
VEXTERNC Gem* Vfetk_getGem | ( | Vfetk * | thee | ) |
VEXTERNC double* Vfetk_getSolution | ( | Vfetk * | thee, |
int * | length | ||
) |
Create an array containing the solution (electrostatic potential in units of ) at the finest mesh level.
thee | Vfetk object with solution |
length | Ste to length of the newly created solution array |
VEXTERNC Vrc_Codes Vfetk_loadGem | ( | Vfetk * | thee, |
Gem * | gm | ||
) |
Load a Gem geometry manager object into Vfetk.
thee | Destination |
gm | Geometry manager source |
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).
If we have an external mesh, load that external mesh from the provided socket. If we specify a non-external mesh type, we generate a new mesh cube based on templates. We then create and store a new Vcsm object in our Vfetk structure, which will carry the mesh data.
thee | Vfetk object to load into |
center | Center for mesh (if constructed) |
length | Mesh lengths (if constructed) |
meshType | Type of mesh to load |
sock | Socket for external mesh data (NULL otherwise) |
VEXTERNC unsigned long int Vfetk_memChk | ( | Vfetk * | thee | ) |
VEXTERNC void Vfetk_PDE_bisectEdge | ( | int | dim, |
int | dimII, | ||
int | edgeType, | ||
int | chart[], | ||
double | vx[][VAPBS_DIM] | ||
) |
Define the way manifold edges are bisected.
dim | Intrinsic dimension of manifold |
dimII | Embedding dimension of manifold |
edgeType | Type of edge being refined |
chart | Chart for edge vertices, used here as accessibility bitfields |
vx | Edge vertex coordindates |
VEXTERNC PDE* Vfetk_PDE_ctor | ( | Vfetk * | fetk | ) |
VEXTERNC int Vfetk_PDE_ctor2 | ( | PDE * | thee, |
Vfetk * | fetk | ||
) |
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.
thee | PDE object |
type | Vertex type |
chart | Chart for point coordinates |
txq | Point coordinates |
user | Vertex object pointer |
F | Set to delta function value |
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. This is the functional linearization of the strong form integrated with a test function to give:
where denotes the functional derivation of the mobile ion term.
thee | The PDE object |
key | Integrand to evaluate (0 = interior weak form, 1 = boundary weak form) |
W | Trial function value at current point |
dW | Trial function gradient at current point |
V | Test function value at current point |
dV | Test function gradient |
VEXTERNC void Vfetk_PDE_dtor | ( | PDE ** | thee | ) |
VEXTERNC void Vfetk_PDE_dtor2 | ( | PDE * | thee | ) |
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 ion term and
is the source charge distribution term (for PBE) or the induced surface charge distribution (for RPBE). For an interior-boundary (simplex face) point, this is:
where is the normal to the simplex face and the term represents the jump in dielectric displacement across the face. There is no outer-boundary contribution for this problem.
This function is not thread-safe
This function is not implemented (sets error to zero)
thee | The PDE object |
key | Type of point (0 = interior, 1 = boundary, 2 = interior boundary |
F | Set to value of residual |
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:
where denotes the mobile ion term.
thee | The PDE object |
key | Integrand to evaluate (0 = interior weak form, 1 = boundary weak form |
V | Test function at current point |
dV | Test function derivative at current point |
VEXTERNC void Vfetk_PDE_initAssemble | ( | PDE * | thee, |
int | ip[], | ||
double | rp[] | ||
) |
VEXTERNC void Vfetk_PDE_initElement | ( | PDE * | thee, |
int | elementType, | ||
int | chart, | ||
double | tvx[][VAPBS_DIM], | ||
void * | data | ||
) |
Do once-per-element initialization.
thee | PDE object |
elementType | Material type (not used) |
chart | Chart in which the vertex coordinates are provided, used here as a bitfield to store molecular accessibility |
tvx | Vertex coordinates |
data | Simplex pointer (hack) |
VEXTERNC void Vfetk_PDE_initFace | ( | PDE * | thee, |
int | faceType, | ||
int | chart, | ||
double | tnvec[] | ||
) |
Do once-per-face initialization.
thee | THe PDE object |
faceType | Simplex face type (interior or various boundary types) |
chart | Chart in which the vertex coordinates are provided, used here as a bitfield for molecular accessibility |
tnvec | Coordinates of outward normal vector for face |
VEXTERNC void Vfetk_PDE_initPoint | ( | PDE * | thee, |
int | pointType, | ||
int | chart, | ||
double | txq[], | ||
double | tU[], | ||
double | tdU[][VAPBS_DIM] | ||
) |
Do once-per-point initialization.
This function is not thread-safe
This function uses pre-defined boudnary definitions for the molecular surface.
thee | The PDE object |
pointType | The type of point – interior or various faces |
chart | The chart in which the point coordinates are provided, used here as bitfield for molecular accessibility |
txq | Point coordinates |
tU | Solution value at point |
tdU | Solution derivative at point |
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 electrolyte where is the output from Vpbe_getZmagic.
thee | The PDE object |
key | What to evluate: interior (0) or boundary (1)? |
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.
dim | Intrinsic dimension of manifold |
dimII | Embedding dimension of manifold |
vertexType | Type of vertex |
chart | Chart for vertex coordinates |
vx | Vertex coordinates |
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 refinement at the dielectric boundary and (for non-regularized PBE) the charges.
dim | Intrinsic manifold dimension |
dimII | Embedding manifold dimension |
simplexType | Type of simplex being refined |
faceType | Types of faces in simplex |
vertexType | Types of vertices in simplex |
chart | Charts for vertex coordinates |
vx | Vertex coordinates |
simplex | Simplex pointer |
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.
dim | Intrinsic manifold dimension |
dimII | Embedding manifold dimension |
objType | ??? |
chart | Charts of vertices' coordinates |
vx | Vertices' coordinates |
dimV | Number of vertices |
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, at all quadrature points on the master simplex element.
key | Basis type to evaluate (0 = trial, 1 = test, 2 = trialB, 3 = testB) |
dim | Spatial dimension |
comp | Which component of elliptic system to produce basis for |
pdkey | Basis partial differential equation evaluation key:
|
xq | Set to quad pt coordinate |
basis | Set to all basis functions evaluated at all quadrature pts |
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, at all quadrature points on the master simplex element.
* The basis ordering is important. For a fixed quadrature * point iq, you must follow the following ordering in p[iq][], * based on how you specify the degrees of freedom in dof[]: * * <v_0 vDF_0>, <v_1 vDF_0>, ..., <v_{nv} vDF_0> * <v_0 vDF_1>, <v_1 vDF_1>, ..., <v_{nv} vDF_1> * ... * <v_0 vDF_{nvDF}>, <v_0 vDF_{nvDF}>, ..., <v_{nv} vDF_{nvDF}> * * <e_0 eDF_0>, <e_1 eDF_0>, ..., <e_{ne} eDF_0> * <e_0 eDF_1>, <e_1 eDF_1>, ..., <e_{ne} eDF_1> * ... * <e_0 eDF_{neDF}>, <e_1 eDF_{neDF}>, ..., <e_{ne} eDF_{neDF}> * * <f_0 fDF_0>, <f_1 fDF_0>, ..., <f_{nf} fDF_0> * <f_0 fDF_1>, <f_1 fDF_1>, ..., <f_{nf} fDF_1> * ... * <f_0 fDF_{nfDF}>, <f_1 fDF_{nfDF}>, ..., <f_{nf} fDF_{nfDF}> * * <s_0 sDF_0>, <s_1 sDF_0>, ..., <s_{ns} sDF_0> * <s_0 sDF_1>, <s_1 sDF_1>, ..., <s_{ns} sDF_1> * ... * <s_0 sDF_{nsDF}>, <s_1 sDF_{nsDF}>, ..., <s_{ns} sDF_{nsDF}> * * For example, linear elements in R^3, with one degree of freedom at each * * vertex, would use the following ordering: * * <v_0 vDF_0>, <v_1 vDF_0>, <v_2 vDF_0>, <v_3 vDF_0> * * Quadratic elements in R^2, with one degree of freedom at each vertex and * edge, would use the following ordering: * * <v_0 vDF_0>, <v_1 vDF_0>, <v_2 vDF_0> * <e_0 eDF_0>, <e_1 eDF_0>, <e_2 eDF_0> * * You can use different trial and test spaces for each component of the * elliptic system, thereby allowing for the use of Petrov-Galerkin methods. * You MUST then tag the bilinear form symmetry entries as nonsymmetric in * your PDE constructor to reflect that DF(u)(w,v) will be different from * DF(u)(v,w), even if your form acts symmetrically when the same basis is * used for w and v. * * You can also use different trial spaces for each component of the elliptic * system, and different test spaces for each component of the elliptic * system. This allows you to e.g. use a basis which is vertex-based for * one component, and a basis which is edge-based for another. This is * useful in fluid mechanics, eletromagnetics, or simply to play around with * different elements. * * This function is called by MC to build new master elements whenever it * reads in a new mesh. Therefore, this function does not have to be all * that fast, and e.g. could involve symbolic computation. *
key | Basis type to evaluate (0 = trial, 1 = test, 2 = trialB, 3 = testB) |
dim | Spatial dimension |
comp | Which component of elliptic system to produce basis for? |
ndof | Set to the number of degrees of freedom |
dof | Set to degree of freedom per v/e/f/s |
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.
This function is hard-coded to call only multiple-sphere Debye-Hü functions.
This function is not thread-safe.
thee | PDE object |
type | Vertex boundary type |
chart | Chart for point coordinates |
txq | Point coordinates |
F | Set to boundary values |
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.
The signature here doesn't match what's in mc's src/pde/mc/pde.h, which g++ seems to dislike for GAMer integration. Trying a change of function signature to match to see if that makes g++ happy. Also see vfetk.h for similar signature change. - P. Ellis 11-8-2011
thee | PDE object |
type | Point type |
chart | Chart for point coordinates |
txq | Point coordinates |
F | Set to value at point |
VEXTERNC double Vfetk_qfEnergy | ( | Vfetk * | thee, |
int | color | ||
) |
Get the "fixed charge" contribution to the electrostatic energy.
Using the solution at the finest mesh level, get the electrostatic energy due to the interaction of the fixed charges with the potential:
and return the result in units of . Clearly, no self-interaction terms are removed. A factor a 1/2 has to be included to convert this to a real energy.
thee | Vfetk object |
color | Partition restriction for energy evaluation, only used if non-negative |
thee | The Vfetk object |
color | Partition restriction for energy evaluation, only used if non-negative |
VEXTERNC void Vfetk_readMesh | ( | Vfetk * | thee, |
int | skey, | ||
Vio * | sock | ||
) |
Read in mesh and initialize associated internal structures.
thee | THe Vfetk object |
skey | The sock format key (0 = MCSF simplex format) |
sock | Socket object ready for reading |
VEXTERNC void Vfetk_setAtomColors | ( | Vfetk * | thee | ) |
Transfer color (partition ID) information frmo a partitioned mesh to the atoms.
Transfer color information from partitioned mesh to the atoms. In the case that a charge is shared between two partitions, the partition color of the first simplex is selected. Due to the arbitrary nature of this selection, THIS METHOD SHOULD ONLY BE USED IMMEDIATELY AFTER PARTITIONING!!!
thee | THe Vfetk object |
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.
thee | Vfetk object |
vec | FEtk Bvec vector to use |
format | Format for data |
iodev | Output device type (FILE/BUFF/UNIX/INET) |
iofmt | Output device format (ASCII/XDR) |
thost | Output hostname (for sockets) |
fname | Output FILE/BUFF/UNIX/INET name |
thee | The Vfetk object |
iodev | Output device type (FILE = file, BUFF = buffer, UNIX = unix pipe, INET = network socket) |
iofmt | Output device format (ASCII = ascii/plaintext, XDR = xdr) |
thost | Output hostname for sockets |
fname | Output filename for other |
vec | Data vector |
format | Data format |