module FMM3D

Wrappers for the Flatiron Institute's FMM3D library.

wrappers

All N-body codes return output in an FMMVals structure. See documentation of N-body codes for details.

N-body interactions with the Helmholtz kernel

• hfmm3d: $O(N)$ fast mutlipole code
• h3ddir: $O(N^2)$ direct code

N-body interactions with the Laplace kernel

• lfmm3d: $O(N)$ fast mutlipole code
• l3ddir: $O(N^2)$ direct code

N-body interactions with Stokes kernels

• stfmm3d: $O(N)$ fast mutlipole code
• st3ddir: $O(N^2)$ direct code

N-body interactions with an electromagnetic kernel

• emfmm3d: $O(N)$ fast mutlipole code
• em3ddir: $O(N^2)$ direct code

lower level routines

For a list of lower level routines see the documentation for lower_level_routs

FMMVals()

Return type for FMM and direct computation function calls. Fields are nothing on return unless requested. See individual FMM/direct computation function documentation for specifics.

function fj, fjder = besseljs3d(nterms,z;scale=1.0,ifder=0)

This subroutine evaluates the first nterms spherical Bessel functions, and if requested, their derivatives. It incorporates a scaling parameter scale so that

  	fjs_n(z)=j_n(z)/SCALE^n
  	fjder_n(z)=\frac{\partial fjs_n(z)}{\partial z}

Input

• nterms::Integer order of expansion of output array fjs
• z::ComplexF64 argument of the spherical Bessel functions
• scale::Float64 scaling factor
• ifder::Integer1 flag indicating whether to calculate fjder 0 NO 1 YES

OUTPUT:

• fjs::Array{ComplexF64} array of length nterms+1 of scaled Bessel functions.
• fjder::Array{ComplexF64} array of derivatives of scaled Bessel functions, if requested.

    vals = em3ddir(zk,sources,targets;h_current=nothing,e_current=nothing,e_charge=nothing,
                ifEtarg=false,ifdivEtarg=false,ifcurlEtarg=false,
                nd=1,thresh=1e-16)

This function computes N-body electromagnetic interactions in three dimensions where the interaction kernels are the Helmholtz kernel and its curl and gradient (see below). This is the $O(N^2)$ direct code. By convention this code only computes the effect of sources on targets. If the value at sources is also needed, the routine can be called again with targets equal to the source locations.

The Helmholtz Green's function (without the $1/4 \pi$ scaling) is

\[ G_k(x,y) = \frac{e^{ik \|x-y\|}}{\|x-y\|} \, , \]

This routine computes the sum

\[E(x) = \sum_{j=1}^m \nabla \times (G_k(x,x_{j}) h_current_j) + G_k(x,x_{j}) e_current_j + \nabla G_k(x,x_{j}) \e_charge_j\]

where $h_current_{j}$ and $e_current_{j}$ are vector densities $\e_charge_{j}$ is a scalar density, and $x_{j}$ are the source locations. When $x=x_{m}$, the term corresponding to $x_{m}$ is dropped from the sum

Input

• zk::ComplexF64 Helmholtz parameter
• sources::Array{Float64} size (3,n) source locations ($x_{m}$)
• targets::Array{Float64} size (3,nt) target locations ($x$)

Optional Keyword Input

• h_current::Array{ComplexF64} size (nd,3,n) or (3,n) vector densities
• e_current::Array{ComplexF64} size (nd,3,n) or (3,n) vector densities
• e_charge::Array{ComplexF64} size (nd,n) or (n) scalar densities
• ifEtarg::Bool E field evaluation flag, if true evaluate E at target points.
• ifdivEtarg::Bool divergence of E field evaluation flag, if true evaluate divergence of E at target points.
• ifcurlEtarg::Bool curl of E field evaluation flag, if true evaluate curl of E at target points.
• nd::Integer number of densities of each type

Note: if all default values are used for optional input, nothing is computed.

Output

vals::FMMVals with the fields

