APBS 3.0.0
Loading...
Searching...
No Matches
vpmg.h File Reference

Contains declarations for class Vpmg. More...

#include "apbscfg.h"
#include "maloc/maloc.h"
#include "generic/vhal.h"
#include "generic/vacc.h"
#include "generic/vcap.h"
#include "generic/vpbe.h"
#include "generic/mgparm.h"
#include "generic/pbeparm.h"
#include "generic/vmatrix.h"
#include "pmgc/mgdrvd.h"
#include "pmgc/newdrvd.h"
#include "pmgc/mgsubd.h"
#include "pmgc/mikpckd.h"
#include "pmgc/matvecd.h"
#include "mg/vpmgp.h"
#include "mg/vgrid.h"
Include dependency graph for vpmg.h:
This graph shows which files directly or indirectly include this file:

Go to the source code of this file.

Data Structures

struct  sVpmg
 Contains public data members for Vpmg class/module. More...
 

Macros

#define VPMGMAXPART   2000
 
#define VCUB(x)   ((x)*(x)*(x))
 
#define VLOG(x)   (log(x))
 
#define IJK(i, j, k)   (((k)*(nx)*(ny))+((j)*(nx))+(i))
 
#define IJKx(j, k, i)   (((i)*(ny)*(nz))+((k)*(ny))+(j))
 
#define IJKy(i, k, j)   (((j)*(nx)*(nz))+((k)*(nx))+(i))
 
#define IJKz(i, j, k)   (((k)*(nx)*(ny))+((j)*(nx))+(i))
 
#define VFCHI(iint, iflt)   (1.5+((double)(iint)-(iflt)))
 

Typedefs

typedef struct sVpmg Vpmg
 Declaration of the Vpmg class as the Vpmg structure.
 

Functions

VEXTERNC unsigned long int Vpmg_memChk (Vpmg *thee)
 Return the memory used by this structure (and its contents) in bytes.
 
VEXTERNC VpmgVpmg_ctor (Vpmgp *parms, Vpbe *pbe, int focusFlag, Vpmg *pmgOLD, MGparm *mgparm, PBEparm_calcEnergy energyFlag)
 Constructor for the Vpmg class (allocates new memory)
 
VEXTERNC int Vpmg_ctor2 (Vpmg *thee, Vpmgp *parms, Vpbe *pbe, int focusFlag, Vpmg *pmgOLD, MGparm *mgparm, PBEparm_calcEnergy energyFlag)
 FORTRAN stub constructor for the Vpmg class (uses previously-allocated memory)
 
VEXTERNC void Vpmg_dtor (Vpmg **thee)
 Object destructor.
 
VEXTERNC void Vpmg_dtor2 (Vpmg *thee)
 FORTRAN stub object destructor.
 
VEXTERNC int Vpmg_fillco (Vpmg *thee, Vsurf_Meth surfMeth, double splineWin, Vchrg_Meth chargeMeth, int useDielXMap, Vgrid *dielXMap, int useDielYMap, Vgrid *dielYMap, int useDielZMap, Vgrid *dielZMap, int useKappaMap, Vgrid *kappaMap, int usePotMap, Vgrid *potMap, int useChargeMap, Vgrid *chargeMap)
 Fill the coefficient arrays prior to solving the equation.
 
VEXTERNC int Vpmg_solve (Vpmg *thee)
 Solve the PBE using PMG.
 
VEXTERNC int Vpmg_solveLaplace (Vpmg *thee)
 Solve Poisson's equation with a homogeneous Laplacian operator using the solvent dielectric constant. This solution is performed by a sine wave decomposition.
 
VEXTERNC double Vpmg_energy (Vpmg *thee, int extFlag)
 Get the total electrostatic energy.
 
VEXTERNC double Vpmg_qfEnergy (Vpmg *thee, int extFlag)
 Get the "fixed charge" contribution to the electrostatic energy.
 
VEXTERNC double Vpmg_qfAtomEnergy (Vpmg *thee, Vatom *atom)
 Get the per-atom "fixed charge" contribution to the electrostatic energy.
 
VEXTERNC double Vpmg_qmEnergy (Vpmg *thee, int extFlag)
 Get the "mobile charge" contribution to the electrostatic energy.
 
