Main Page   Modules   Data Structures   File List   Data Fields   Globals   Related Pages  

Vfetk class

FEtk master class (interface between FEtk and APBS). More...


Files

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 Contains variables used when solving the PDE with FEtk. More...


Typedefs

typedef enum eVfetk_LsolvType Vfetk_LsolvType
 Declare FEMparm_LsolvType 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 sVfetk_LocalVar Vfetk_LocalVar
 Declaration of the Vfetk_LocalVar subclass as the Vfetk_LocalVar structure.

typedef 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_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

Gem * Vfetk_getGem (Vfetk *thee)
 Get a pointer to the Gem (grid manager) object.

AM * Vfetk_getAM (Vfetk *thee)
 Get a pointer to the AM (algebra manager) object.

VpbeVfetk_getVpbe (Vfetk *thee)
 Get a pointer to the Vpbe (PBE manager) object.

VcsmVfetk_getVcsm (Vfetk *thee)
 Get a pointer to the Vcsm (charge-simplex map) object.

int Vfetk_getAtomColor (Vfetk *thee, int iatom)
 Get the partition information for a particular atom.

Vfetk * Vfetk_ctor (Vpbe *pbe, Vhal_PBEType type)
 Constructor for Vfetk object.

int Vfetk_ctor2 (Vfetk *thee, Vpbe *pbe, Vhal_PBEType type)
 FORTRAN stub constructor for Vfetk object.

void Vfetk_dtor (Vfetk **thee)
 Object destructor.

void Vfetk_dtor2 (Vfetk *thee)
 FORTRAN stub object destructor.

double * Vfetk_getSolution (Vfetk *thee, int *length)
 Create an array containing the solution (electrostatic potential in units of ) at the finest mesh level.

void Vfetk_setParameters (Vfetk *thee, PBEparm *pbeparm, FEMparm *feparm)
 Set the parameter objects.

double Vfetk_energy (Vfetk *thee, int color, int nonlin)
 Return the total electrostatic energy.

double Vfetk_dqmEnergy (Vfetk *thee, int color)
 Get the "mobile charge" and "polarization" contributions to the electrostatic energy.

double Vfetk_qfEnergy (Vfetk *thee, int color)
 Get the "fixed charge" contribution to the electrostatic energy.

double Vfetk_lnDet (Vfetk *thee, int color, int oflag, int mflag)
 Calculate the log determinant of the specified operator.

int Vfetk_memChk (Vfetk *thee)
 Return the memory used by this structure (and its contents) in bytes.

void Vfetk_setAtomColors (Vfetk *thee)
 Transfer color (partition ID) information frmo a partitioned mesh to the atoms.

void Bmat_printHB (Bmat *thee, char *fname)
 Writes a Bmat to disk in Harwell-Boeing sparse matrix format.

int Bmat_choleskyFactor (Bmat *thee, int flag)
 Assembles the Cholesky factorization of a Bmat.

double Mat_lnDetDiag (Mat *thee)
 Returns the log(abs(det(D))) of a diagonal Mat, D.

int Vfetk_genCube (Vfetk *thee, double center[3], double length[3])
 Generates a simple cubic tetrahedral mesh.

PDE * Vfetk_PDE_ctor (Vfetk *fetk)
 Constructs the FEtk PDE object.

int Vfetk_PDE_ctor2 (PDE *thee, Vfetk *fetk)
 Intializes the FEtk PDE object.

void Vfetk_PDE_dtor (PDE **thee)
 Destroys FEtk PDE object.

void Vfetk_PDE_dtor2 (PDE *thee)
 FORTRAN stub: destroys FEtk PDE object.

void Vfetk_PDE_initAssemble (PDE *thee, int ip[], double rp[])
 Do once-per-assembly initialization.

void Vfetk_PDE_initElement (PDE *thee, int elementType, int chart, double tvx[][3], void *data)
 Do once-per-element initialization.