• vals.Etarg::Array{ComplexF64} size (nd,3,nt) or (3,nt) E field at target locations if requested
• vals.divEtarg::Array{ComplexF64} size (nd,nt) or (nt) divergence of E field at target locations if requested
• vals.curlEtarg::Array{ComplexF64} size (nd,3,nt) or (3,nt) curl of E field at target locations if requested

    vals = emfmm3d(eps,zk,sources;h_current=nothing,e_current=nothing,e_charge=nothing,
                ifE=false,ifdivE=false,ifcurlE=false,
                ifEtarg=false,ifdivEtarg=false,ifcurlEtarg=false,
                nd=1,targets=nothing)

This function computes N-body electromagnetic interactions in three dimensions where the interaction kernels are the Helmholtz kernel and its curl and gradient (see below). This is the $O(N)$ fast multipole code which computes the interactions to the requested precision.

The Helmholtz Green's function (without the $1/4 \pi$ scaling) is

\[ G_k(x,y) = \frac{e^{ik \|x-y\|}}{\|x-y\|} \, , \]

This routine computes the sum

\[E(x) = \sum_{j=1}^m \nabla \times (G_k(x,x_{j}) h_current_j) + G_k(x,x_{j}) e_current_j + \nabla G_k(x,x_{j}) \e_charge_j\]

where $h_current_{j}$ and $e_current_{j}$ are vector densities $\e_charge_{j}$ is a scalar density, and $x_{j}$ are the source locations. When $x=x_{m}$, the term corresponding to $x_{m}$ is dropped from the sum

Input

• eps::Float64 precision requested
• zk::ComplexF64 Helmholtz parameter
• sources::Array{Float64} size (3,n) source locations ($x_{m}$)

Optional Keyword Input

• targets::Array{Float64} size (3,nt) target locations ($x$)
• h_current::Array{ComplexF64} size (nd,3,n) or (3,n) vector densities
• e_current::Array{ComplexF64} size (nd,3,n) or (3,n) vector densities
• e_charge::Array{ComplexF64} size (nd,n) or (n) scalar densities
• ifE::Bool E field evaluation flag, if true evaluate E at source points.
• ifdivE::Bool divergence of E field evaluation flag, if true evaluate divergence of E at source points.
• ifcurlE::Bool curl of E field evaluation flag, if true evaluate curl of E at source points.
• ifEtarg::Bool E field evaluation flag, if true evaluate E at target points.
• ifdivEtarg::Bool divergence of E field evaluation flag, if true evaluate divergence of E at target points.
• ifcurlEtarg::Bool curl of E field evaluation flag, if true evaluate curl of E at target points.
• nd::Integer number of densities of each type

Note: if all default values are used for optional input, nothing is computed.

Output

vals::FMMVals with the fields

• vals.E::Array{ComplexF64} size (nd,3,n) or (3,n) E field at source locations if requested
• vals.divE::Array{ComplexF64} size (nd,n) or (n) divergence of E field at source locations if requested
• vals.curlE::Array{ComplexF64} size (nd,3,n) or (3,n) curl of E field at source locations if requested
• vals.Etarg::Array{ComplexF64} size (nd,3,nt) or (3,nt) E field at target locations if requested
• vals.divEtarg::Array{ComplexF64} size (nd,nt) or (nt) divergence of E field at target locations if requested
• vals.curlEtarg::Array{ComplexF64} size (nd,3,nt) or (3,nt) curl of E field at target locations if requested
• vals.ier error flag as returned by FMM3D library. A value of 0 indicates a successful call.

Non-zero values may indicate insufficient memory available. See the documentation for the FMM3D library. If not set (nothing), then FMM3D library was never called.

    anyfail, n, nt, ifh_current, ife_current, ife_charge = (
        emfmm3dinputcheck(sources,h_current,e_current,e_charge,
                        targets,ifE,ifdivE,ifcurlE,
                        ifEtarg,ifdivEtarg,
                        ifcurlEtarg,nd) )