VEXTERNC double Vpmg_dielEnergy (Vpmg *thee, int extFlag)
 Get the "polarization" contribution to the electrostatic energy.
 
VEXTERNC double Vpmg_dielGradNorm (Vpmg *thee)
 Get the integral of the gradient of the dielectric function.
 
VEXTERNC int Vpmg_force (Vpmg *thee, double *force, int atomID, Vsurf_Meth srfm, Vchrg_Meth chgm)
 Calculate the total force on the specified atom in units of k_B T/AA.
 
VEXTERNC int Vpmg_qfForce (Vpmg *thee, double *force, int atomID, Vchrg_Meth chgm)
 Calculate the "charge-field" force on the specified atom in units of k_B T/AA.
 
VEXTERNC int Vpmg_dbForce (Vpmg *thee, double *dbForce, int atomID, Vsurf_Meth srfm)
 Calculate the dielectric boundary forces on the specified atom in units of k_B T/AA.
 
VEXTERNC int Vpmg_ibForce (Vpmg *thee, double *force, int atomID, Vsurf_Meth srfm)
 Calculate the osmotic pressure on the specified atom in units of k_B T/AA.
 
VEXTERNC void Vpmg_setPart (Vpmg *thee, double lowerCorner[3], double upperCorner[3], int bflags[6])
 Set partition information which restricts the calculation of observables to a (rectangular) subset of the problem domain.
 
VEXTERNC void Vpmg_unsetPart (Vpmg *thee)
 Remove partition restrictions.
 
VEXTERNC int Vpmg_fillArray (Vpmg *thee, double *vec, Vdata_Type type, double parm, Vhal_PBEType pbetype, PBEparm *pbeparm)
 Fill the specified array with accessibility values.
 
VPUBLIC void Vpmg_fieldSpline4 (Vpmg *thee, int atomID, double field[3])
 Computes the field at an atomic center using a stencil based on the first derivative of a 5th order B-spline.
 
VEXTERNC double Vpmg_qfPermanentMultipoleEnergy (Vpmg *thee, int atomID)
 Computes the permanent multipole electrostatic hydration energy (the polarization component of the hydration energy currently computed in TINKER).
 
VEXTERNC void Vpmg_qfPermanentMultipoleForce (Vpmg *thee, int atomID, double force[3], double torque[3])
 Computes the q-Phi Force for permanent multipoles based on 5th order B-splines.
 
VEXTERNC void Vpmg_ibPermanentMultipoleForce (Vpmg *thee, int atomID, double force[3])
 Compute the ionic boundary force for permanent multipoles.
 
VEXTERNC void Vpmg_dbPermanentMultipoleForce (Vpmg *thee, int atomID, double force[3])
 Compute the dielectric boundary force for permanent multipoles.
 
VEXTERNC void Vpmg_qfDirectPolForce (Vpmg *thee, Vgrid *perm, Vgrid *induced, int atomID, double force[3], double torque[3])
 q-Phi direct polarization force between permanent multipoles and induced dipoles, which are induced by the sum of the permanent intramolecular field and the permanent reaction field.
 
VEXTERNC void Vpmg_qfNLDirectPolForce (Vpmg *thee, Vgrid *perm, Vgrid *nlInduced, int atomID, double force[3], double torque[3])
 q-Phi direct polarization force between permanent multipoles and non-local induced dipoles based on 5th Order B-Splines. Keep in mind that the "non-local" induced dipooles are just a mathematical quantity that result from differentiation of the AMOEBA polarization energy.
 
VEXTERNC void Vpmg_ibDirectPolForce (Vpmg *thee, Vgrid *perm, Vgrid *induced, int atomID, double force[3])
 Ionic boundary direct polarization force between permanent multipoles and induced dipoles, which are induced by the sum of the permanent intramolecular field and the permanent reaction field.
 
VEXTERNC void Vpmg_ibNLDirectPolForce (Vpmg *thee, Vgrid *perm, Vgrid *nlInduced, int atomID, double force[3])
 Ionic boundary direct polarization force between permanent multipoles and non-local induced dipoles based on 5th order Keep in mind that the "non-local" induced dipooles are just a mathematical quantity that result from differentiation of the AMOEBA polarization energy.
 