void Vfetk_PDE_initFace (PDE *thee, int faceType, int chart, double tnvec[])
 Do once-per-face initialization.

void Vfetk_PDE_initPoint (PDE *thee, int pointType, int chart, double txq[], double tU[], double tdU[][3])
 Do once-per-point initialization.

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.

double Vfetk_PDE_Fu_v (PDE *thee, int key, double V[], double dV[][3])
 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.

double Vfetk_PDE_DFu_wv (PDE *thee, int key, double W[], double dW[][3], double V[], double dV[][3])
 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.

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.

void Vfetk_PDE_u_D (PDE *thee, int type, int chart, double txq[], double F[])
 Evaluate the Dirichlet boundary condition at the given point.

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.

void Vfetk_PDE_bisectEdge (int dim, int dimII, int edgeType, int chart[], double vx[][3])
 Define the way manifold edges are bisected.

void Vfetk_PDE_mapBoundary (int dim, int dimII, int vertexType, int chart, double vx[3])
 Map a boundary point to some pre-defined shape.

int Vfetk_PDE_markSimplex (int dim, int dimII, int simplexType, int faceType[4], int vertexType[4], int chart[], double vx[][3], 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.

void Vfetk_PDE_oneChart (int dim, int dimII, int objType, int chart[], double vx[][3], int dimV)
 Unify the chart for different coordinate systems -- a no-op for us.

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.

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).

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.

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.

void Vfetk_readMesh (Vfetk *thee, int skey, Vio *sock)
 Read in mesh and initialize associated internal structures.

void Vfetk_dumpLocalVar ()
 Debugging routine to print out local variables used by PDE object.

int Vfetk_fillArray (Vfetk *thee, Bvec *vec, Vdata_Type type)
 Fill an array with the specified data.

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.


Detailed Description

FEtk master class (interface between FEtk and APBS).


Enumeration Type Documentation

enum eVfetk_GuessType
 

Initial guess type.

Note:
Do not change these values; they correspond to settings in FEtk
Enumeration values:
VGT_ZERO  Zero initial guess
VGT_DIRI  Dirichlet boundary condition initial guess
VGT_PREV  Previous level initial guess

enum eVfetk_LsolvType
 

Linear solver type.

Note:
Do not change these values; they correspond to settings in FEtk
Enumeration values:
VLT_SLU  SuperLU direct solve
VLT_MG  Multigrid
VLT_CG  Conjugate gradient
VLT_BCG  BiCGStab

enum eVfetk_NsolvType
 

Non-linear solver type.

Note:
Do not change these values; they correspond to settings in FEtk
Enumeration values:
VNT_NEW  Newton solver
VNT_INC  Incremental
VNT_ARC  Psuedo-arclength

enum eVfetk_PrecType
 

Preconditioner type.

Note:
Do not change these values; they correspond to settings in FEtk
Enumeration values:
VPT_IDEN  Identity matrix
VPT_DIAG  Diagonal scaling
VPT_MG  Multigrid


Function Documentation

int Bmat_choleskyFactor Bmat *    thee,
int    flag
 

Assembles the Cholesky factorization of a Bmat.

Author:
Stephen Bond
Note:
This is a friend function of Bmat
Parameters:
thee  Bmat object
flag  Type of factor to be assembled:
  • 0: Full Cholesky factor
  • 1: Diagonal of the Cholesky factor

void Bmat_printHB Bmat *    thee,
char *    fname
 

Writes a Bmat to disk in Harwell-Boeing sparse matrix format.

Author:
Stephen Bond
Note:
This is a friend function of Bmat
Parameters:
thee  Bmat object
fname  char Output filename
Bug:
Hardwired to only handle the single block symmetric case.

double Mat_lnDetDiag Mat *    thee
 

Returns the log(abs(det(D))) of a diagonal Mat, D.