Input checking for electromagnetics routines. Checks the sizes of these arrays and the compatibility of the various flags, nd, and provided arrays. If something is off, a warning is issued and anyfail is set to true. Non-fatal mistakes result in a warning but anyfail remains false.

Output:

• anyfail - boolean is true if any checks fail, false otherwise
• n - integer, number of sources
• nt - integer, number of targets
• ifh_current - integer, is 1 if there is an h_current vector, 0 otherwise
• ife_current - integer, is 1 if there is a e_current vector, 0 otherwise
• ife_charge - integer, is 1 if there is a e_charge density, 0 otherwise

    vals = h3ddir(zk,sources,targets;charges=nothing,
                    dipvecs=nothing,pgt=0,nd=1,
                    thresh=1e-16)

This function computes the N-body Helmholtz interactions in three dimensions where the interaction kernel is given by $e^{ikr}/r$ and its gradients. This is the $O(N^2)$ direct evaluation code. By convention this code only computes the effect of sources on targets. If the value at sources is also needed, the routine can be called again with targets equal to the source locations.

\[ 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) \, , \]

where $c_{j}$ are the charge densities, $v_{j}$ are the dipole orientation vectors, and $x_{j}$ are the source locations. When $x=x_{m}$, the term corresponding to $x_{m}$ is dropped from the sum

Input

• zk::ComplexF64 Helmholtz parameter
• sources::Array{Float64} size (3,n) source locations ($x_{j}$)
• targets::Array{Float64} size (3,nt) target locations ($x$)

Optional Keyword Input

• charges::Array{ComplexF64} size (nd,n) or (n) charge densities (c_{j})
• dipvecs::Array{ComplexF64} size (nd,3,n) or (3,n)) dipole orientation vectors ($v_{j}$)
• pgt::Integer target eval flag.
  • Potential at targets evaluated if pgt == 1.
  • Potential and gradient at targets evaluated if pgt == 2
  • Potential, gradient, and Hessian at targets evaluated if pgt == 3
• nd::Integer number of densities
• thresh::Float64 threshold for ignoring interactions when $\|x-x_{j}\| \leq thresh$

Note: if all default values are used for optional input, nothing is computed.

Output

vals::FMMVals with the fields

Note that the Hessian is returned as a length 6 vector at each point with the second derivatives in the order: $\partial_{xx}$, $\partial_{yy}$, $\partial_{zz}$, $\partial_{xy}$, $\partial_{xz}$, $\partial_{yz}$.

• vals.pottarg::Array{ComplexF64} size (nd,nt) or (nt) potential at target locations if requested
• vals.gradtarg::Array{ComplexF64} size (nd,3,nt) or (3,nt) gradient at target locations if requested
• vals.hesstarg::Array{ComplexF64} size (nd,6,nt) or (6,nt) gradient at target locations if requested

    vals = hfmm3d(eps,zk,sources;charges=nothing,dipvecs=nothing,
                  targets=nothing,pg=0,pgt=0,nd=1)

This function computes the N-body Helmholtz interactions in three dimensions where the interaction kernel is given by $e^{ikr}/r$ and its gradients. This is the $O(N)$ fast multipole code which computes the interactions to the requested precision.

\[ 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) \, , \]

where $c_{j}$ are the charge densities, $v_{j}$ are the dipole orientation vectors, and $x_{j}$ are the source locations. When $x=x_{m}$, the term corresponding to $x_{m}$ is dropped from the sum

Input

• eps::Float64 precision requested
• zk::ComplexF64 Helmholtz parameter
• sources::Array{Float64} size (3,n) source locations ($x_{j}$)

Optional Keyword Input