VEXTERNC void Vpmg_dbDirectPolForce (Vpmg *thee, Vgrid *perm, Vgrid *induced, int atomID, double force[3])
 Dielectric boundary direct polarization force between permanent multipoles and induced dipoles, which are induced by the sum of the permanent intramolecular field and the permanent reaction field.
 
VEXTERNC void Vpmg_dbNLDirectPolForce (Vpmg *thee, Vgrid *perm, Vgrid *nlInduced, int atomID, double force[3])
 Dielectric bounday direct polarization force between permanent multipoles and non-local induced dipoles. Keep in mind that the "non-local" induced dipooles are just a mathematical quantity that result from differentiation of the AMOEBA polarization energy.
 
VEXTERNC void Vpmg_qfMutualPolForce (Vpmg *thee, Vgrid *induced, Vgrid *nlInduced, int atomID, double force[3])
 Mutual polarization force for induced dipoles based on 5th order B-Splines. This force arises due to self-consistent convergence of the solute induced dipoles and reaction field.
 
VEXTERNC void Vpmg_ibMutualPolForce (Vpmg *thee, Vgrid *induced, Vgrid *nlInduced, int atomID, double force[3])
 Ionic boundary mutual polarization force for induced dipoles based on 5th order B-Splines. This force arises due to self-consistent convergence of the solute induced dipoles and reaction field.
 
VEXTERNC void Vpmg_dbMutualPolForce (Vpmg *thee, Vgrid *induced, Vgrid *nlInduced, int atomID, double force[3])
 Dielectric boundary mutual polarization force for induced dipoles based on 5th order B-Splines. This force arises due to self-consistent convergence of the solute induced dipoles and reaction field.
 
VEXTERNC void Vpmg_printColComp (Vpmg *thee, char path[72], char title[72], char mxtype[3], int flag)
 Print out a column-compressed sparse matrix in Harwell-Boeing format.
 
VPRIVATE void bcolcomp (int *iparm, double *rparm, int *iwork, double *rwork, double *values, int *rowind, int *colptr, int *flag)
 Build a column-compressed matrix in Harwell-Boeing format.
 
VPRIVATE void bcolcomp2 (int *iparm, double *rparm, int *nx, int *ny, int *nz, int *iz, int *ipc, double *rpc, double *ac, double *cc, double *values, int *rowind, int *colptr, int *flag)
 Build a column-compressed matrix in Harwell-Boeing format.
 
VPRIVATE void bcolcomp3 (int *nx, int *ny, int *nz, int *ipc, double *rpc, double *ac, double *cc, double *values, int *rowind, int *colptr, int *flag)
 Build a column-compressed matrix in Harwell-Boeing format.
 
VPRIVATE void bcolcomp4 (int *nx, int *ny, int *nz, int *ipc, double *rpc, double *oC, double *cc, double *oE, double *oN, double *uC, double *values, int *rowind, int *colptr, int *flag)
 Build a column-compressed matrix in Harwell-Boeing format.
 
VPRIVATE void pcolcomp (int *nrow, int *ncol, int *nnzero, double *values, int *rowind, int *colptr, char *path, char *title, char *mxtype)
 Print a column-compressed matrix in Harwell-Boeing format.
 
VPRIVATE double bspline2 (double x)
 Evaluate a cubic B-spline.
 
VPRIVATE double dbspline2 (double x)
 Evaluate a cubic B-spline derivative.
 
VPRIVATE double VFCHI4 (int i, double f)
 Return 2.5 plus difference of i - f.
 
VPRIVATE double bspline4 (double x)
 Evaluate a 5th Order B-Spline (4th order polynomial)
 
VPRIVATE double dbspline4 (double x)
 Evaluate a 5th Order B-Spline derivative (4th order polynomial)
 
VPRIVATE double d2bspline4 (double x)
 Evaluate the 2nd derivative of a 5th Order B-Spline.
 
VPRIVATE double d3bspline4 (double x)
 Evaluate the 3rd derivative of a 5th Order B-Spline.
 
VPRIVATE double Vpmg_polarizEnergy (Vpmg *thee, int extFlag)
 Determines energy from polarizeable charge and interaction with fixed charges according to Rocchia et al.
 