Author:
Stephen Bond
Note:
This is a friend function of Mat
Parameters:
thee  Mat object
Bug:
Only works for Mat's of type RLN, CLN, or DRC.

Vfetk* Vfetk_ctor Vpbe   pbe,
Vhal_PBEType    type
 

Constructor for Vfetk object.

Author:
Nathan Baker
Parameters:
pbe  Vpbe (PBE manager) object
type  Version of PBE to solve
Returns:
Pointer to newly allocated Vfetk object
Note:
This sets up the Gem, AM, and Aprx FEtk objects but does not create a mesh. The easiest way to create a mesh is to then call Vfetk_genCube

int Vfetk_ctor2 Vfetk   thee,
Vpbe   pbe,
Vhal_PBEType    type
 

FORTRAN stub constructor for Vfetk object.

Author:
Nathan Baker
Parameters:
thee  Vfetk obeject memory address
pbe  Vpbe (PBE manager) object
type  Version of PBE to solve
Returns:
1 if successful, 0 otherwise
Note:
This sets up the Gem, AM, and Aprx FEtk objects but does not create a mesh. The easiest way to create a mesh is to then call Vfetk_genCube

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 .

Author:
Nathan Baker
Parameters:
thee  Vfetk object
color  Partition restriction for energy evaluation, only used if non-negative
Returns:
The "mobile charge" and "polarization" contributions to the electrostatic energy in units of .

void Vfetk_dtor Vfetk **    thee
 

Object destructor.

Author:
Nathan Baker
Parameters:
thee  Pointer to memory location of object to be destroyed

void Vfetk_dtor2 Vfetk   thee
 

FORTRAN stub object destructor.

Author:
Nathan Baker
Parameters:
thee  Pointer to object to be destroyed

void Vfetk_dumpLocalVar  
 

Debugging routine to print out local variables used by PDE object.

Author:
Nathan Baker
Bug:
This function is not thread-safe

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.

Author:
Nathan Baker
Parameters:
thee  Vfetk object
color  Partition restriction; if non-negative, energy calculation is restricted to the specified partition.
nonlin  If 1, the NPBE energy functional is used, if 0 then the LPBE energy functional is used.
Returns:
Total electrostatic energy in units of .

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).

Author:
Nathan Baker
Parameters:
simps  List of parent and children simplices
num  Number of simplices in list
Bug:
This function is not thread-safe.

int Vfetk_fillArray Vfetk   thee,
Bvec *    vec,
Vdata_Type    type
 

Fill an array with the specified data.

Author:
Nathan Baker
Parameters:
thee  Vfetk object
vec  FEtk Bvec vector to use
type  Type of data to write
Note:
This function is thread-safe
Bug:
Several values of type are not implemented
Returns:
1 if successful, 0 otherwise

int Vfetk_genCube Vfetk   thee,
double    center[3],
double    length[3]
 

Generates a simple cubic tetrahedral mesh.