• charges::Array{ComplexF64} size (nd,n) or (n) charge densities (c_{j})
• dipvecs::Array{ComplexF64} size (nd,3,n) or (3,n)) dipole orientation vectors ($v_{j}$)
• targets::Array{Float64} size (3,nt) target locations ($x$)
• pg::Integer source eval flag.
  • Do not compute any quantities at sources if pg == 0
  • Potential ($u$) at sources evaluated if pg == 1.
  • Potential and gradient ($\nabla u$) at sources evaluated if pg == 2
  • Potential, gradient, and Hessian ($\nabla \nabla u$) at sources evaluated if pg == 3
• pgt::Integer target eval flag.
  • Do not compute any quantities at targets if pgt == 0
  • Potential at targets evaluated if pgt == 1.
  • Potential and gradient at targets evaluated if pgt == 2
  • Potential, gradient, and Hessian at targets evaluated if pgt == 3
• nd::Integer number of densities

Note: if all default values are used for optional input, nothing is computed.

Output

vals::FMMVals with the fields

Note that the Hessian is returned as a length 6 vector at each point with the second derivatives in the order: $\partial_{xx}$, $\partial_{yy}$, $\partial_{zz}$, $\partial_{xy}$, $\partial_{xz}$, $\partial_{yz}$.

• vals.pot::Array{ComplexF64} size (nd,n) or (n) potential at source locations if requested
• vals.grad::Array{ComplexF64} size (nd,3,n) or (3,n) gradient at source locations if requested
• vals.hess::Array{ComplexF64} size (nd,6,n) or (6,n) Hessian at source locations if requested
• vals.pottarg::Array{ComplexF64} size (nd,nt) or (nt) potential at target locations if requested
• vals.gradtarg::Array{ComplexF64} size (nd,3,nt) or (3,nt) gradient at target locations if requested
• vals.hesstarg::Array{ComplexF64} size (nd,6,nt) or (6,nt) Hessian at target locations if requested
• vals.ier error flag as returned by FMM3D library. A value of 0 indicates a successful call.

Non-zero values may indicate insufficient memory available. See the documentation for the FMM3D library. If not set (nothing), then FMM3D library was never called.

    vals = l3ddir(sources,targets;charges=nothing,
                dipvecs=nothing,pgt=0,nd=1,
                thresh=1e-16)

This function computes the N-body Laplace interactions in three dimensions where the interaction kernel is given by $1/r$ and its gradients. This is the $O(N^2)$ direct evaluation code. By convention this code only computes the effect of sources on targets. If the value at sources is also needed, the routine can be called again with targets equal to the source locations.

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

where $c_{j}$ are the charge densities, $v_{j}$ are the dipole orientation vectors, and $x_{j}$ are the source locations. When $x=x_{m}$, the term corresponding to $x_{m}$ is dropped from the sum

Input

• sources::Array{Float64} size (3,n) source locations ($x_{j}$)
• targets::Array{Float64} size (3,nt) target locations ($x$)

Optional Keyword Input

• charges::Array{ComplexF64} size (nd,n) or (n) charge densities (c_{j})
• dipvecs::Array{ComplexF64} size (nd,3,n) or (3,n)) dipole orientation vectors ($v_{j}$)
• pgt::Integer target eval flag.
  • Potential at targets evaluated if pgt == 1.
  • Potential and gradient at targets evaluated if pgt == 2
• nd::Integer number of densities
• thresh::Float64 threshold for ignoring interactions when $\|x-x_{j}\| \leq thresh$

Note: if all default values are used for optional input, nothing is computed.

Output

vals::FMMVals with the fields

• vals.pottarg::Array{ComplexF64} size (nd,nt) or (nt) potential at target locations if requested
• vals.gradtarg::Array{ComplexF64} size (nd,3,nt) or (3,nt) gradient at target locations if requested

    vals = lfmm3d(eps,sources;charges=nothing,dipvecs=nothing,
                  targets=nothing,pg=0,pgt=0,nd=1)

This function computes the N-body Laplace interactions in three dimensions where the interaction kernel is given by $$1/r$$ and its gradients. This is the $O(N)$ fast multipole code which computes the interactions to the requested precision.

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

