Fortran and C interfaces

Laplace FMM

The Laplace FMM evaluates the following potential and its gradient

\[u(x) = \sum_{j=1}^{N} \frac{c_{j}}{\|x-x_{j}\|} - v_{j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|} \right) \, .\]

Here \(x_{j}\) are the source locations, \(c_{j}\) are the charge strengths and \(v_{j}\) are the dipole strengths, and the collection of \(x\) at which the potential and its gradient are evaluated are referred to as the evalution points.

There are 18 different Fortran wrappers for the Laplace FMM to account for collection of evaluation points (sources only, targets only, sources+targets), interaction kernel (charges only, dipoles only, charges + dipoles), output request (potential, potential+gradient).

For example, the subroutine to evaluate the potential and gradient, at a collection of targets \(t_{i}\) due to a collection of charges is:

lfmm3d_t_c_g

In general, the subroutine names take the following form:

lfmm3d_<eval-pts>_<int-ker>_<out>
  • <eval-pts>: evaluation points. Collection of x where \(u\) and its gradient is to be evaluated

    • s: Evaluate \(u\) and its gradient at the source locations \(x_{i}\)

    • t: Evaluate \(u\) and its gradient at \(t_{i}\), a collection of target locations specified by the user.

    • st: Evaluate \(u\) and its gradient at both source and target locations \(x_{i}\) and \(t_{i}\).

  • <int-ker>: kernel of interaction (charges/dipoles/both). The charge interactions are given by \(c_{j}/\|x-x_{j}\| \), and the dipole interactions are given by \(-v_{j} \cdot \nabla (1/\|x-x_{j}\|)\)

    • c: charges

    • d: dipoles

    • cd: charges + dipoles

  • <out>: Flag for evaluating potential or potential + gradient

    • p: on output only \(u\) is evaluated

    • g: on output both \(u\) and its gradient are evaluated

    • h: on output \(u\), its gradient and its hessian are evaluated

These are all the single density routines. To get a vectorized version of any of the routines use:

<subroutine name>_vec

Note

For the vectorized subroutines, the charge strengths, dipole strengths, potentials, and gradients are interleaved as opposed to provided in a sequential manner. For example for three sets of charge strengths, they should be stored as \(c_{1,1}, c_{2,1}, c_{3,1}, c_{1,2}, c_{2,2},c_{3,2} \ldots c_{1,N}, c_{2,N}, c_{3,N}\).

Example drivers:

  • examples/lfmm3d_example.f. The corresponding makefile is examples/lfmm3d_example.make

  • examples/lfmm3d_vec_example.f. The corresponding makefile is examples/lfmm3d_vec_example.make

List of interfaces

lfmm3d_s_c_p

  • Evaluation points: Sources

  • Interaction kernel: Charges

  • Outputs requested: Potential


subroutine lfmm3d_s_c_p(eps,nsource,source,charge,pot,ier)

This subroutine evaluates the potential

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{1}{\|x-x_{j}\|}\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double precision(nsource)

    Charge strengths, \(c_{j}\)