VPRIVATE double Vpmg_qfEnergyPoint (Vpmg *thee, int extFlag)
 Calculates charge-potential energy using summation over delta function positions (i.e. something like an Linf norm)
 
VPRIVATE double Vpmg_qfEnergyVolume (Vpmg *thee, int extFlag)
 Calculates charge-potential energy as integral over a volume.
 
VPRIVATE void Vpmg_splineSelect (int srfm, Vacc *acc, double *gpos, double win, double infrad, Vatom *atom, double *force)
 Selects a spline based surface method from either VSM_SPLINE, VSM_SPLINE5 or VSM_SPLINE7.
 
VPRIVATE void bcfl1 (double size, double *apos, double charge, double xkappa, double pre1, double *gxcf, double *gycf, double *gzcf, double *xf, double *yf, double *zf, int nx, int ny, int nz)
 Increment all boundary points by pre1*(charge/d)*(exp(-xkappa*(d-size))/(1+xkappa*size) to add the effect of the Debye-Huckel potential due to a single charge.
 
VPRIVATE void multipolebc (double r, double kappa, double eps_p, double eps_w, double rad, double tsr[3])
 This routine serves bcfl2. It returns (in tsr) the contraction independent portion of the Debye-Huckel potential tensor for a spherical ion with a central charge, dipole and quadrupole. See the code for an in depth description.
 
VPRIVATE void bcCalc (Vpmg *thee)
 Fill boundary condition arrays.
 
VPRIVATE void fillcoCoef (Vpmg *thee)
 Top-level driver to fill all operator coefficient arrays.
 
VPRIVATE void fillcoCoefMap (Vpmg *thee)
 Fill operator coefficient arrays from pre-calculated maps.
 
VPRIVATE void fillcoCoefMol (Vpmg *thee)
 Fill operator coefficient arrays from a molecular surface calculation.
 
VPRIVATE void fillcoCoefMolIon (Vpmg *thee)
 Fill ion (nonlinear) operator coefficient array from a molecular surface calculation.
 
VPRIVATE void fillcoCoefMolDiel (Vpmg *thee)
 Fill differential operator coefficient arrays from a molecular surface calculation.
 
VPRIVATE void fillcoCoefMolDielNoSmooth (Vpmg *thee)
 Fill differential operator coefficient arrays from a molecular surface calculation without smoothing.
 
VPRIVATE void fillcoCoefMolDielSmooth (Vpmg *thee)
 Fill differential operator coefficient arrays from a molecular surface calculation with smoothing.
 
VPRIVATE void fillcoCoefSpline (Vpmg *thee)
 Fill operator coefficient arrays from a spline-based surface calculation.
 
VPRIVATE void fillcoCoefSpline3 (Vpmg *thee)
 Fill operator coefficient arrays from a 5th order polynomial based surface calculation.
 
VPRIVATE void fillcoCoefSpline4 (Vpmg *thee)
 Fill operator coefficient arrays from a 7th order polynomial based surface calculation.
 
VPRIVATE Vrc_Codes fillcoCharge (Vpmg *thee)
 Top-level driver to fill source term charge array.
 
VPRIVATE Vrc_Codes fillcoChargeMap (Vpmg *thee)
 Fill source term charge array from a pre-calculated map.
 
VPRIVATE void fillcoChargeSpline1 (Vpmg *thee)
 Fill source term charge array from linear interpolation.
 
VPRIVATE void fillcoChargeSpline2 (Vpmg *thee)
 Fill source term charge array from cubic spline interpolation.
 
VPRIVATE void fillcoPermanentMultipole (Vpmg *thee)
 Fill source term charge array for the use of permanent multipoles.
 
VPRIVATE void fillcoInducedDipole (Vpmg *thee)
 Fill source term charge array for use of induced dipoles.
 
VPRIVATE void fillcoNLInducedDipole (Vpmg *thee)
 Fill source term charge array for non-local induced dipoles.
 
VPRIVATE void qfForceSpline1 (Vpmg *thee, double *force, int atomID)
 Charge-field force due to a linear spline charge function.
 
VPRIVATE void qfForceSpline2 (Vpmg *thee, double *force, int atomID)
 Charge-field force due to a cubic spline charge function.
 
VPRIVATE void qfForceSpline4 (Vpmg *thee, double *force, int atomID)
 Charge-field force due to a quintic spline charge function.
 
VPRIVATE void zlapSolve (Vpmg *thee, double **solution, double **source, double **work1)
 Calculate the solution to Poisson's equation with a simple Laplacian operator and zero-valued Dirichlet boundary conditions. Store the solution in thee->u.
 
VPRIVATE void markSphere (double rtot, double *tpos, int nx, int ny, int nz, double hx, double hy, double hzed, double xmin, double ymin, double zmin, double *array, double markVal)
 Mark the grid points inside a sphere with a particular value. This marks by resetting the the grid points inside the sphere to the specified value.
 
VPRIVATE double Vpmg_qmEnergySMPBE (Vpmg *thee, int extFlag)
 Vpmg_qmEnergy for SMPBE.
 
VPRIVATE double Vpmg_qmEnergyNONLIN (Vpmg *thee, int extFlag)
 

Detailed Description

Contains declarations for class Vpmg.

Version
$Id$
Author
Nathan A. Baker
Attention
*
* APBS -- Adaptive Poisson-Boltzmann Solver
*
*  Nathan A. Baker (nathan.baker@pnnl.gov)
*  Pacific Northwest National Laboratory
*
*  Additional contributing authors listed in the code documentation.
*
* Copyright (c) 2010-2020 Battelle Memorial Institute. Developed at the
* Pacific Northwest National Laboratory, operated by Battelle Memorial
* Institute, Pacific Northwest Division for the U.S. Department of Energy.
*
* Portions Copyright (c) 2002-2010, Washington University in St. Louis.
* Portions Copyright (c) 2002-2010, Nathan A. Baker.
* Portions Copyright (c) 1999-2002, The Regents of the University of
* California.
* Portions Copyright (c) 1995, Michael Holst.
* All rights reserved.
*
* Redistribution and use in source and binary forms, with or without
* modification, are permitted provided that the following conditions are met:
*
* -  Redistributions of source code must retain the above copyright notice, this
* list of conditions and the following disclaimer.
*
* - Redistributions in binary form must reproduce the above copyright notice,
* this list of conditions and the following disclaimer in the documentation
* and/or other materials provided with the distribution.
*
* - Neither the name of Washington University in St. Louis nor the names of its
* contributors may be used to endorse or promote products derived from this
* software without specific prior written permission.
*
* THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
* "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
* LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
* A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR
* CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
* EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
* PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
* PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
* LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
* NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
* SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
*
* Neither the name of the developer nor the names of its contributors may be
* used to endorse or promote products derived from this software without
* specific prior written permission.
*
* THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
* AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
* IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
* ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE
* LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
* CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
* SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
* INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
* CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
* ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF
* THE POSSIBILITY OF SUCH DAMAGE.
*
* 

Definition in file vpmg.h.

Macro Definition Documentation

◆ IJK

#define IJK ( i,
j,
k )   (((k)*(nx)*(ny))+((j)*(nx))+(i))