where $c_{j}$ are the charge densities, $v_{j}$ are the dipole orientation vectors, and $x_{j}$ are the source locations. When $x=x_{m}$, the term corresponding to $x_{m}$ is dropped from the sum

Input

• eps::Float64 precision requested
• sources::Array{Float64} size (3,n) source locations ($x_{j}$)

Optional Keyword Input

• charges::Array{Float64} size (nd,n) or (n) charge densities (c_{j})
• dipvecs::Array{Float64} size (nd,3,n) or (3,n)) dipole orientation vectors ($v_{j}$)
• targets::Array{Float64} size (3,nt) target locations ($x$)
• pg::Integer source eval flag.
  • Potential ($u$) at sources evaluated if pg == 1.
  • Potential and gradient ($\nabla u$) at sources evaluated if pg == 2
• pgt::Integer target eval flag.
  • Potential at targets evaluated if pgt == 1.
  • Potential and gradient at targets evaluated if pgt == 2
• nd::Integer number of densities

Note: if all default values are used for optional input, nothing is computed.

Output

vals::FMMVals with the fields

Note that the Hessian is returned as a length 6 vector at each point with the second derivatives in the order: $\partial_{xx}$, $\partial_{yy}$, $\partial_{zz}$, $\partial_{xy}$, $\partial_{xz}$, $\partial_{yz}$.

• vals.pot::Array{Float64} size (nd,n) potential at source locations if requested
• vals.grad::Array{Float64} size (nd,3,n) gradient at source locations if requested
• vals.hess::Array{Float64} size (nd,6,n) Hessian at source locations if requested
• vals.pottarg::Array{Float64} size (nd,nt) potential at target locations if requested
• vals.gradtarg::Array{Float64} size (nd,3,nt) gradient at target locations if requested
• vals.hesstarg::Array{Float64} size (nd,6,nt) Hessian at target locations if requested
• vals.ier error flag as returned by FMM3D library. A value of 0 indicates a successful call.

Non-zero values may indicate insufficient memory available. See the documentation for the FMM3D library. If not set (nothing), then FMM3D library was never called.

Available lower level routines

• besseljs3d spherical Bessel function eval

    anyfail, n, nt, ifcharge, ifdipole = (
        scalarfmm3dinputcheck(sources,charges,dipvecs,
                                targets,pg,pgt,nd) )

Input checking for Helmholtz and Laplace routines. Checks the sizes of these arrays and the compatibility of the various flags, nd, and provided arrays. If something is off, a warning is issued and anyfail is set to true. Non-fatal mistakes result in a warning but anyfail remains false.

Output:

• anyfail - boolean is true if any checks fail, false otherwise
• n - integer, number of sources
• nt - integer, number of targets
• ifcharge - integer, is 1 if there are charges, 0 otherwise
• ifdipole - integer, is 1 if there are dipoles, 0 otherwise

    vals = st3ddir(sources,targets;stoklet=nothing,strslet=nothing,
                   strsvec=nothing,ppregt=0,nd=1,thresh=1e-16)

This function computes the N-body Stokes interactions in three dimensions where the interaction kernels are the Stokeslet and stresslet (see below). This is the $O(N^2)$ direct evaluation code. By convention this code only computes the effect of sources on targets. If the value at sources is also needed, the routine can be called again with targets equal to the source locations.

We take the following conventions for the Stokes kernels:

For a source $y$ and target $x$, let $r_i = x_i-y_i$ and let $r = \sqrt{r_1^2 + r_2^2 + r_3^2}$

The Stokeslet, $G_{ij}$, and its associated pressure tensor, $P_j$, (without the $1/4\pi$ scaling) are

\[G_{ij}(x,y) = (r_i r_j)/(2r^3) + \delta_{ij}/(2r) \; , \quad P_j(x,y) = r_j/r^3\]

The (Type I) stresslet, $T_{ijk}$, and its associated pressure tensor, $\Pi_{jk}$, (without the $1/4\pi$ scaling) are