Output arguments:

  • pot: double precision(nsource)

    Potential at source locations, \(u(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_s_c_p_vec(nd,eps,nsource,source,charge,pot,ier)

This subroutine evaluates the potential

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{1}{\|x-x_{j}\|}\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double precision(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

Output arguments:

  • pot: double precision(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

lfmm3d_s_c_g

  • Evaluation points: Sources

  • Interaction kernel: Charges

  • Outputs requested: Potential and Gradient


subroutine lfmm3d_s_c_g(eps,nsource,source,charge,pot,grad,ier)

This subroutine evaluates the potential and its gradient

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{1}{\|x-x_{j}\|}\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double precision(nsource)

    Charge strengths, \(c_{j}\)

Output arguments:

  • pot: double precision(nsource)

    Potential at source locations, \(u(x_{j})\)

  • grad: double precision(3,nsource)

    Gradient at source locations, \(\nabla u(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_s_c_g_vec(nd,eps,nsource,source,charge,pot,grad,ier)

This subroutine evaluates the potential and its gradient

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{1}{\|x-x_{j}\|}\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double precision(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

Output arguments:

  • pot: double precision(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • grad: double precision(nd,3,nsource)

    Gradient at source locations, \(\nabla u_{\ell}(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

lfmm3d_s_c_h

  • Evaluation points: Sources

  • Interaction kernel: Charges

  • Outputs requested: Potential, Gradient and Hessians


subroutine lfmm3d_s_c_h(eps,nsource,source,charge,pot,grad,hess,ier)

This subroutine evaluates the potential, its gradient and its hessian

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{1}{\|x-x_{j}\|}\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double precision(nsource)

    Charge strengths, \(c_{j}\)

Output arguments:

  • pot: double precision(nsource)

    Potential at source locations, \(u(x_{j})\)

  • grad: double precision(3,nsource)

    Gradient at source locations, \(\nabla u(x_{j})\)

  • hess: double precision(6,nsource)

    Hessian at source locations, \(\nabla \nabla u(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_s_c_h_vec(nd,eps,nsource,source,charge,pot,grad,hess,ier)

This subroutine evaluates the potential, its gradient and its hessian

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{1}{\|x-x_{j}\|}\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double precision(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

Output arguments:

  • pot: double precision(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • grad: double precision(nd,3,nsource)

    Gradient at source locations, \(\nabla u_{\ell}(x_{j})\)

  • hess: double precision(nd,6,nsource)

    Gradient at source locations, \(\nabla \nabla u_{\ell}(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

lfmm3d_s_d_p

  • Evaluation points: Sources

  • Interaction kernel: Dipoles

  • Outputs requested: Potential


subroutine lfmm3d_s_d_p(eps,nsource,source,dipvec,pot,ier)

This subroutine evaluates the potential

\[u(x) = -\sum_{j=1}^{N} v_{j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • dipvec: double precision(3,nsource)

    Dipole strengths, \(v_{j}\)

Output arguments:

  • pot: double precision(nsource)

    Potential at source locations, \(u(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_s_d_p_vec(nd,eps,nsource,source,dipvec,pot,ier)

This subroutine evaluates the potential

\[u_{\ell}(x) = -\sum_{j=1}^{N} v_{\ell,j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • dipvec: double precision(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pot: double precision(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

lfmm3d_s_d_g

  • Evaluation points: Sources

  • Interaction kernel: Dipoles

  • Outputs requested: Potential and Gradient


subroutine lfmm3d_s_d_g(eps,nsource,source,dipvec,pot,grad,ier)

This subroutine evaluates the potential and its gradient

\[u(x) = -\sum_{j=1}^{N} v_{j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • dipvec: double precision(3,nsource)

    Dipole strengths, \(v_{j}\)

Output arguments:

  • pot: double precision(nsource)

    Potential at source locations, \(u(x_{j})\)

  • grad: double precision(3,nsource)

    Gradient at source locations, \(\nabla u(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_s_d_g_vec(nd,eps,nsource,source,dipvec,pot,grad,ier)

This subroutine evaluates the potential and its gradient

\[u_{\ell}(x) = -\sum_{j=1}^{N} v_{\ell,j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • dipvec: double precision(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pot: double precision(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • grad: double precision(nd,3,nsource)

    Gradient at source locations, \(\nabla u_{\ell}(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

lfmm3d_s_d_h

  • Evaluation points: Sources

  • Interaction kernel: Dipoles

  • Outputs requested: Potential, Gradient and Hessians


subroutine lfmm3d_s_d_h(eps,nsource,source,dipvec,pot,grad,hess,ier)

This subroutine evaluates the potential, its gradient and its hessian

\[u(x) = -\sum_{j=1}^{N} v_{j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • dipvec: double precision(3,nsource)

    Dipole strengths, \(v_{j}\)

Output arguments:

  • pot: double precision(nsource)

    Potential at source locations, \(u(x_{j})\)

  • grad: double precision(3,nsource)

    Gradient at source locations, \(\nabla u(x_{j})\)

  • hess: double precision(6,nsource)

    Hessian at source locations, \(\nabla \nabla u(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_s_d_h_vec(nd,eps,nsource,source,dipvec,pot,grad,hess,ier)

This subroutine evaluates the potential, its gradient and its hessian

\[u_{\ell}(x) = -\sum_{j=1}^{N} v_{\ell,j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • dipvec: double precision(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pot: double precision(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • grad: double precision(nd,3,nsource)

    Gradient at source locations, \(\nabla u_{\ell}(x_{j})\)

  • hess: double precision(nd,6,nsource)

    Gradient at source locations, \(\nabla \nabla u_{\ell}(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

lfmm3d_s_cd_p

  • Evaluation points: Sources

  • Interaction kernel: Charges and Dipoles

  • Outputs requested: Potential


subroutine lfmm3d_s_cd_p(eps,nsource,source,charge,dipvec,pot,ier)

This subroutine evaluates the potential

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{1}{\|x-x_{j}\|} - v_{j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double precision(nsource)

    Charge strengths, \(c_{j}\)

  • dipvec: double precision(3,nsource)

    Dipole strengths, \(v_{j}\)

Output arguments:

  • pot: double precision(nsource)

    Potential at source locations, \(u(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_s_cd_p_vec(nd,eps,nsource,source,charge,dipvec,pot,ier)

This subroutine evaluates the potential

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{1}{\|x-x_{j}\|} - v_{\ell,j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double precision(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

  • dipvec: double precision(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pot: double precision(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

lfmm3d_s_cd_g

  • Evaluation points: Sources

  • Interaction kernel: Charges and Dipoles

  • Outputs requested: Potential and Gradient


subroutine lfmm3d_s_cd_g(eps,nsource,source,charge,dipvec,pot,grad,ier)

This subroutine evaluates the potential and its gradient

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{1}{\|x-x_{j}\|} - v_{j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double precision(nsource)

    Charge strengths, \(c_{j}\)

  • dipvec: double precision(3,nsource)

    Dipole strengths, \(v_{j}\)

Output arguments:

  • pot: double precision(nsource)

    Potential at source locations, \(u(x_{j})\)

  • grad: double precision(3,nsource)

    Gradient at source locations, \(\nabla u(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_s_cd_g_vec(nd,eps,nsource,source,charge,dipvec,pot,grad,ier)

This subroutine evaluates the potential and its gradient

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{1}{\|x-x_{j}\|} - v_{\ell,j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double precision(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

  • dipvec: double precision(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pot: double precision(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • grad: double precision(nd,3,nsource)

    Gradient at source locations, \(\nabla u_{\ell}(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

lfmm3d_s_cd_h

  • Evaluation points: Sources

  • Interaction kernel: Charges and Dipoles

  • Outputs requested: Potential, Gradient and Hessians


subroutine lfmm3d_s_cd_h(eps,nsource,source,charge,dipvec,pot,grad,hess,ier)

This subroutine evaluates the potential, its gradient and its hessian

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{1}{\|x-x_{j}\|} - v_{j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double precision(nsource)

    Charge strengths, \(c_{j}\)

  • dipvec: double precision(3,nsource)

    Dipole strengths, \(v_{j}\)

Output arguments:

  • pot: double precision(nsource)

    Potential at source locations, \(u(x_{j})\)

  • grad: double precision(3,nsource)

    Gradient at source locations, \(\nabla u(x_{j})\)

  • hess: double precision(6,nsource)

    Hessian at source locations, \(\nabla \nabla u(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_s_cd_h_vec(nd,eps,nsource,source,charge,dipvec,pot,grad,hess,ier)

This subroutine evaluates the potential, its gradient and its hessian

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{1}{\|x-x_{j}\|} - v_{\ell,j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double precision(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

  • dipvec: double precision(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pot: double precision(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • grad: double precision(nd,3,nsource)

    Gradient at source locations, \(\nabla u_{\ell}(x_{j})\)

  • hess: double precision(nd,6,nsource)

    Gradient at source locations, \(\nabla \nabla u_{\ell}(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

lfmm3d_t_c_p

  • Evaluation points: Targets

  • Interaction kernel: Charges

  • Outputs requested: Potential


subroutine lfmm3d_t_c_p(eps,nsource,source,charge,ntarg,targ,pottarg,ier)

This subroutine evaluates the potential

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{1}{\|x-x_{j}\|}\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double precision(nsource)

    Charge strengths, \(c_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pottarg: double precision(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_t_c_p_vec(nd,eps,nsource,source,charge,ntarg,targ,pottarg,ier)

This subroutine evaluates the potential

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{1}{\|x-x_{j}\|}\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double precision(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

Output arguments:

  • pottarg: double precision(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

lfmm3d_t_c_g

  • Evaluation points: Targets

  • Interaction kernel: Charges

  • Outputs requested: Potential and Gradient


subroutine lfmm3d_t_c_g(eps,nsource,source,charge,ntarg,targ,pottarg,gradtarg,ier)

This subroutine evaluates the potential and its gradient

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{1}{\|x-x_{j}\|}\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double precision(nsource)

    Charge strengths, \(c_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pottarg: double precision(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • gradtarg: double precision(3,ntarg)

    Gradient at target locations, \(\nabla u(t_{i})\)

  • hesstarg: double precision(6,ntarg)

    Hessian at target locations, \(\nabla \nabla u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_t_c_g_vec(nd,eps,nsource,source,charge,ntarg,targ,pottarg,gradtarg,ier)

This subroutine evaluates the potential and its gradient

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{1}{\|x-x_{j}\|}\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double precision(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

Output arguments:

  • pottarg: double precision(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • gradtarg: double precision(nd,3,ntarg)

    Gradient at target locations, \(\nabla u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

lfmm3d_t_c_h

  • Evaluation points: Targets

  • Interaction kernel: Charges

  • Outputs requested: Potential, Gradient and Hessians


subroutine lfmm3d_t_c_h(eps,nsource,source,charge,ntarg,targ,pottarg,gradtarg,hesstarg,ier)

This subroutine evaluates the potential, its gradient and its hessian

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{1}{\|x-x_{j}\|}\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double precision(nsource)

    Charge strengths, \(c_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pottarg: double precision(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • gradtarg: double precision(3,ntarg)

    Gradient at target locations, \(\nabla u(t_{i})\)

  • hesstarg: double precision(6,ntarg)

    Hessian at target locations, \(\nabla \nabla u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_t_c_h_vec(nd,eps,nsource,source,charge,ntarg,targ,pottarg,gradtarg,hesstarg,ier)

This subroutine evaluates the potential, its gradient and its hessian

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{1}{\|x-x_{j}\|}\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double precision(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

Output arguments:

  • pottarg: double precision(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • gradtarg: double precision(nd,3,ntarg)

    Gradient at target locations, \(\nabla u_{\ell}(t_{i})\)

  • hesstarg: double precision(nd,3,ntarg)

    Hessian at target locations, \(\nabla \nabla u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

lfmm3d_t_d_p

  • Evaluation points: Targets

  • Interaction kernel: Dipoles

  • Outputs requested: Potential


subroutine lfmm3d_t_d_p(eps,nsource,source,dipvec,ntarg,targ,pottarg,ier)

This subroutine evaluates the potential

\[u(x) = -\sum_{j=1}^{N} v_{j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • dipvec: double precision(3,nsource)

    Dipole strengths, \(v_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pottarg: double precision(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_t_d_p_vec(nd,eps,nsource,source,dipvec,ntarg,targ,pottarg,ier)

This subroutine evaluates the potential

\[u_{\ell}(x) = -\sum_{j=1}^{N} v_{\ell,j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • dipvec: double precision(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pottarg: double precision(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

lfmm3d_t_d_g

  • Evaluation points: Targets

  • Interaction kernel: Dipoles

  • Outputs requested: Potential and Gradient


subroutine lfmm3d_t_d_g(eps,nsource,source,dipvec,ntarg,targ,pottarg,gradtarg,ier)

This subroutine evaluates the potential and its gradient

\[u(x) = -\sum_{j=1}^{N} v_{j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • dipvec: double precision(3,nsource)

    Dipole strengths, \(v_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pottarg: double precision(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • gradtarg: double precision(3,ntarg)

    Gradient at target locations, \(\nabla u(t_{i})\)

  • hesstarg: double precision(6,ntarg)

    Hessian at target locations, \(\nabla \nabla u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_t_d_g_vec(nd,eps,nsource,source,dipvec,ntarg,targ,pottarg,gradtarg,ier)

This subroutine evaluates the potential and its gradient

\[u_{\ell}(x) = -\sum_{j=1}^{N} v_{\ell,j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • dipvec: double precision(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pottarg: double precision(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • gradtarg: double precision(nd,3,ntarg)

    Gradient at target locations, \(\nabla u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

lfmm3d_t_d_h

  • Evaluation points: Targets

  • Interaction kernel: Dipoles

  • Outputs requested: Potential, Gradient and Hessians


subroutine lfmm3d_t_d_h(eps,nsource,source,dipvec,ntarg,targ,pottarg,gradtarg,hesstarg,ier)

This subroutine evaluates the potential, its gradient and its hessian

\[u(x) = -\sum_{j=1}^{N} v_{j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • dipvec: double precision(3,nsource)

    Dipole strengths, \(v_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pottarg: double precision(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • gradtarg: double precision(3,ntarg)

    Gradient at target locations, \(\nabla u(t_{i})\)

  • hesstarg: double precision(6,ntarg)

    Hessian at target locations, \(\nabla \nabla u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_t_d_h_vec(nd,eps,nsource,source,dipvec,ntarg,targ,pottarg,gradtarg,hesstarg,ier)

This subroutine evaluates the potential, its gradient and its hessian

\[u_{\ell}(x) = -\sum_{j=1}^{N} v_{\ell,j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • dipvec: double precision(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pottarg: double precision(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • gradtarg: double precision(nd,3,ntarg)

    Gradient at target locations, \(\nabla u_{\ell}(t_{i})\)

  • hesstarg: double precision(nd,3,ntarg)

    Hessian at target locations, \(\nabla \nabla u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

lfmm3d_t_cd_p

  • Evaluation points: Targets

  • Interaction kernel: Charges and Dipoles

  • Outputs requested: Potential


subroutine lfmm3d_t_cd_p(eps,nsource,source,charge,dipvec,ntarg,targ,pottarg,ier)

This subroutine evaluates the potential

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{1}{\|x-x_{j}\|} - v_{j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double precision(nsource)

    Charge strengths, \(c_{j}\)

  • dipvec: double precision(3,nsource)

    Dipole strengths, \(v_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pottarg: double precision(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_t_cd_p_vec(nd,eps,nsource,source,charge,dipvec,ntarg,targ,pottarg,ier)

This subroutine evaluates the potential

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{1}{\|x-x_{j}\|} - v_{\ell,j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double precision(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

  • dipvec: double precision(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pottarg: double precision(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

lfmm3d_t_cd_g

  • Evaluation points: Targets

  • Interaction kernel: Charges and Dipoles

  • Outputs requested: Potential and Gradient


subroutine lfmm3d_t_cd_g(eps,nsource,source,charge,dipvec,ntarg,targ,pottarg,gradtarg,ier)

This subroutine evaluates the potential and its gradient

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{1}{\|x-x_{j}\|} - v_{j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double precision(nsource)

    Charge strengths, \(c_{j}\)

  • dipvec: double precision(3,nsource)

    Dipole strengths, \(v_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pottarg: double precision(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • gradtarg: double precision(3,ntarg)

    Gradient at target locations, \(\nabla u(t_{i})\)

  • hesstarg: double precision(6,ntarg)

    Hessian at target locations, \(\nabla \nabla u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_t_cd_g_vec(nd,eps,nsource,source,charge,dipvec,ntarg,targ,pottarg,gradtarg,ier)

This subroutine evaluates the potential and its gradient

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{1}{\|x-x_{j}\|} - v_{\ell,j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double precision(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

  • dipvec: double precision(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pottarg: double precision(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • gradtarg: double precision(nd,3,ntarg)

    Gradient at target locations, \(\nabla u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

lfmm3d_t_cd_h

  • Evaluation points: Targets

  • Interaction kernel: Charges and Dipoles

  • Outputs requested: Potential, Gradient and Hessians


subroutine lfmm3d_t_cd_h(eps,nsource,source,charge,dipvec,ntarg,targ,pottarg,gradtarg,hesstarg,ier)

This subroutine evaluates the potential, its gradient and its hessian

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{1}{\|x-x_{j}\|} - v_{j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double precision(nsource)

    Charge strengths, \(c_{j}\)

  • dipvec: double precision(3,nsource)

    Dipole strengths, \(v_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pottarg: double precision(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • gradtarg: double precision(3,ntarg)

    Gradient at target locations, \(\nabla u(t_{i})\)

  • hesstarg: double precision(6,ntarg)

    Hessian at target locations, \(\nabla \nabla u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_t_cd_h_vec(nd,eps,nsource,source,charge,dipvec,ntarg,targ,pottarg,gradtarg,hesstarg,ier)

This subroutine evaluates the potential, its gradient and its hessian

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{1}{\|x-x_{j}\|} - v_{\ell,j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double precision(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

  • dipvec: double precision(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pottarg: double precision(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • gradtarg: double precision(nd,3,ntarg)

    Gradient at target locations, \(\nabla u_{\ell}(t_{i})\)

  • hesstarg: double precision(nd,3,ntarg)

    Hessian at target locations, \(\nabla \nabla u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

lfmm3d_st_c_p

  • Evaluation points: Sources and Targets

  • Interaction kernel: Charges

  • Outputs requested: Potential


subroutine lfmm3d_st_c_p(eps,nsource,source,charge,pot,ntarg,targ,pottarg,ier)

This subroutine evaluates the potential

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{1}{\|x-x_{j}\|}\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double precision(nsource)

    Charge strengths, \(c_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pot: double precision(nsource)

    Potential at source locations, \(u(x_{j})\)

  • pottarg: double precision(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_st_c_p_vec(nd,eps,nsource,source,charge,pot,ntarg,targ,pottarg,ier)

This subroutine evaluates the potential

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{1}{\|x-x_{j}\|}\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double precision(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

Output arguments:

  • pot: double precision(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • pottarg: double precision(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

lfmm3d_st_c_g

  • Evaluation points: Sources and Targets

  • Interaction kernel: Charges

  • Outputs requested: Potential and Gradient


subroutine lfmm3d_st_c_g(eps,nsource,source,charge,pot,grad,ntarg,targ,pottarg,gradtarg,ier)

This subroutine evaluates the potential and its gradient

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{1}{\|x-x_{j}\|}\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double precision(nsource)

    Charge strengths, \(c_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pot: double precision(nsource)

    Potential at source locations, \(u(x_{j})\)

  • grad: double precision(3,nsource)

    Gradient at source locations, \(\nabla u(x_{j})\)

  • pottarg: double precision(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • gradtarg: double precision(3,ntarg)

    Gradient at target locations, \(\nabla u(t_{i})\)

  • hesstarg: double precision(6,ntarg)

    Hessian at target locations, \(\nabla \nabla u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_st_c_g_vec(nd,eps,nsource,source,charge,pot,grad,ntarg,targ,pottarg,gradtarg,ier)

This subroutine evaluates the potential and its gradient

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{1}{\|x-x_{j}\|}\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double precision(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

Output arguments:

  • pot: double precision(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • grad: double precision(nd,3,nsource)

    Gradient at source locations, \(\nabla u_{\ell}(x_{j})\)

  • pottarg: double precision(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • gradtarg: double precision(nd,3,ntarg)

    Gradient at target locations, \(\nabla u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

lfmm3d_st_c_h

  • Evaluation points: Sources and Targets

  • Interaction kernel: Charges

  • Outputs requested: Potential, Gradient and Hessians


subroutine lfmm3d_st_c_h(eps,nsource,source,charge,pot,grad,hess,ntarg,targ,pottarg,gradtarg,hesstarg,ier)

This subroutine evaluates the potential, its gradient and its hessian

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{1}{\|x-x_{j}\|}\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double precision(nsource)

    Charge strengths, \(c_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pot: double precision(nsource)

    Potential at source locations, \(u(x_{j})\)

  • grad: double precision(3,nsource)

    Gradient at source locations, \(\nabla u(x_{j})\)

  • hess: double precision(6,nsource)

    Hessian at source locations, \(\nabla \nabla u(x_{j})\)

  • pottarg: double precision(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • gradtarg: double precision(3,ntarg)

    Gradient at target locations, \(\nabla u(t_{i})\)

  • hesstarg: double precision(6,ntarg)

    Hessian at target locations, \(\nabla \nabla u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_st_c_h_vec(nd,eps,nsource,source,charge,pot,grad,hess,ntarg,targ,pottarg,gradtarg,hesstarg,ier)

This subroutine evaluates the potential, its gradient and its hessian

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{1}{\|x-x_{j}\|}\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double precision(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

Output arguments:

  • pot: double precision(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • grad: double precision(nd,3,nsource)

    Gradient at source locations, \(\nabla u_{\ell}(x_{j})\)

  • hess: double precision(nd,6,nsource)

    Gradient at source locations, \(\nabla \nabla u_{\ell}(x_{j})\)

  • pottarg: double precision(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • gradtarg: double precision(nd,3,ntarg)

    Gradient at target locations, \(\nabla u_{\ell}(t_{i})\)

  • hesstarg: double precision(nd,3,ntarg)

    Hessian at target locations, \(\nabla \nabla u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

lfmm3d_st_d_p

  • Evaluation points: Sources and Targets

  • Interaction kernel: Dipoles

  • Outputs requested: Potential


subroutine lfmm3d_st_d_p(eps,nsource,source,dipvec,pot,ntarg,targ,pottarg,ier)

This subroutine evaluates the potential

\[u(x) = -\sum_{j=1}^{N} v_{j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • dipvec: double precision(3,nsource)

    Dipole strengths, \(v_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pot: double precision(nsource)

    Potential at source locations, \(u(x_{j})\)

  • pottarg: double precision(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_st_d_p_vec(nd,eps,nsource,source,dipvec,pot,ntarg,targ,pottarg,ier)

This subroutine evaluates the potential

\[u_{\ell}(x) = -\sum_{j=1}^{N} v_{\ell,j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • dipvec: double precision(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pot: double precision(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • pottarg: double precision(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

lfmm3d_st_d_g

  • Evaluation points: Sources and Targets

  • Interaction kernel: Dipoles

  • Outputs requested: Potential and Gradient


subroutine lfmm3d_st_d_g(eps,nsource,source,dipvec,pot,grad,ntarg,targ,pottarg,gradtarg,ier)

This subroutine evaluates the potential and its gradient

\[u(x) = -\sum_{j=1}^{N} v_{j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • dipvec: double precision(3,nsource)

    Dipole strengths, \(v_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pot: double precision(nsource)

    Potential at source locations, \(u(x_{j})\)

  • grad: double precision(3,nsource)

    Gradient at source locations, \(\nabla u(x_{j})\)

  • pottarg: double precision(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • gradtarg: double precision(3,ntarg)

    Gradient at target locations, \(\nabla u(t_{i})\)

  • hesstarg: double precision(6,ntarg)

    Hessian at target locations, \(\nabla \nabla u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_st_d_g_vec(nd,eps,nsource,source,dipvec,pot,grad,ntarg,targ,pottarg,gradtarg,ier)

This subroutine evaluates the potential and its gradient

\[u_{\ell}(x) = -\sum_{j=1}^{N} v_{\ell,j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • dipvec: double precision(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pot: double precision(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • grad: double precision(nd,3,nsource)

    Gradient at source locations, \(\nabla u_{\ell}(x_{j})\)

  • pottarg: double precision(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • gradtarg: double precision(nd,3,ntarg)

    Gradient at target locations, \(\nabla u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

lfmm3d_st_d_h

  • Evaluation points: Sources and Targets

  • Interaction kernel: Dipoles

  • Outputs requested: Potential, Gradient and Hessians


subroutine lfmm3d_st_d_h(eps,nsource,source,dipvec,pot,grad,hess,ntarg,targ,pottarg,gradtarg,hesstarg,ier)

This subroutine evaluates the potential, its gradient and its hessian

\[u(x) = -\sum_{j=1}^{N} v_{j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • dipvec: double precision(3,nsource)

    Dipole strengths, \(v_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pot: double precision(nsource)

    Potential at source locations, \(u(x_{j})\)

  • grad: double precision(3,nsource)

    Gradient at source locations, \(\nabla u(x_{j})\)

  • hess: double precision(6,nsource)

    Hessian at source locations, \(\nabla \nabla u(x_{j})\)

  • pottarg: double precision(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • gradtarg: double precision(3,ntarg)

    Gradient at target locations, \(\nabla u(t_{i})\)

  • hesstarg: double precision(6,ntarg)

    Hessian at target locations, \(\nabla \nabla u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_st_d_h_vec(nd,eps,nsource,source,dipvec,pot,grad,hess,ntarg,targ,pottarg,gradtarg,hesstarg,ier)

This subroutine evaluates the potential, its gradient and its hessian

\[u_{\ell}(x) = -\sum_{j=1}^{N} v_{\ell,j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • dipvec: double precision(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pot: double precision(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • grad: double precision(nd,3,nsource)

    Gradient at source locations, \(\nabla u_{\ell}(x_{j})\)

  • hess: double precision(nd,6,nsource)

    Gradient at source locations, \(\nabla \nabla u_{\ell}(x_{j})\)

  • pottarg: double precision(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • gradtarg: double precision(nd,3,ntarg)

    Gradient at target locations, \(\nabla u_{\ell}(t_{i})\)

  • hesstarg: double precision(nd,3,ntarg)

    Hessian at target locations, \(\nabla \nabla u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

lfmm3d_st_cd_p

  • Evaluation points: Sources and Targets

  • Interaction kernel: Charges and Dipoles

  • Outputs requested: Potential


subroutine lfmm3d_st_cd_p(eps,nsource,source,charge,dipvec,pot,ntarg,targ,pottarg,ier)

This subroutine evaluates the potential

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{1}{\|x-x_{j}\|} - v_{j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double precision(nsource)

    Charge strengths, \(c_{j}\)

  • dipvec: double precision(3,nsource)

    Dipole strengths, \(v_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pot: double precision(nsource)

    Potential at source locations, \(u(x_{j})\)

  • pottarg: double precision(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_st_cd_p_vec(nd,eps,nsource,source,charge,dipvec,pot,ntarg,targ,pottarg,ier)

This subroutine evaluates the potential

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{1}{\|x-x_{j}\|} - v_{\ell,j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double precision(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

  • dipvec: double precision(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pot: double precision(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • pottarg: double precision(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

lfmm3d_st_cd_g

  • Evaluation points: Sources and Targets

  • Interaction kernel: Charges and Dipoles

  • Outputs requested: Potential and Gradient


subroutine lfmm3d_st_cd_g(eps,nsource,source,charge,dipvec,pot,grad,ntarg,targ,pottarg,gradtarg,ier)

This subroutine evaluates the potential and its gradient

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{1}{\|x-x_{j}\|} - v_{j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double precision(nsource)

    Charge strengths, \(c_{j}\)

  • dipvec: double precision(3,nsource)

    Dipole strengths, \(v_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pot: double precision(nsource)

    Potential at source locations, \(u(x_{j})\)

  • grad: double precision(3,nsource)

    Gradient at source locations, \(\nabla u(x_{j})\)

  • pottarg: double precision(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • gradtarg: double precision(3,ntarg)

    Gradient at target locations, \(\nabla u(t_{i})\)

  • hesstarg: double precision(6,ntarg)

    Hessian at target locations, \(\nabla \nabla u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_st_cd_g_vec(nd,eps,nsource,source,charge,dipvec,pot,grad,ntarg,targ,pottarg,gradtarg,ier)

This subroutine evaluates the potential and its gradient

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{1}{\|x-x_{j}\|} - v_{\ell,j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double precision(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

  • dipvec: double precision(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pot: double precision(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • grad: double precision(nd,3,nsource)

    Gradient at source locations, \(\nabla u_{\ell}(x_{j})\)

  • pottarg: double precision(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • gradtarg: double precision(nd,3,ntarg)

    Gradient at target locations, \(\nabla u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

lfmm3d_st_cd_h

  • Evaluation points: Sources and Targets

  • Interaction kernel: Charges and Dipoles

  • Outputs requested: Potential, Gradient and Hessians


subroutine lfmm3d_st_cd_h(eps,nsource,source,charge,dipvec,pot,grad,hess,ntarg,targ,pottarg,gradtarg,hesstarg,ier)

This subroutine evaluates the potential, its gradient and its hessian

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{1}{\|x-x_{j}\|} - v_{j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double precision(nsource)

    Charge strengths, \(c_{j}\)

  • dipvec: double precision(3,nsource)

    Dipole strengths, \(v_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pot: double precision(nsource)

    Potential at source locations, \(u(x_{j})\)

  • grad: double precision(3,nsource)

    Gradient at source locations, \(\nabla u(x_{j})\)

  • hess: double precision(6,nsource)

    Hessian at source locations, \(\nabla \nabla u(x_{j})\)

  • pottarg: double precision(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • gradtarg: double precision(3,ntarg)

    Gradient at target locations, \(\nabla u(t_{i})\)

  • hesstarg: double precision(6,ntarg)

    Hessian at target locations, \(\nabla \nabla u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine lfmm3d_st_cd_h_vec(nd,eps,nsource,source,charge,dipvec,pot,grad,hess,ntarg,targ,pottarg,gradtarg,hesstarg,ier)

This subroutine evaluates the potential, its gradient and its hessian

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{1}{\|x-x_{j}\|} - v_{\ell,j} \cdot \nabla \left( \frac{1}{\|x-x_{j}\|}\right)\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double precision(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

  • dipvec: double precision(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pot: double precision(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • grad: double precision(nd,3,nsource)

    Gradient at source locations, \(\nabla u_{\ell}(x_{j})\)

  • hess: double precision(nd,6,nsource)

    Gradient at source locations, \(\nabla \nabla u_{\ell}(x_{j})\)

  • pottarg: double precision(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • gradtarg: double precision(nd,3,ntarg)

    Gradient at target locations, \(\nabla u_{\ell}(t_{i})\)

  • hesstarg: double precision(nd,3,ntarg)

    Hessian at target locations, \(\nabla \nabla u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

Helmholtz FMM

The Helmholtz FMM evaluates the following potential and its gradient

\[u(x) = \sum_{j=1}^{N} \frac{c_{j} e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|} - v_{j} \cdot \nabla \left( \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|} \right) \, .\]

Here \(x_{j}\) are the source locations, \(c_{j}\) are the charge strengths and \(v_{j}\) are the dipole strengths, and the collection of \(x\) at which the potential and its gradient are evaluated are referred to as the evalution points.

There are 18 different Fortran wrappers for the Helmholtz FMM to account for collection of evaluation points (sources only, targets only, sources+targets), interaction kernel (charges only, dipoles only, charges + dipoles), output request (potential, potential+gradient).

For example, the subroutine to evaluate the potential and gradient, at a collection of targets \(t_{i}\) due to a collection of charges is:

hfmm3d_t_c_g

In general, the subroutine names take the following form:

hfmm3d_<eval-pts>_<int-ker>_<out>
  • <eval-pts>: evaluation points. Collection of x where \(u\) and its gradient is to be evaluated

    • s: Evaluate \(u\) and its gradient at the source locations \(x_{i}\)

    • t: Evaluate \(u\) and its gradient at \(t_{i}\), a collection of target locations specified by the user.

    • st: Evaluate \(u\) and its gradient at both source and target locations \(x_{i}\) and \(t_{i}\).

  • <int-ker>: kernel of interaction (charges/dipoles/both). The charge interactions are given by \(c_{j}/\|x-x_{j}\| \), and the dipole interactions are given by \(-v_{j} \cdot \nabla (1/\|x-x_{j}\|)\)

    • c: charges

    • d: dipoles

    • cd: charges + dipoles

  • <out>: Flag for evaluating potential or potential + gradient

    • p: on output only \(u\) is evaluated

    • g: on output both \(u\) and its gradient are evaluated

These are all the single density routines. To get a vectorized version of any of the routines use:

<subroutine name>_vec

Note

For the vectorized subroutines, the charge strengths, dipole strengths, potentials, and gradients are interleaved as opposed to provided in a sequential manner. For example for three sets of charge strengths, they should be stored as \(c_{1,1}, c_{2,1}, c_{3,1}, c_{1,2}, c_{2,2},c_{3,2} \ldots c_{1,N}, c_{2,N}, c_{3,N}\).

Example drivers:

  • examples/hfmm3d_example.f. The corresponding makefile is examples/hfmm3d_example.make

  • examples/hfmm3d_vec_example.f. The corresponding makefile is examples/hfmm3d_vec_example.make

List of interfaces

hfmm3d_s_c_p

  • Evaluation points: Sources

  • Interaction kernel: Charges

  • Outputs requested: Potential


subroutine hfmm3d_s_c_p(eps,zk,nsource,source,charge,pot,ier)

This subroutine evaluates the potential

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{e^{ik\|x- x_{j}\|}}{\|x-x_{j}\|}\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • zk: double complex

    Helmholtz parameter, k

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double complex(nsource)

    Charge strengths, \(c_{j}\)

Output arguments:

  • pot: double complex(nsource)

    Potential at source locations, \(u(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine hfmm3d_s_c_p_vec(nd,eps,zk,nsource,source,charge,pot,ier)

This subroutine evaluates the potential

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{e^{ik\|x- x_{j}\|}}{\|x-x_{j}\|}\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double complex(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

Output arguments:

  • pot: double complex(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

hfmm3d_s_c_g

  • Evaluation points: Sources

  • Interaction kernel: Charges

  • Outputs requested: Potential and Gradient


subroutine hfmm3d_s_c_g(eps,zk,nsource,source,charge,pot,grad,ier)

This subroutine evaluates the potential and its gradient

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{e^{ik\|x- x_{j}\|}}{\|x-x_{j}\|}\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • zk: double complex

    Helmholtz parameter, k

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double complex(nsource)

    Charge strengths, \(c_{j}\)

Output arguments:

  • pot: double complex(nsource)

    Potential at source locations, \(u(x_{j})\)

  • grad: double complex(3,nsource)

    Gradient at source locations, \(\nabla u(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine hfmm3d_s_c_g_vec(nd,eps,zk,nsource,source,charge,pot,grad,ier)

This subroutine evaluates the potential and its gradient

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{e^{ik\|x- x_{j}\|}}{\|x-x_{j}\|}\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double complex(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

Output arguments:

  • pot: double complex(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • grad: double complex(nd,3,nsource)

    Gradient at source locations, \(\nabla u_{\ell}(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

hfmm3d_s_d_p

  • Evaluation points: Sources

  • Interaction kernel: Dipoles

  • Outputs requested: Potential


subroutine hfmm3d_s_d_p(eps,zk,nsource,source,dipvec,pot,ier)

This subroutine evaluates the potential

\[u(x) = -\sum_{j=1}^{N} v_{j} \cdot \nabla \left( \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|}\right)\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • zk: double complex

    Helmholtz parameter, k

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • dipvec: double complex(3,nsource)

    Dipole strengths, \(v_{j}\)

Output arguments:

  • pot: double complex(nsource)

    Potential at source locations, \(u(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine hfmm3d_s_d_p_vec(nd,eps,zk,nsource,source,dipvec,pot,ier)

This subroutine evaluates the potential

\[u_{\ell}(x) = -\sum_{j=1}^{N} v_{\ell,j} \cdot \nabla \left( \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|}\right)\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • dipvec: double complex(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pot: double complex(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

hfmm3d_s_d_g

  • Evaluation points: Sources

  • Interaction kernel: Dipoles

  • Outputs requested: Potential and Gradient


subroutine hfmm3d_s_d_g(eps,zk,nsource,source,dipvec,pot,grad,ier)

This subroutine evaluates the potential and its gradient

\[u(x) = -\sum_{j=1}^{N} v_{j} \cdot \nabla \left( \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|}\right)\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • zk: double complex

    Helmholtz parameter, k

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • dipvec: double complex(3,nsource)

    Dipole strengths, \(v_{j}\)

Output arguments:

  • pot: double complex(nsource)

    Potential at source locations, \(u(x_{j})\)

  • grad: double complex(3,nsource)

    Gradient at source locations, \(\nabla u(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine hfmm3d_s_d_g_vec(nd,eps,zk,nsource,source,dipvec,pot,grad,ier)

This subroutine evaluates the potential and its gradient

\[u_{\ell}(x) = -\sum_{j=1}^{N} v_{\ell,j} \cdot \nabla \left( \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|}\right)\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • dipvec: double complex(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pot: double complex(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • grad: double complex(nd,3,nsource)

    Gradient at source locations, \(\nabla u_{\ell}(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

hfmm3d_s_cd_p

  • Evaluation points: Sources

  • Interaction kernel: Charges and Dipoles

  • Outputs requested: Potential


subroutine hfmm3d_s_cd_p(eps,zk,nsource,source,charge,dipvec,pot,ier)

This subroutine evaluates the potential

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{e^{ik\|x- x_{j}\|}}{\|x-x_{j}\|} - v_{j} \cdot \nabla \left( \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|}\right)\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • zk: double complex

    Helmholtz parameter, k

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double complex(nsource)

    Charge strengths, \(c_{j}\)

  • dipvec: double complex(3,nsource)

    Dipole strengths, \(v_{j}\)

Output arguments:

  • pot: double complex(nsource)

    Potential at source locations, \(u(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine hfmm3d_s_cd_p_vec(nd,eps,zk,nsource,source,charge,dipvec,pot,ier)

This subroutine evaluates the potential

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{e^{ik\|x- x_{j}\|}}{\|x-x_{j}\|} - v_{\ell,j} \cdot \nabla \left( \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|}\right)\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double complex(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

  • dipvec: double complex(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pot: double complex(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

hfmm3d_s_cd_g

  • Evaluation points: Sources

  • Interaction kernel: Charges and Dipoles

  • Outputs requested: Potential and Gradient


subroutine hfmm3d_s_cd_g(eps,zk,nsource,source,charge,dipvec,pot,grad,ier)

This subroutine evaluates the potential and its gradient

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{e^{ik\|x- x_{j}\|}}{\|x-x_{j}\|} - v_{j} \cdot \nabla \left( \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|}\right)\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • zk: double complex

    Helmholtz parameter, k

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double complex(nsource)

    Charge strengths, \(c_{j}\)

  • dipvec: double complex(3,nsource)

    Dipole strengths, \(v_{j}\)

Output arguments:

  • pot: double complex(nsource)

    Potential at source locations, \(u(x_{j})\)

  • grad: double complex(3,nsource)

    Gradient at source locations, \(\nabla u(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine hfmm3d_s_cd_g_vec(nd,eps,zk,nsource,source,charge,dipvec,pot,grad,ier)

This subroutine evaluates the potential and its gradient

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{e^{ik\|x- x_{j}\|}}{\|x-x_{j}\|} - v_{\ell,j} \cdot \nabla \left( \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|}\right)\]

at the source locations \(x=x_{j}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double complex(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

  • dipvec: double complex(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pot: double complex(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • grad: double complex(nd,3,nsource)

    Gradient at source locations, \(\nabla u_{\ell}(x_{j})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

hfmm3d_t_c_p

  • Evaluation points: Targets

  • Interaction kernel: Charges

  • Outputs requested: Potential


subroutine hfmm3d_t_c_p(eps,zk,nsource,source,charge,ntarg,targ,pottarg,ier)

This subroutine evaluates the potential

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{e^{ik\|x- x_{j}\|}}{\|x-x_{j}\|}\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • zk: double complex

    Helmholtz parameter, k

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double complex(nsource)

    Charge strengths, \(c_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pottarg: double complex(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine hfmm3d_t_c_p_vec(nd,eps,zk,nsource,source,charge,ntarg,targ,pottarg,ier)

This subroutine evaluates the potential

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{e^{ik\|x- x_{j}\|}}{\|x-x_{j}\|}\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double complex(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

Output arguments:

  • pottarg: double complex(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

hfmm3d_t_c_g

  • Evaluation points: Targets

  • Interaction kernel: Charges

  • Outputs requested: Potential and Gradient


subroutine hfmm3d_t_c_g(eps,zk,nsource,source,charge,ntarg,targ,pottarg,gradtarg,ier)

This subroutine evaluates the potential and its gradient

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{e^{ik\|x- x_{j}\|}}{\|x-x_{j}\|}\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • zk: double complex

    Helmholtz parameter, k

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double complex(nsource)

    Charge strengths, \(c_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pottarg: double complex(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • gradtarg: double complex(3,ntarg)

    Gradient at target locations, \(\nabla u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine hfmm3d_t_c_g_vec(nd,eps,zk,nsource,source,charge,ntarg,targ,pottarg,gradtarg,ier)

This subroutine evaluates the potential and its gradient

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{e^{ik\|x- x_{j}\|}}{\|x-x_{j}\|}\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double complex(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

Output arguments:

  • pottarg: double complex(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • gradtarg: double complex(nd,3,ntarg)

    Gradient at target locations, \(\nabla u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

hfmm3d_t_d_p

  • Evaluation points: Targets

  • Interaction kernel: Dipoles

  • Outputs requested: Potential


subroutine hfmm3d_t_d_p(eps,zk,nsource,source,dipvec,ntarg,targ,pottarg,ier)

This subroutine evaluates the potential

\[u(x) = -\sum_{j=1}^{N} v_{j} \cdot \nabla \left( \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|}\right)\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • zk: double complex

    Helmholtz parameter, k

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • dipvec: double complex(3,nsource)

    Dipole strengths, \(v_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pottarg: double complex(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine hfmm3d_t_d_p_vec(nd,eps,zk,nsource,source,dipvec,ntarg,targ,pottarg,ier)

This subroutine evaluates the potential

\[u_{\ell}(x) = -\sum_{j=1}^{N} v_{\ell,j} \cdot \nabla \left( \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|}\right)\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • dipvec: double complex(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pottarg: double complex(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

hfmm3d_t_d_g

  • Evaluation points: Targets

  • Interaction kernel: Dipoles

  • Outputs requested: Potential and Gradient


subroutine hfmm3d_t_d_g(eps,zk,nsource,source,dipvec,ntarg,targ,pottarg,gradtarg,ier)

This subroutine evaluates the potential and its gradient

\[u(x) = -\sum_{j=1}^{N} v_{j} \cdot \nabla \left( \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|}\right)\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • zk: double complex

    Helmholtz parameter, k

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • dipvec: double complex(3,nsource)

    Dipole strengths, \(v_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pottarg: double complex(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • gradtarg: double complex(3,ntarg)

    Gradient at target locations, \(\nabla u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine hfmm3d_t_d_g_vec(nd,eps,zk,nsource,source,dipvec,ntarg,targ,pottarg,gradtarg,ier)

This subroutine evaluates the potential and its gradient

\[u_{\ell}(x) = -\sum_{j=1}^{N} v_{\ell,j} \cdot \nabla \left( \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|}\right)\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • dipvec: double complex(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pottarg: double complex(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • gradtarg: double complex(nd,3,ntarg)

    Gradient at target locations, \(\nabla u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

hfmm3d_t_cd_p

  • Evaluation points: Targets

  • Interaction kernel: Charges and Dipoles

  • Outputs requested: Potential


subroutine hfmm3d_t_cd_p(eps,zk,nsource,source,charge,dipvec,ntarg,targ,pottarg,ier)

This subroutine evaluates the potential

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{e^{ik\|x- x_{j}\|}}{\|x-x_{j}\|} - v_{j} \cdot \nabla \left( \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|}\right)\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • zk: double complex

    Helmholtz parameter, k

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double complex(nsource)

    Charge strengths, \(c_{j}\)

  • dipvec: double complex(3,nsource)

    Dipole strengths, \(v_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pottarg: double complex(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine hfmm3d_t_cd_p_vec(nd,eps,zk,nsource,source,charge,dipvec,ntarg,targ,pottarg,ier)

This subroutine evaluates the potential

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{e^{ik\|x- x_{j}\|}}{\|x-x_{j}\|} - v_{\ell,j} \cdot \nabla \left( \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|}\right)\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double complex(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

  • dipvec: double complex(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pottarg: double complex(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

hfmm3d_t_cd_g

  • Evaluation points: Targets

  • Interaction kernel: Charges and Dipoles

  • Outputs requested: Potential and Gradient


subroutine hfmm3d_t_cd_g(eps,zk,nsource,source,charge,dipvec,ntarg,targ,pottarg,gradtarg,ier)

This subroutine evaluates the potential and its gradient

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{e^{ik\|x- x_{j}\|}}{\|x-x_{j}\|} - v_{j} \cdot \nabla \left( \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|}\right)\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • zk: double complex

    Helmholtz parameter, k

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double complex(nsource)

    Charge strengths, \(c_{j}\)

  • dipvec: double complex(3,nsource)

    Dipole strengths, \(v_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pottarg: double complex(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • gradtarg: double complex(3,ntarg)

    Gradient at target locations, \(\nabla u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine hfmm3d_t_cd_g_vec(nd,eps,zk,nsource,source,charge,dipvec,ntarg,targ,pottarg,gradtarg,ier)

This subroutine evaluates the potential and its gradient

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{e^{ik\|x- x_{j}\|}}{\|x-x_{j}\|} - v_{\ell,j} \cdot \nabla \left( \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|}\right)\]

at the target locations \(x=t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double complex(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

  • dipvec: double complex(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pottarg: double complex(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • gradtarg: double complex(nd,3,ntarg)

    Gradient at target locations, \(\nabla u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

hfmm3d_st_c_p

  • Evaluation points: Sources and Targets

  • Interaction kernel: Charges

  • Outputs requested: Potential


subroutine hfmm3d_st_c_p(eps,zk,nsource,source,charge,pot,ntarg,targ,pottarg,ier)

This subroutine evaluates the potential

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{e^{ik\|x- x_{j}\|}}{\|x-x_{j}\|}\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • zk: double complex

    Helmholtz parameter, k

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double complex(nsource)

    Charge strengths, \(c_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pot: double complex(nsource)

    Potential at source locations, \(u(x_{j})\)

  • pottarg: double complex(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine hfmm3d_st_c_p_vec(nd,eps,zk,nsource,source,charge,pot,ntarg,targ,pottarg,ier)

This subroutine evaluates the potential

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{e^{ik\|x- x_{j}\|}}{\|x-x_{j}\|}\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double complex(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

Output arguments:

  • pot: double complex(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • pottarg: double complex(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

hfmm3d_st_c_g

  • Evaluation points: Sources and Targets

  • Interaction kernel: Charges

  • Outputs requested: Potential and Gradient


subroutine hfmm3d_st_c_g(eps,zk,nsource,source,charge,pot,grad,ntarg,targ,pottarg,gradtarg,ier)

This subroutine evaluates the potential and its gradient

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{e^{ik\|x- x_{j}\|}}{\|x-x_{j}\|}\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • zk: double complex

    Helmholtz parameter, k

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double complex(nsource)

    Charge strengths, \(c_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pot: double complex(nsource)

    Potential at source locations, \(u(x_{j})\)

  • grad: double complex(3,nsource)

    Gradient at source locations, \(\nabla u(x_{j})\)

  • pottarg: double complex(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • gradtarg: double complex(3,ntarg)

    Gradient at target locations, \(\nabla u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine hfmm3d_st_c_g_vec(nd,eps,zk,nsource,source,charge,pot,grad,ntarg,targ,pottarg,gradtarg,ier)

This subroutine evaluates the potential and its gradient

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{e^{ik\|x- x_{j}\|}}{\|x-x_{j}\|}\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double complex(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

Output arguments:

  • pot: double complex(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • grad: double complex(nd,3,nsource)

    Gradient at source locations, \(\nabla u_{\ell}(x_{j})\)

  • pottarg: double complex(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • gradtarg: double complex(nd,3,ntarg)

    Gradient at target locations, \(\nabla u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

hfmm3d_st_d_p

  • Evaluation points: Sources and Targets

  • Interaction kernel: Dipoles

  • Outputs requested: Potential


subroutine hfmm3d_st_d_p(eps,zk,nsource,source,dipvec,pot,ntarg,targ,pottarg,ier)

This subroutine evaluates the potential

\[u(x) = -\sum_{j=1}^{N} v_{j} \cdot \nabla \left( \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|}\right)\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • zk: double complex

    Helmholtz parameter, k

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • dipvec: double complex(3,nsource)

    Dipole strengths, \(v_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pot: double complex(nsource)

    Potential at source locations, \(u(x_{j})\)

  • pottarg: double complex(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine hfmm3d_st_d_p_vec(nd,eps,zk,nsource,source,dipvec,pot,ntarg,targ,pottarg,ier)

This subroutine evaluates the potential

\[u_{\ell}(x) = -\sum_{j=1}^{N} v_{\ell,j} \cdot \nabla \left( \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|}\right)\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • dipvec: double complex(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pot: double complex(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • pottarg: double complex(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

hfmm3d_st_d_g

  • Evaluation points: Sources and Targets

  • Interaction kernel: Dipoles

  • Outputs requested: Potential and Gradient


subroutine hfmm3d_st_d_g(eps,zk,nsource,source,dipvec,pot,grad,ntarg,targ,pottarg,gradtarg,ier)

This subroutine evaluates the potential and its gradient

\[u(x) = -\sum_{j=1}^{N} v_{j} \cdot \nabla \left( \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|}\right)\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • zk: double complex

    Helmholtz parameter, k

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • dipvec: double complex(3,nsource)

    Dipole strengths, \(v_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pot: double complex(nsource)

    Potential at source locations, \(u(x_{j})\)

  • grad: double complex(3,nsource)

    Gradient at source locations, \(\nabla u(x_{j})\)

  • pottarg: double complex(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • gradtarg: double complex(3,ntarg)

    Gradient at target locations, \(\nabla u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine hfmm3d_st_d_g_vec(nd,eps,zk,nsource,source,dipvec,pot,grad,ntarg,targ,pottarg,gradtarg,ier)

This subroutine evaluates the potential and its gradient

\[u_{\ell}(x) = -\sum_{j=1}^{N} v_{\ell,j} \cdot \nabla \left( \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|}\right)\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • dipvec: double complex(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pot: double complex(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • grad: double complex(nd,3,nsource)

    Gradient at source locations, \(\nabla u_{\ell}(x_{j})\)

  • pottarg: double complex(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • gradtarg: double complex(nd,3,ntarg)

    Gradient at target locations, \(\nabla u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

hfmm3d_st_cd_p

  • Evaluation points: Sources and Targets

  • Interaction kernel: Charges and Dipoles

  • Outputs requested: Potential


subroutine hfmm3d_st_cd_p(eps,zk,nsource,source,charge,dipvec,pot,ntarg,targ,pottarg,ier)

This subroutine evaluates the potential

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{e^{ik\|x- x_{j}\|}}{\|x-x_{j}\|} - v_{j} \cdot \nabla \left( \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|}\right)\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • zk: double complex

    Helmholtz parameter, k

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double complex(nsource)

    Charge strengths, \(c_{j}\)

  • dipvec: double complex(3,nsource)

    Dipole strengths, \(v_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pot: double complex(nsource)

    Potential at source locations, \(u(x_{j})\)

  • pottarg: double complex(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine hfmm3d_st_cd_p_vec(nd,eps,zk,nsource,source,charge,dipvec,pot,ntarg,targ,pottarg,ier)

This subroutine evaluates the potential

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{e^{ik\|x- x_{j}\|}}{\|x-x_{j}\|} - v_{\ell,j} \cdot \nabla \left( \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|}\right)\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double complex(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

  • dipvec: double complex(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pot: double complex(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • pottarg: double complex(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

hfmm3d_st_cd_g

  • Evaluation points: Sources and Targets

  • Interaction kernel: Charges and Dipoles

  • Outputs requested: Potential and Gradient


subroutine hfmm3d_st_cd_g(eps,zk,nsource,source,charge,dipvec,pot,grad,ntarg,targ,pottarg,gradtarg,ier)

This subroutine evaluates the potential and its gradient

\[u(x) = \sum_{j=1}^{N} c_{j} \frac{e^{ik\|x- x_{j}\|}}{\|x-x_{j}\|} - v_{j} \cdot \nabla \left( \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|}\right)\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • eps: double precision

    precision requested

  • zk: double complex

    Helmholtz parameter, k

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • charge: double complex(nsource)

    Charge strengths, \(c_{j}\)

  • dipvec: double complex(3,nsource)

    Dipole strengths, \(v_{j}\)

  • ntarg: integer

    Number of targets

  • targ: double precision(3,ntarg)

    Target locations, \(t_{i}\)

Output arguments:

  • pot: double complex(nsource)

    Potential at source locations, \(u(x_{j})\)

  • grad: double complex(3,nsource)

    Gradient at source locations, \(\nabla u(x_{j})\)

  • pottarg: double complex(ntarg)

    Potential at target locations, \(u(t_{i})\)

  • gradtarg: double complex(3,ntarg)

    Gradient at target locations, \(\nabla u(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory


Vectorized version:

subroutine hfmm3d_st_cd_g_vec(nd,eps,zk,nsource,source,charge,dipvec,pot,grad,ntarg,targ,pottarg,gradtarg,ier)

This subroutine evaluates the potential and its gradient

\[u_{\ell}(x) = \sum_{j=1}^{N} c_{\ell,j} \frac{e^{ik\|x- x_{j}\|}}{\|x-x_{j}\|} - v_{\ell,j} \cdot \nabla \left( \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|}\right)\]

at the source and target locations \(x=x_{j},t_{i}\). When \(x=x_{j}\), the term corresponding to \(x_{j}\) is dropped from the sum.

Input arguments:

  • nd: integer

    number of densities

  • charge: double complex(nd,nsource)

    Charge strengths, \(c_{\ell,j}\)

  • dipvec: double complex(nd,3,nsource)

    Dipole strengths, \(v_{\ell,j}\)

Output arguments:

  • pot: double complex(nd,nsource)

    Potential at source locations, \(u_{\ell}(x_{j})\)

  • grad: double complex(nd,3,nsource)

    Gradient at source locations, \(\nabla u_{\ell}(x_{j})\)

  • pottarg: double complex(nd,ntarg)

    Potential at target locations, \(u_{\ell}(t_{i})\)

  • gradtarg: double complex(nd,3,ntarg)

    Gradient at target locations, \(\nabla u_{\ell}(t_{i})\)

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

Stokes FMM

Let \(\mathcal{G}^{\textrm{stok}}(x,y)\) denote the Stokeslet given by

\[\begin{split}\mathcal{G}^{\textrm{stok}}(x,y)=\frac{1}{2 \|x-y\|^3} \begin{bmatrix} (x_{1}-y_{1})^2 + \|x-y \|^2 & (x_{1}-y_{1})(x_{2}-y_{2}) & (x_{1}-y_{1})(x_{3}-y_{3}) \\ (x_{2}-y_{2})(x_{1}-y_{1}) & (x_{2}-y_{2})^2 + \|x-y \|^2 & (x_{2}-y_{2})(x_{3}-y_{3}) \\ (x_{3}-y_{3})(x_{1}-y_{1}) & (x_{3}-y_{3})(x_{2}-y_{2}) & (x_{3}-y_{3})^2 + \|x-y \|^2 \end{bmatrix} \, ,\end{split}\]

and \(\mathcal{T}^{\textrm{stok}}(x,y)\) denote the Stresslet whose action on a vector \(v\) is given by

\[\begin{split}v\cdot \mathcal{T}^{\textrm{stok}}(x,y) = \frac{3 v \cdot (x-y)}{\|x-y \|^5} \begin{bmatrix} (x_{1}-y_{1})^2 & (x_{1}-y_{1})(x_{2}-y_{2}) & (x_{1}-y_{1})(x_{3}-y_{3}) \\ (x_{2}-y_{2})(x_{1}-y_{1}) & (x_{2}-y_{2})^2 & (x_{2}-y_{2})(x_{3}-y_{3}) \\ (x_{3}-y_{3})(x_{1}-y_{1}) & (x_{3}-y_{3})(x_{2}-y_{2}) & (x_{3}-y_{3})^2 \end{bmatrix} \, .\end{split}\]

The Stokes FMM evaluates the following velocity, its gradient and the associated pressure

\[u(x) = \sum_{m=1}^{N} \mathcal{G}^{\textrm{stok}}(x,x_{j}) \sigma_{j} + \nu_{j} \cdot \mathcal{T}^{\textrm{stok}}(x,x_{j}) \cdot \mu_{j} \, .\]

Here \(x_{j}\) are the source locations, \(\sigma_{j}\) are the Stokeslet densities, \(\nu_{j}\) are the stresslet orientation vectors, \(\mu_{j}\) are the stresslet densities, and rhw xollwxrion of \(x\) at which the velocity and its gradient are evaluated are referred to as the evaluation points.

Unlike the Laplace and Helmholtz FMM, currently we have only the guru interface for the Stokes FMM (for both the single density and the vectorized density cases) with appropriate flags for including or excluding the stokeslet/stresslet term in the interaction, and flags for computing velocity/velocity and pressure/velocity, pressure, and gradients at the evaluation points.

subroutine stfmm3d(nd,eps,nsource,source,ifstoklet,stoklet,ifstrslet,strslet,strsvec,ifppreg,pot,pre,grad,ntarg,targ,ifppregtarg,pottarg,pretarg,gradtarg,ier)

Input arguments:

  • nd: integer

    Number of densities

  • eps: double precision

    Precision requested

  • nsource: integer

    Number of sources

  • source: double precision(3,nsource)

    Source locations, \(x_{j}\)

  • ifstoklet: integer

    Flag for including Stokeslet (\(\sigma_{j}\)) term in interaction kernel Stokeslet term will be included if ifstoklet = 1

  • stoklet: double precision(nd,3,nsource)

    Stokeslet strengths, \(\sigma_{j}\)

  • ifstrslet: integer

    Flag for including Stresslet (\(\mu_{j},\nu_{j}\)) term in interaction kernel Stresslet term will be included if ifstrslet = 1

  • strslet: double precision(nd,3,nsource)

    Stresslet strengths, \(\mu_{j}\)

  • strsvec: double precision(nd,3,nsource)

    Stresslet orientation vectors, \(\nu_{j}\)

  • ifppreg: integer
    Flag for computing velocity, pressure and/or gradients at source locations
    ifppreg = 1, compute velocity
    ifppreg = 2, compute velocity and pressure
    ifppreg = 3, compute veloicty, pressure and gradient
  • ntarg: integer

    Number of targets

  • targets: double precision (3,ntarg)

    Target locations \(x\)

  • ifppregtarg: integer
    Flag for computing velocity, pressure and/or gradients at target locations
    ifppregtarg = 1, compute velocity
    ifppregtarg = 2, compute velocity and pressure
    ifppregtarg = 3, compute veloicty, pressure and gradient

Output arguments:

  • pot: double precision (nd,3,nsource)

    Velocity at source locations if requested

  • pre: double precision (nd,nsource)

    Pressure at source locations if requested

  • grad: double precision (nd,3,3,nsource)

    Gradient at source locations if requested

  • pottarg: double precision (nd,3,ntarg)

    Velocity at target locations if requested

  • pretarg: double precision (nd,ntarg)

    Pressure at target locations if requested

  • gradtarg: double precision (nd,3,3,ntarg)

    Gradient at target locations if requested

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

Example drivers:

  • examples/stfmm3d_example.f. The corresponding makefile is examples/stfmm3d_example.make

Maxwell FMM

The Maxwell FMM evaluates the following field, its curl, and its divergence

\[E(x) = \sum_{j=1}^{N} \nabla \times \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|} M_{j} + \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|} J_{j} + \nabla \frac{e^{ik\|x-x_{j}\|}}{\|x-x_{j}\|} \rho_{j} \, .\]

Here \(x_{j}\) are the source locations, \(M_{j}\) are the magnetic current densities, \(J_{j}\) are the electric current densities, \(\rho_{j}\) are the electric charge densities, and the collection of \(x\) at which the field, its curl and its divergence are evalauated are referred to as the evaluation points.

Unlike the Laplace and Helmholtz FMM, currently we have only the guru interface for the Maxwell FMM (for both the single density and the vectorized density cases) with appropriate flags for including or excluding the magnetic current/electric current/electric charge term in the interaction, and flags for computing field/curl/divergence at the evaluation points.

subroutine emfmm3d(nd,eps,zk,ns,source,ifh_current,h_current,ife_current,e_current,ife_charge,e_charge,nt,targets,ifE,E,ifcurlE,curlE,ifdivE,divE,ier)

Input arguments:

  • nd: integer

    Number of densities

  • eps: double precision

    Precision requested

  • zk: double complex

    Wave number \(k\)

  • ns: integer

    Number of sources

  • source: double precision(3,ns)

    Source locations, \(x_{j}\)

  • ifh_current: integer

    Flag for including magnetic current (\(M_{j}\)) term in interaction kernel. Magnetic current term will be included if ifh_current = 1

  • h_current: double complex(nd,3,ns)

    Magnetic currents, \(M_{j}\)

  • ife_current: integer

    Flag for including electric current (\(J_{j}\)) term in interaction kernel. Electric current term will be included if ife_current = 1

  • e_current: double complex(nd,3,ns)

    Electric currents, \(J_{j}\)

  • ife_charge: integer

    Flag for including electric charge (\(\rho_{j}\)) term in interaction kernel. Electric charge term will be included if ife_charge = 1

  • e_charge: double complex(nd,ns)

    Electric charges, \(\rho_{j}\)

  • nt: integer

    Number of targets

  • targets: double precision (3,nt)

    Target locations \(x\)

  • ifE: integer

    Flag for computing field. The field \(E\) will be returned if ifE = 1

  • ifcurlE: integer

    Flag for computing curl of field. \(\nabla \times E\) will be returned if ifcurlE = 1

  • ifdivE: integer

    Flag for computing divergence of field. \(\nabla \cdot E\) will be returned if ifdivE = 1

Output arguments:

  • E: double complex (nd,3,nt)

    Field at the evaluation points if requested

  • curlE: double complex (nd,3,nt)

    Curl of field at the evaluation points if requested

  • divE: double complex (nd,nt)

    Divergence of field at the evaluation points if requested

  • ier: integer

    Error flag; ier=0 implies successful execution, and ier=4/8 implies insufficient memory

Example drivers:

  • examples/emfmm3d_example.f. The corresponding makefile is examples/emfmm3d_example.make

List of interfaces

C interfaces

All of the above fortran routines can be called from c by including the header utils.h and lfmm3d_c.h for Laplace FMMs or hfmm3d_c.h for Helmholtz FMMs.

For example, the subroutine to evaluate the potential and gradient, at a collection of targets \(t_{i}\) due to a collection of Helmholtz charges is:

hfmm3d_t_c_g_

In general, to call a fortran subroutine from c use:

"<fortran subroutine name>"_("<calling sequence>")

Note

All the variables in the calling sequence must be passed as pointers from c.

Note

For the vectorized subroutines, the charge strengths, dipole strengths, potentials, and gradients are interleaved as opposed to provided in a sequential manner. For example for three sets of charge strengths, they should be stored as \(c_{1,1}, c_{2,1}, c_{3,1}, c_{1,2}, c_{2,2},c_{3,2} \ldots c_{1,N}, c_{2,N}, c_{3,N}\).

Example drivers:

  • Laplace:

    • c/lfmm3d_example.c. The corresponding makefile is c/lfmm3d_example.make

    • c/lfmm3d_vec_example.c. The corresponding makefile is c/lfmm3d_vec_example.make

  • Helmholtz:

    • c/hfmm3d_example.c. The corresponding makefile is c/hfmm3d_example.make

    • c/hfmm3d_vec_example.c. The corresponding makefile is c/hfmm3d_vec_example.make

The Maxwell and Stokes interfaces are currently unavailable in C, and will be made available soon.