Definition at line 1395 of file vpmg.h.

◆ IJKx

#define IJKx ( j,
k,
i )   (((i)*(ny)*(nz))+((k)*(ny))+(j))

Definition at line 1396 of file vpmg.h.

◆ IJKy

#define IJKy ( i,
k,
j )   (((j)*(nx)*(nz))+((k)*(nx))+(i))

Definition at line 1397 of file vpmg.h.

◆ IJKz

#define IJKz ( i,
j,
k )   (((k)*(nx)*(ny))+((j)*(nx))+(i))

Definition at line 1398 of file vpmg.h.

◆ VCUB

#define VCUB ( x)    ((x)*(x)*(x))

Definition at line 1392 of file vpmg.h.

◆ VFCHI

#define VFCHI ( iint,
iflt )   (1.5+((double)(iint)-(iflt)))

Definition at line 1399 of file vpmg.h.

◆ VLOG

#define VLOG ( x)    (log(x))

Definition at line 1393 of file vpmg.h.

◆ VPMGMAXPART

#define VPMGMAXPART   2000

Definition at line 105 of file vpmg.h.

Function Documentation

◆ bcCalc()

VPRIVATE void bcCalc ( Vpmg * thee)

Fill boundary condition arrays.

Author
Nathan Baker

Definition at line 4382 of file vpmg.c.