Author:
Nathan Baker (based on Mike Holst's Gem_makeCube code)
Parameters:
thee  Vfetk object
center  desired mesh center (in Å)
length  desired mesh length (in Å)
Returns:
1 if successful, 0 otherwise

AM* Vfetk_getAM Vfetk   thee
 

Get a pointer to the AM (algebra manager) object.

Author:
Nathan Baker
Parameters:
thee  Vfetk object
Returns:
Pointer to the AM (algebra manager) object

int Vfetk_getAtomColor Vfetk   thee,
int    iatom
 

Get the partition information for a particular atom.

Author:
Nathan Baker
Note:
Friend function of Vatom
Parameters:
thee  Vfetk object
iatom  Valist atom ID
Returns:
Partition ID

Gem* Vfetk_getGem Vfetk   thee
 

Get a pointer to the Gem (grid manager) object.

Author:
Nathan Baker
Parameters:
thee  Vfetk object
Returns:
Pointer to the Gem (grid manager) object

double* Vfetk_getSolution Vfetk   thee,
int *    length
 

Create an array containing the solution (electrostatic potential in units of ) at the finest mesh level.

Author:
Nathan Baker and Michael Holst
Note:
The user is responsible for destroying the newly created array
Parameters:
thee  Vfetk object
length  Set to the length of the newly created array
Returns:
Newly created array of length "length" (see above); the user is responsible for destruction

Vcsm* Vfetk_getVcsm Vfetk   thee
 

Get a pointer to the Vcsm (charge-simplex map) object.

Author:
Nathan Baker
Parameters:
thee  Vfetk object
Returns:
Pointer to the Vcsm (charge-simplex map) object

Vpbe* Vfetk_getVpbe Vfetk   thee
 

Get a pointer to the Vpbe (PBE manager) object.

Author:
Nathan Baker
Parameters:
thee  Vfetk object
Returns:
Pointer to the Vpbe (PBE manager) object

double Vfetk_lnDet Vfetk   thee,
int    color,
int    oflag,
int    mflag
 

Calculate the log determinant of the specified operator.

Author:
Nathan Baker and Stephen Bond
Note:
  • Only works with symmetric positive definite matrices
  • Uses LU or Recycled Cholesky factorization and can be very memory- and time-consuming for large matrices.
Parameters:
thee  Vfetk object
color  Partition to evaluate over (ignored if <0)
oflag  Operator to evaluate:
  • 0: Helmholtz operator (NPBE tangent operator evaluated at zero solution)
  • 1: Response function (NPBE tangent operator evaluated at NPBE solution)
mflag  Method to use:
  • 0: Full nonsymmetric SuperLU factor with ROW/COL reordering
  • 1: Recycled symmetric Cholesky factor with no reordering
Returns:
The log determinant of the specified operator
Bug:
color argument ignored

int Vfetk_memChk Vfetk   thee
 

Return the memory used by this structure (and its contents) in bytes.

Author:
Nathan Baker
Parameters:
thee  Vfetk object
Returns:
The memory used by this structure and its contents in bytes

void Vfetk_PDE_bisectEdge int    dim,
int    dimII,
int    edgeType,
int    chart[],
double    vx[][3]
 

Define the way manifold edges are bisected.

Author:
Nathan Baker and Mike Holst
Parameters:
dim  intrinsic dimension of the manifold
dimII  imbedding dimension of the manifold
edgeType  edge type being refined
chart  Accessibility bitfields
  • chart[0] Accessibility bitfield for first vertex of edge
  • chart[1] Accessibility bitfield for second vertex of edge
  • chart[2] Set to accessibility bitfield for new vertex
vx  Vertex coordinates
  • vx[0] First vertex coordinates
  • vx[1] Second vertex coordinates
  • vx[2] Set to new vertex coordinates
Note:
This function is thread-safe.

PDE* Vfetk_PDE_ctor Vfetk   fetk
 

Constructs the FEtk PDE object.

Author:
Nathan Baker
Parameters:
fetk  Vfetk (parent) object
Returns:
Newly-allocated PDE object
Bug:
Not thread-safe

int Vfetk_PDE_ctor2 PDE *    thee,
Vfetk   fetk
 

Intializes the FEtk PDE object.

Author:
Nathan Baker (with code by Mike Holst)
Parameters:
thee  newly-allocated PDE object
fetk  Vfetk (parent) object
Returns:
1 if successful, 0 otherwise
Bug:
Not thread-safe

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.

Author:
Nathan Baker
Parameters:
thee  PDE object
chart  Chart for point coordinates
type  Vertex type
txq  Point coordinates
user  User-defined data (cast to vertex of interest)
F  set to delta function value
Bug:
This function is not thread-safe

double Vfetk_PDE_DFu_wv PDE *    thee,
int    key,
double    W[],
double    dW[][3],
double    V[],
double    dV[][3]
 

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.

Author:
Nathan Baker and Mike Holst
Parameters:
thee  PDE object
key  Integrand to evaluate
  • 0 => Interior form (see above)
  • 1 => Boundary form (not required)
W  trial function at current point
dW  trial function gradient at current point
V  test function at current point
dV  test function gradient at current point
Returns:
Integrand value
Bug:
This function is not thread-safe

void Vfetk_PDE_dtor PDE **    thee
 

Destroys FEtk PDE object.

Author:
Nathan Baker
Parameters:
thee  Pointer to PDE object memory
Note:
Thread-safe

void Vfetk_PDE_dtor2 PDE *    thee
 

FORTRAN stub: destroys FEtk PDE object.

Author:
Nathan Baker
Parameters:
thee  Pointer to PDE object
Note:
Thread-safe

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.

Author:
Nathan Baker
Parameters:
thee  PDE object
key  type of point (0 = interior, 1 = boundary, 2 = interior boundary)
F  value of residual
Bug:
This function is not thread-safe

This function is not implemented (sets error to zero)

double Vfetk_PDE_Fu_v PDE *    thee,
int    key,
double    V[],
double    dV[][3]
 

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.

Author:
Nathan Baker and Mike Holst
Parameters:
thee  PDE object
key  Integrand to evaluate
  • 0 => Interior weak form (see above)
  • 1 => Boundary weak form (not required)
V  test function at current point
dV  test function gradient at current point
Returns:
Integrand value
Bug:
This function is not thread-safe

void Vfetk_PDE_initAssemble PDE *    thee,
int    ip[],
double    rp[]
 

Do once-per-assembly initialization.

Author:
Nathan Baker and Mike Holst
Parameters:
thee  PDE object
ip  Integer parameter array (not used)
rp  Double parameter array (not used)
Note:
Thread-safe

void Vfetk_PDE_initElement PDE *    thee,
int    elementType,
int    chart,
double    tvx[][3],
void *    data
 

Do once-per-element initialization.

Author:
Nathan Baker and Mike Holst
Parameters:
thee  PDE object
elementType  Material type of element
chart  Chart in which vertex coordinates are provided (used here as bitfield for molecular accessibiliity)
tvx  Vertex coordinates
data  Simplex pointer
Todo:
Jump term is not implemented
Bug:
This function is not thread-safe

void Vfetk_PDE_initFace PDE *    thee,
int    faceType,
int    chart,
double    tnvec[]
 

Do once-per-face initialization.

Author:
Nathan Baker and Mike Holst
Parameters:
thee  PDE object
faceType  (interior or various boundary types)
chart  Chart in which vertex coordinates are provided (used here as bitfield for molecular accessibiliity)
tnvec  Coordinates of outward normal vector for face
Bug:
This function is not thread-safe

void Vfetk_PDE_initPoint PDE *    thee,
int    pointType,
int    chart,
double    txq[],
double    tU[],
double    tdU[][3]
 

Do once-per-point initialization.

Author:
Nathan Baker
Parameters:
thee  PDE object
pointType  (interior or various boundary types)
chart  Chart in which vertex coordinates are provided (used here as bitfield for molecular accessibiliity)
txq  Coordinates of point
tU  Solution value at point
tdU  Solution gradient at point
Bug:
This function is not thread-safe

This function uses pre-defined boudnary definitions for the molecular surface.

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.

Author:
Nathan Baker
Parameters:
thee  The PDE object
key  interior (0) or boundary (1)
Returns:
Energy value (in kT)
Bug:
This function is not thread-safe.

void Vfetk_PDE_mapBoundary int    dim,
int    dimII,
int    vertexType,
int    chart,
double    vx[3]
 

Map a boundary point to some pre-defined shape.

Author:
Nathan Baker and Mike Holst
Parameters:
dim  intrinsic dimension of the manifold
dimII  imbedding dimension of the manifold
vertexType  type of boundary vertex
chart  Vertex chart
vx  Vertex coordinates
Note:
This function is thread-safe and is a no-op

int Vfetk_PDE_markSimplex int    dim,
int    dimII,
int    simplexType,
int    faceType[4],
int    vertexType[4],
int    chart[],
double    vx[][3],
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.

Author:
Nathan Baker
Parameters:
dim  intrinsic dimension of the manifold
dimII  imbedding dimension of the manifold
simplexType  type of simplex
faceType  types of simplex faces
vertexType  types of simplex vertices
chart  Vertices' charts
vx  Vertex coordinates
simplex  Simplex pointer
Returns:
1 if mark simplex for refinement, 0 otherwise
Bug:
This function is not thread-safe

void Vfetk_PDE_oneChart int    dim,
int    dimII,
int    objType,
int    chart[],
double    vx[][3],
int    dimV
 

Unify the chart for different coordinate systems -- a no-op for us.

Author:
Nathan Baker
Parameters:
dim  Dimension of manifold
dimII  Imbedding dimension
objType  Don't know
chart  Charts of vertices
vx  Vertex coordinates
dimV  Number of vertices
Note:
Thread-safe; a no-op

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.

Author:
Mike Holst
Parameters:
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:
  • 0 = evaluate basis(x,y,z)
  • 1 = evaluate basis_x(x,y,z)
  • 2 = evaluate basis_y(x,y,z)
  • 3 = evaluate basis_z(x,y,z)
  • 4 = evaluate basis_xx(x,y,z)
  • 5 = evaluate basis_yy(x,y,z)
  • 6 = evaluate basis_zz(x,y,z)
  • 7 = etc...
xq  set to quadrature point coordinate
basis  set to all basis functions evaluated at all quadrature points

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.

Author:
Mike Holst
Parameters:
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 number of degrees of freedom
dof  set to degree of freedom per v/e/f/s
Note:
 *   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.
 *   

void Vfetk_PDE_u_D PDE *    thee,
int    type,
int    chart,
double    txq[],
double    F[]
 

Evaluate the Dirichlet boundary condition at the given point.

Author:
Nathan Baker
Parameters:
thee  PDE object
type  Vertex type (i.e., different types of boundary points)
chart  Chart in which vertex coordinates are provided
txq  Vertex coordinates
F  Set to boundary condition value
Bug:
This function is hard-coded to call only multiple-sphere Debye-Hü functions.

This function is not thread-safe.

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.

Author:
Nathan Baker
Parameters:
thee  PDE object
type  Vertex type (i.e., different types of boundary points)
chart  Chart in which vertex coordinates are provided
txq  Vertex coordinates
F  Set to true solution value
Note:
This function only returns zero.
Bug:
This function is not thread-safe.

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.

Author:
Nathan Baker
Parameters:
thee  Vfetk object
color  Partition restriction for energy evaluation, only used if non-negative
Returns:
The fixed charge electrostatic energy in units of .

void Vfetk_readMesh Vfetk   thee,
int    skey,
Vio *    sock
 

Read in mesh and initialize associated internal structures.

Author:
Nathan Baker
Parameters:
thee  Vfetk object
sock  Socket ready for reading
skey  Format key (0 => simplex format)
Note:
See also:
Vfetk_genCube

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!!!

Warning:
This function should only be used immediately after mesh partitioning
Author:
Nathan Baker
Note:
This is a friend function of Vcsm
Parameters:
thee  Vfetk object

void Vfetk_setParameters Vfetk   thee,
PBEparm   pbeparm,
FEMparm   feparm
 

Set the parameter objects.

Author:
Nathan Baker
Parameters:
thee  Vfetk object
pbeparm  PBE generic parameters
feparm  FEM-specific parameters

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.

Author:
Nathan Baker
Parameters:
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
Note:
This function is thread-safe
Bug:
Some values of format are not implemented
Returns:
1 if successful, 0 otherwise


Generated on Tue Feb 17 09:27:33 2004 for APBS by doxygen1.2.18