\[ T_{ijk}(x,y) = -3 r_i r_j r_k/ r^5 \; , \quad PI_{jk} = -2 \delta_{jk} + 6 r_j r_k/r^5 \]

The output of this routine gives the velocity

\[ u_i(x) = \sum_{m=1}^{N} G_{ij}(x,x_m) \sigma^{(m)}_j + T_{ijk}(x,x_m) \mu^{(m)}_j \nu^{(m)}_k \, , \]

and, optionally, the pressure

\[ p(x) = \sum_{m=1}^{N} P_{j}(x,x_m) \sigma^{(m)}_j + \Pi{jk}(x,x_m) \mu^{(m)}_j \nu^{(m)}_k \, , \]

where $\sigma^{(m)}$ is a Stokeslet density 3-vector, $\mu^{(m)}$ and $\nu^{(m)}$ are the stresslet density and orientation 3-vectors, and $x_{m}$ are the source locations. When $x=x_{m}$, the term corresponding to $x_{m}$ is dropped from the sum. The gradient of the velocity can also be computed by request.

Input

• sources::Array{Float64} size (3,n) source locations ($x_{m}$)
• targets::Array{Float64} size (3,nt) target locations ($x$)

Optional Keyword Input

• stoklet::Array{Float64} size (nd,3,n) or (3,n) Stokeslet densities ($\sigma$)
• strslet::Array{Float64} size (nd,3,n) or (3,n) Stresslet densities (``\mu'') if either strsvec or strslet is provided, the other must be as well
• strsvec::Array{Float64} size (nd,3,n) or (3,n) Stresslet orientations ($\nu$) if either strsvec or strslet is provided, the other must be as well
• ppregt::Integer source eval flag.
  • Velocity ($u$) at targets evaluated if ppregt == 1.
  • Velocity and pressure ($p$) at targets evaluated if ppregt == 2
  • Velocity, pressure, and velocity gradient ($\nabla u$)
  at targets evaluated if ppregt == 3
• nd::Integer number of densities of each type

Note: if all default values are used for optional input, nothing is computed.

Output

vals::FMMVals with the fields

• vals.pottarg::Array{Float64} size (nd,3,nt) or (3,nt) velocity at target locations if requested
• vals.gradtarg::Array{Float64} size (nd,3,3,nt) or (3,nt) gradient of velocity at target locations if requested.

gradtarg[l,i,j,k] is the i-th component of the gradient of the j-th component of the velocity for the l-th density at the k-th target location

• vals.pretarg::Array{Float64} size (nd,nt) or (nt) pressure at target locations if requested

    vals = stfmm3d(eps,sources;stoklet=nothing,strslet=nothing,
                   strsvec=nothing,targets=nothing,ppreg=0,
                   ppregt=0,nd=1)

This function computes the N-body Stokes interactions in three dimensions where the interaction kernels are the Stokeslet and stresslet (see below). This is the $O(N)$ fast multipole code which computes the interactions to the requested precision.

We take the following conventions for the Stokes kernels:

For a source $y$ and target $x$, let $r_i = x_i-y_i$ and let $r = \sqrt{r_1^2 + r_2^2 + r_3^2}$

The Stokeslet, $G_{ij}$, and its associated pressure tensor, $P_j$, (without the $1/4\pi$ scaling) are

\[G_{ij}(x,y) = (r_i r_j)/(2r^3) + \delta_{ij}/(2r) \; , \quad P_j(x,y) = r_j/r^3\]

The (Type I) stresslet, $T_{ijk}$, and its associated pressure tensor, $\Pi_{jk}$, (without the $1/4\pi$ scaling) are

\[ T_{ijk}(x,y) = -3 r_i r_j r_k/ r^5 \; , \quad PI_{jk} = -2 \delta_{jk} + 6 r_j r_k/r^5 \]

The output of this routine gives the velocity