◆ bcfl1()

VPRIVATE void bcfl1 ( double size,
double * apos,
double charge,
double xkappa,
double pre1,
double * gxcf,
double * gycf,
double * gzcf,
double * xf,
double * yf,
double * zf,
int nx,
int ny,
int nz )

Increment all boundary points by pre1*(charge/d)*(exp(-xkappa*(d-size))/(1+xkappa*size) to add the effect of the Debye-Huckel potential due to a single charge.

Author
Nathan Baker
Parameters
aposSize of the ion
chargePosition of the ion
xkappaCharge of the ion
pre1Exponential screening factor
gxcfUnit- and dielectric-dependent prefactor
gycfSet to x-boundary values
gzcfSet to y-boundary values
xfSet to z-boundary values
yfBoundary point x-coordinates
zfBoundary point y-coordinates
nxBoundary point z-coordinates
nyNumber of grid points in x-direction
nzNumber of grid points in y-direction Number of grid points in y-direction

Definition at line 2564 of file vpmg.c.

◆ bspline2()

VPRIVATE double bspline2 ( double x)

Evaluate a cubic B-spline.

Author
Nathan Baker
Returns
Cubic B-spline value
Parameters
xPosition

Definition at line 5496 of file vpmg.c.

◆ bspline4()

VPRIVATE double bspline4 ( double x)

Evaluate a 5th Order B-Spline (4th order polynomial)

Author
: Michael Schnieders
Returns
5th Order B-Spline
Parameters
xPosition

Definition at line 7136 of file vpmg.c.

◆ d2bspline4()

VPRIVATE double d2bspline4 ( double x)

Evaluate the 2nd derivative of a 5th Order B-Spline.

Author
: Michael Schnieders
Returns
2nd derivative of a 5th Order B-Spline
Parameters
xPosition

Definition at line 7202 of file vpmg.c.

◆ d3bspline4()

VPRIVATE double d3bspline4 ( double x)

Evaluate the 3rd derivative of a 5th Order B-Spline.

Author
: Michael Schnieders
Returns
3rd derivative of a 5th Order B-Spline
Parameters
xPosition

Definition at line 7229 of file vpmg.c.

◆ dbspline2()

VPRIVATE double dbspline2 ( double x)

Evaluate a cubic B-spline derivative.

Author
Nathan Baker
Returns
Cubic B-spline derivative
Parameters
xPosition

Definition at line 5512 of file vpmg.c.

◆ dbspline4()

VPRIVATE double dbspline4 ( double x)

Evaluate a 5th Order B-Spline derivative (4th order polynomial)

Author
: Michael Schnieders
Returns
5th Order B-Spline derivative
Parameters
xPosition

Definition at line 7170 of file vpmg.c.

◆ fillcoCharge()

VPRIVATE Vrc_Codes fillcoCharge ( Vpmg * thee)

Top-level driver to fill source term charge array.

Returns
Success/failure status
Author
Nathan Baker

Definition at line 5287 of file vpmg.c.

◆ fillcoChargeMap()

VPRIVATE Vrc_Codes fillcoChargeMap ( Vpmg * thee)

Fill source term charge array from a pre-calculated map.

Returns
Success/failure status
Author
Nathan Baker

Definition at line 5343 of file vpmg.c.

◆ fillcoChargeSpline1()

VPRIVATE void fillcoChargeSpline1 ( Vpmg * thee)

Fill source term charge array from linear interpolation.

Author
Nathan Baker

Definition at line 5391 of file vpmg.c.

◆ fillcoChargeSpline2()

VPRIVATE void fillcoChargeSpline2 ( Vpmg * thee)

Fill source term charge array from cubic spline interpolation.

Author
Nathan Baker

Definition at line 5528 of file vpmg.c.

◆ fillcoCoef()

VPRIVATE void fillcoCoef ( Vpmg * thee)

Top-level driver to fill all operator coefficient arrays.

Author
Nathan Baker

Definition at line 5247 of file vpmg.c.

◆ fillcoCoefMap()

VPRIVATE void fillcoCoefMap ( Vpmg * thee)

Fill operator coefficient arrays from pre-calculated maps.

Author
Nathan Baker

Definition at line 4489 of file vpmg.c.

◆ fillcoCoefMol()

VPRIVATE void fillcoCoefMol ( Vpmg * thee)

Fill operator coefficient arrays from a molecular surface calculation.

Author
Nathan Baker

Definition at line 4612 of file vpmg.c.

◆ fillcoCoefMolDiel()

VPRIVATE void fillcoCoefMolDiel ( Vpmg * thee)

Fill differential operator coefficient arrays from a molecular surface calculation.

Author
Nathan Baker

Definition at line 4726 of file vpmg.c.

◆ fillcoCoefMolDielNoSmooth()

VPRIVATE void fillcoCoefMolDielNoSmooth ( Vpmg * thee)

Fill differential operator coefficient arrays from a molecular surface calculation without smoothing.

Author
Nathan Baker

Definition at line 4737 of file vpmg.c.

◆ fillcoCoefMolDielSmooth()

VPRIVATE void fillcoCoefMolDielSmooth ( Vpmg * thee)

Fill differential operator coefficient arrays from a molecular surface calculation with smoothing.

Molecular surface, dielectric smoothing following an implementation of Bruccoleri, et al. J Comput Chem 18 268-276 (1997).

This algorithm uses a 9 point harmonic smoothing technique - the point in question and all grid points 1/sqrt(2) grid spacings away.

Note
This uses thee->a1cf, thee->a2cf, thee->a3cf as temporary storage.
Author
Todd Dolinsky

Definition at line 4891 of file vpmg.c.

◆ fillcoCoefMolIon()

VPRIVATE void fillcoCoefMolIon ( Vpmg * thee)

Fill ion (nonlinear) operator coefficient array from a molecular surface calculation.

Author
Nathan Baker

Definition at line 4628 of file vpmg.c.

◆ fillcoCoefSpline()

VPRIVATE void fillcoCoefSpline ( Vpmg * thee)

Fill operator coefficient arrays from a spline-based surface calculation.

Author
Nathan Baker

Definition at line 5022 of file vpmg.c.

◆ fillcoCoefSpline3()

VPRIVATE void fillcoCoefSpline3 ( Vpmg * thee)

Fill operator coefficient arrays from a 5th order polynomial based surface calculation.

Author
Michael Schnieders

Definition at line 10430 of file vpmg.c.

◆ fillcoCoefSpline4()

VPRIVATE void fillcoCoefSpline4 ( Vpmg * thee)

Fill operator coefficient arrays from a 7th order polynomial based surface calculation.

Author
Michael Schnieders

Definition at line 9939 of file vpmg.c.

◆ fillcoInducedDipole()

VPRIVATE void fillcoInducedDipole ( Vpmg * thee)

Fill source term charge array for use of induced dipoles.

Author
Michael Schnieders

◆ fillcoNLInducedDipole()

VPRIVATE void fillcoNLInducedDipole ( Vpmg * thee)

Fill source term charge array for non-local induced dipoles.

Author
Michael Schnieders

◆ fillcoPermanentMultipole()

VPRIVATE void fillcoPermanentMultipole ( Vpmg * thee)

Fill source term charge array for the use of permanent multipoles.

Author
Michael Schnieders

Definition at line 7240 of file vpmg.c.

◆ markSphere()

VPRIVATE void markSphere ( double rtot,
double * tpos,
int nx,
int ny,
int nz,
double hx,
double hy,
double hzed,
double xmin,
double ymin,
double zmin,
double * array,
double markVal )

Mark the grid points inside a sphere with a particular value. This marks by resetting the the grid points inside the sphere to the specified value.

Author
Nathan Baker
Parameters
tposSphere radius
nxSphere position
nyNumber of grid points
nzNumber of grid points
hxNumber of grid points
hyGrid spacing
hzedGrid spacing
xminGrid spacing
yminGrid lower corner
zminGrid lower corner
arrayGrid lower corner
markValGrid values Value to mark with

Definition at line 6849 of file vpmg.c.

◆ multipolebc()

VPRIVATE void multipolebc ( double r,
double kappa,
double eps_p,
double eps_w,
double rad,
double tsr[3] )

This routine serves bcfl2. It returns (in tsr) the contraction independent portion of the Debye-Huckel potential tensor for a spherical ion with a central charge, dipole and quadrupole. See the code for an in depth description.

Author
Michael Schnieders
Parameters
kappaDistance to the boundary
eps_pExponential screening factor
eps_wSolute dielectric
radSolvent dielectric
tsrRadius of the sphere Contraction-independent portion of each tensor

Definition at line 3487 of file vpmg.c.

◆ qfForceSpline1()

VPRIVATE void qfForceSpline1 ( Vpmg * thee,
double * force,
int atomID )

Charge-field force due to a linear spline charge function.

Author
Nathan Baker
Parameters
atomIDSet to force Valist atom ID

Definition at line 6311 of file vpmg.c.

◆ qfForceSpline2()

VPRIVATE void qfForceSpline2 ( Vpmg * thee,
double * force,
int atomID )

Charge-field force due to a cubic spline charge function.

Author
Nathan Baker
Parameters
atomIDSet to force Valist atom ID

Definition at line 6448 of file vpmg.c.

◆ qfForceSpline4()

VPRIVATE void qfForceSpline4 ( Vpmg * thee,
double * force,
int atomID )

Charge-field force due to a quintic spline charge function.

Author
Michael Schnieders
Parameters
atomIDSet to force Valist atom ID

Definition at line 6561 of file vpmg.c.

◆ VFCHI4()

VPRIVATE double VFCHI4 ( int i,
double f )

Return 2.5 plus difference of i - f.

Author
Michael Schnieders
Returns
(2.5+((double)(i)-(f)))

Definition at line 7132 of file vpmg.c.

◆ Vpmg_polarizEnergy()

VPRIVATE double Vpmg_polarizEnergy ( Vpmg * thee,
int extFlag )

Determines energy from polarizeable charge and interaction with fixed charges according to Rocchia et al.

Author
Nathan Baker
Returns
Energy in kT
Parameters
extFlagIf 1, add external energy contributions to result

Definition at line 1148 of file vpmg.c.

◆ Vpmg_qfEnergyPoint()

VPRIVATE double Vpmg_qfEnergyPoint ( Vpmg * thee,
int extFlag )

Calculates charge-potential energy using summation over delta function positions (i.e. something like an Linf norm)

Author
Nathan Baker
Returns
Energy in kT
Parameters
extFlagIf 1, add external energy contributions to result

Definition at line 1704 of file vpmg.c.

◆ Vpmg_qfEnergyVolume()

VPRIVATE double Vpmg_qfEnergyVolume ( Vpmg * thee,
int extFlag )

Calculates charge-potential energy as integral over a volume.

Author
Nathan Baker
Returns
Energy in kT
Parameters
extFlagIf 1, add external energy contributions to result

Definition at line 1861 of file vpmg.c.

◆ Vpmg_qmEnergyNONLIN()

VPRIVATE double Vpmg_qmEnergyNONLIN ( Vpmg * thee,
int extFlag )

Definition at line 1401 of file vpmg.c.

◆ Vpmg_qmEnergySMPBE()

VPRIVATE double Vpmg_qmEnergySMPBE ( Vpmg * thee,
int extFlag )

Vpmg_qmEnergy for SMPBE.

Author
Vincent Chu

Definition at line 1490 of file vpmg.c.

◆ Vpmg_splineSelect()

VPRIVATE void Vpmg_splineSelect ( int srfm,
Vacc * acc,
double * gpos,
double win,
double infrad,
Vatom * atom,
double * force )

Selects a spline based surface method from either VSM_SPLINE, VSM_SPLINE5 or VSM_SPLINE7.

Author
David Gohara
Parameters
accSurface method, currently VSM_SPLINE, VSM_SPLINE5, or VSM_SPLINE7
gposAccessibility object
winPosition array -> array[3]
infradSpline window
atomInflation radius
forceAtom object Force array -> array[3]

Definition at line 1893 of file vpmg.c.

◆ zlapSolve()

VPRIVATE void zlapSolve ( Vpmg * thee,
double ** solution,
double ** source,
double ** work1 )

Calculate the solution to Poisson's equation with a simple Laplacian operator and zero-valued Dirichlet boundary conditions. Store the solution in thee->u.

Author
Nathan Baker
Note
Vpmg_fillco must be called first
Parameters
sourceSolution term vector
work1Source term vector Work vector

Definition at line 6898 of file vpmg.c.