\[ u_i(x) = \sum_{m=1}^{N} G_{ij}(x,x_m) \sigma^{(m)}_j + T_{ijk}(x,x_m) \mu^{(m)}_j \nu^{(m)}_k \, , \]

and, optionally, the pressure

\[ p(x) = \sum_{m=1}^{N} P_{j}(x,x_m) \sigma^{(m)}_j + \Pi{jk}(x,x_m) \mu^{(m)}_j \nu^{(m)}_k \, , \]

where $\sigma^{(m)}$ is a Stokeslet density 3-vector, $\mu^{(m)}$ and $\nu^{(m)}$ are the stresslet density and orientation 3-vectors, and $x_{m}$ are the source locations. When $x=x_{m}$, the term corresponding to $x_{m}$ is dropped from the sum. The gradient of the velocity can also be computed by request.

Input

• eps::Float64 precision requested
• sources::Array{Float64} size (3,n) source locations ($x_{m}$)

Optional Keyword Input

• stoklet::Array{Float64} size (nd,3,n) or (3,n) Stokeslet densities ($\sigma$)
• strslet::Array{Float64} size (nd,3,n) or (3,n) Stresslet densities (``\mu'') if either strsvec or strslet is provided, the other must be as well
• strsvec::Array{Float64} size (nd,3,n) or (3,n) Stresslet orientations ($\nu$) if either strsvec or strslet is provided, the other must be as well
• targets::Array{Float64} size (3,nt) target locations ($x$)
• ppreg::Integer source eval flag.
  • Velocity ($u$) at sources evaluated if ppreg == 1.
  • Velocity and pressure ($p$) at sources evaluated if ppreg == 2
  • Velocity, pressure, and velocity gradient ($\nabla u$)
  at sources evaluated if ppreg == 3
• ppregt::Integer source eval flag.
  • Velocity ($u$) at targets evaluated if ppregt == 1.
  • Velocity and pressure ($p$) at targets evaluated if ppregt == 2
  • Velocity, pressure, and velocity gradient ($\nabla u$)
  at targets evaluated if ppregt == 3
• nd::Integer number of densities of each type

Note: if all default values are used for optional input, nothing is computed.

Output

vals::FMMVals with the fields

• vals.pot::Array{Float64} size (nd,3,n) or (3,n) velocity at source locations if requested
• vals.grad::Array{Float64} size (nd,3,3,n) or (3,n) gradient of velocity at source locations if requested.

grad[l,i,j,k] is the i-th component of the gradient of the j-th component of the velocity for the l-th density at the k-th source location

• vals.pre::Array{Float64} size (nd,n) or (n) pressure at source locations if requested
• vals.pottarg::Array{Float64} size (nd,3,nt) or (3,nt) velocity at target locations if requested
• vals.gradtarg::Array{Float64} size (nd,3,3,nt) or (3,nt) gradient of velocity at target locations if requested.

gradtarg[l,i,j,k] is the i-th component of the gradient of the j-th component of the velocity for the l-th density at the k-th target location

• vals.pretarg::Array{Float64} size (nd,nt) or (nt) pressure at target locations if requested
• vals.ier error flag as returned by FMM3D library. A value of 0 indicates a successful call.

Non-zero values may indicate insufficient memory available. See the documentation for the FMM3D library. If not set (nothing), then FMM3D library was never called.

    anyfail, n, nt, ifstoklet, ifstrslet = (
        stfmm3dinputcheck(sources,stoklet,strslet,strsvec,
                                targets,ppreg,ppregt,nd) )

Input checking for Stokes routines. Checks the sizes of these arrays and the compatibility of the various flags, nd, and provided arrays. If something is off, a warning is issued and anyfail is set to true. Non-fatal mistakes result in a warning but anyfail remains false.

Output:

• anyfail - boolean is true if any checks fail, false otherwise
• n - integer, number of sources
• nt - integer, number of targets
• ifstoklet - integer, is 1 if there are Stokeslets, 0 otherwise
• ifstrslet - integer, is 1 if there are (type I) stresslets, 0 otherwise