Skip to content

Commit

Permalink
Merge pull request #63 from flatironinstitute/stokes-rotlet-doublet
Browse files Browse the repository at this point in the history 
Stokes update
  • Loading branch information
@lu1and10
lu1and10 authored Nov 21, 2024
2 parents 3766fef + 2cb4aef commit c0dbe29
Show file tree
Hide file tree
Showing 25 changed files with 2,097 additions and 649 deletions.
67 changes: 62 additions & 5 deletions docs/fortran-c.rst
View file
Original file line number Diff line number Diff line change
Expand Up @@ -330,7 +330,7 @@ denote the Stokeslet given by
(x_{3}-y_{3})^2 + \|x-y \|^2
\end{bmatrix} \, ,
and $\mathcal{T}^{\textrm{stok}}(x,y)$ denote the Stresslet whose action on
let $\mathcal{T}^{\textrm{stok}}(x,y)$ denote the Stresslet whose action on
a vector $v$ is given by

.. math::
Expand All @@ -343,19 +343,60 @@ a vector $v$ is given by
(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{bmatrix} \, ,
let $\mathcal{R}^{\textrm{stok}}(x,y)$ denote the Rotlet whose action on
a vector $v$ is given by

.. math::
v\cdot \mathcal{R}^{\textrm{stok}}(x,y) =
\frac{v \cdot (x-y)}{4\pi\|x-y \|^3}
\begin{bmatrix}
1 & 0 & 0 \\
0 & 1 & 0 \\
0 & 0 & 1
\end{bmatrix} \, ,
and $\mathcal{D}^{\textrm{stok}}(x,y)$ denote the symmetric part of Doublet whose action on
a vector $v$ is given by

.. math::
v\cdot \mathcal{D}^{\textrm{stok}}(x,y) =
\frac{3 v \cdot (x-y)}{4\pi\|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} - \\
\frac{1}{4\pi\|x-y \|^3}
\begin{bmatrix}
v_1(x_{1}-y_{1}) & v_2(x_{1}-y_{1}) & v_3(x_{1}-y_{1}) \\
v_2(x_{2}-y_{2}) & v_2(x_{2}-y_{2}) & v_3(x_{2}-y_{2}) \\
v_3(x_{3}-y_{3}) & v_3(x_{3}-y_{3}) & v_3(x_{3}-y_{3})
\end{bmatrix} \, .
The Stokes FMM evaluates the following velocity, its gradient
and the associated pressure

.. math::
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} \, .
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} +
\nu^{r}_{j} \cdot \mathcal{R}^{\textrm{stok}}(x,x_{j}) \cdot \mu^{r}_{j} - \mu^{r}_{j} \cdot \mathcal{R}^{\textrm{stok}}(x,x_{j}) \cdot \nu^{r}_{j} + \\
\nu^{d}_{j} \cdot \mathcal{R}^{\textrm{stok}}(x,x_{j}) \cdot \mu^{d}_{j} - \mu^{d}_{j} \cdot \mathcal{R}^{\textrm{stok}}(x,x_{j}) \cdot \nu^{d}_{j} +
\nu^{d}_{j} \cdot \mathcal{D}^{\textrm{stok}}(x,x_{j}) \cdot \mu^{d}_{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 the locations $x$
are the stresslet densities,
$\nu^{r}_{j}$ are the rotlet orientation vectors, $\mu^{r}_{j}$
are the rotlet densities,
$\nu^{d}_{j}$ are the doublet orientation vectors, $\mu^{d}_{j}$
are the doublet densities,
and the locations $x$
at which the velocity and its gradient are evaluated are referred to
as the evaluation points.

Expand All @@ -370,7 +411,7 @@ gradients at the evaluation points.

.. code::
subroutine stfmm3d(nd,eps,nsource,source,ifstoklet,stoklet,ifstrslet,strslet,strsvec,ifppreg,pot,pre,grad,ntarg,targ,ifppregtarg,pottarg,pretarg,gradtarg,ier)
subroutine stfmm3d(nd,eps,nsource,source,ifstoklet,stoklet,ifstrslet,strslet,strsvec,ifrotlet,rotlet,rotvec,ifdoublet,doublet,doubvec,ifppreg,pot,pre,grad,ntarg,targ,ifppregtarg,pottarg,pretarg,gradtarg,ier)
Input arguments:

Expand All @@ -396,6 +437,22 @@ Input arguments:
Stresslet strengths, $\mu_{j}$
- strsvec: double precision(nd,3,nsource)
Stresslet orientation vectors, $\nu_{j}$
- ifrotlet: integer
Flag for including Rotlet ($\mu^{r}_{j},\nu^{r}_{j}$) term in interaction kernel
Rotlet term will be included if ifrotlet
= 1
- rotlet: double precision(nd,3,nsource)
Rotlet strengths, $\mu^{r}_{j}$
- rotvec: double precision(nd,3,nsource)
Rotlet orientation vectors, $\nu^{r}_{j}$
- ifdoublet: integer
Flag for including Doublet ($\mu^{d}_{j},\nu^{d}_{j}$) term in interaction kernel
Doublet term will be included if ifdoublet
= 1
- doublet: double precision(nd,3,nsource)
Doublet strengths, $\mu^{d}_{j}$
- doubvec: double precision(nd,3,nsource)
Doublet orientation vectors, $\nu^{d}_{j}$
- ifppreg: integer
| Flag for computing velocity, pressure and/or gradients at source locations
| ifppreg = 1, compute velocity
Expand Down
2 changes: 1 addition & 1 deletion docs/genpdfmanual.sh
View file
Original file line number Diff line number Diff line change
Expand Up @@ -3,4 +3,4 @@
# Barnett 12/6/17

make latexpdf
mv _build/latex/finufft.pdf ../finufft-manual.pdf
mv _build/latex/fmm3d.pdf ../fmm3d_manual.pdf
63 changes: 57 additions & 6 deletions docs/julia.rst
View file
Original file line number Diff line number Diff line change
Expand Up @@ -204,7 +204,7 @@ denote the Stokeslet given by
(x_{3}-y_{3})^2 + \|x-y \|^2
\end{bmatrix} \, ,
and $\mathcal{T}^{\textrm{stok}}(x,y)$ denote the Stresslet whose action on
let $\mathcal{T}^{\textrm{stok}}(x,y)$ denote the Stresslet whose action on
a vector $v$ is given by

.. math::
Expand All @@ -217,27 +217,69 @@ a vector $v$ is given by
(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{bmatrix} \, ,
let $\mathcal{R}^{\textrm{stok}}(x,y)$ denote the Rotlet whose action on
a vector $v$ is given by

.. math::
v\cdot \mathcal{R}^{\textrm{stok}}(x,y) =
\frac{v \cdot (x-y)}{4\pi\|x-y \|^3}
\begin{bmatrix}
1 & 0 & 0 \\
0 & 1 & 0 \\
0 & 0 & 1
\end{bmatrix} \, ,
and $\mathcal{D}^{\textrm{stok}}(x,y)$ denote the symmetric part of Doublet whose action on
a vector $v$ is given by

.. math::
v\cdot \mathcal{D}^{\textrm{stok}}(x,y) =
\frac{3 v \cdot (x-y)}{4\pi\|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} - \\
\frac{1}{4\pi\|x-y \|^3}
\begin{bmatrix}
v_1(x_{1}-y_{1}) & v_2(x_{1}-y_{1}) & v_3(x_{1}-y_{1}) \\
v_2(x_{2}-y_{2}) & v_2(x_{2}-y_{2}) & v_3(x_{2}-y_{2}) \\
v_3(x_{3}-y_{3}) & v_3(x_{3}-y_{3}) & v_3(x_{3}-y_{3})
\end{bmatrix} \, .
This subroutine computes the N-body Stokes
interactions, its gradients and the corresponding pressure
in three dimensions given by

.. math::
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}
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} +
\nu^{r}_{j} \cdot \mathcal{R}^{\textrm{stok}}(x,x_{j}) \cdot \mu^{r}_{j} - \mu^{r}_{j} \cdot \mathcal{R}^{\textrm{stok}}(x,x_{j}) \cdot \nu^{r}_{j} + \\
\nu^{d}_{j} \cdot \mathcal{R}^{\textrm{stok}}(x,x_{j}) \cdot \mu^{d}_{j} - \mu^{d}_{j} \cdot \mathcal{R}^{\textrm{stok}}(x,x_{j}) \cdot \nu^{d}_{j} +
\nu^{d}_{j} \cdot \mathcal{D}^{\textrm{stok}}(x,x_{j}) \cdot \mu^{d}_{j} \, .
where $\sigma_{j}$ are the Stokeslet densities,
$\nu_{j}$ are the stresslet orientation vectors, $\mu_{j}$
are the stresslet densities, and
are the stresslet densities,
$\nu^{r}_{j}$ are the rotlet orientation vectors, $\mu^{r}_{j}$
are the rotlet densities,
$\nu^{d}_{j}$ are the doublet orientation vectors, $\mu^{d}_{j}$
are the doublet densities, and
$x_{j}$ are the source locations.
When $x=x_{j}$, the term corresponding to $x_{j}$ is dropped
from the sum.

.. code:: julia
vals = stfmm3d(eps,sources;stoklet=nothing,strslet=nothing,
strsvec=nothing,targets=nothing,ppreg=0,
strsvec=nothing,rotlet=nothing,rotvec=nothing,
doublet=nothing,doubvec=nothing,targets=nothing,ppreg=0,
ppregt=0,nd=1)
Wrapper for fast multipole implementation for Stokes N-body
Expand All @@ -255,6 +297,14 @@ Args:
stresslet strengths ($mu_{j}$ above)
- strsvec: float(nd,3,n) or float(3,n)
stresslet orientations ($nu_{j}$ above)
- rotlet: float(nd,3,n) or float(3,n)
rotlet strengths ($mu^{r}_{j}$ above)
- rotvec: float(nd,3,n) or float(3,n)
rotlet orientations ($nu^{r}_{j}$ above)
- doublet: float(nd,3,n) or float(3,n)
doublet strengths ($mu^{d}_{j}$ above)
- doubvec: float(nd,3,n) or float(3,n)
doublet orientations ($nu^{d}_{j}$ above)
- targets: float(3,nt)
target locations (x)
- ifppreg: integer
Expand Down Expand Up @@ -286,7 +336,8 @@ target locations.
.. code:: julia
vals = st3ddir(sources,targets;stoklet=nothing,strslet=nothing,
strsvec=nothing,ppregt=0,nd=1,thresh=1e-16)
strsvec=nothing,rotlet=nothing,rotvec=nothing,
doublet=nothing,doubvec=nothing,ppregt=0,nd=1,thresh=1e-16)
------------------------------------------------------------------

Expand Down
67 changes: 60 additions & 7 deletions docs/matlab.rst
View file
Original file line number Diff line number Diff line change
Expand Up @@ -205,7 +205,7 @@ denote the Stokeslet given by
(x_{3}-y_{3})^2 + \|x-y \|^2
\end{bmatrix} \, ,
and $\mathcal{T}^{\textrm{stok}}(x,y)$ denote the Stresslet whose action on
let $\mathcal{T}^{\textrm{stok}}(x,y)$ denote the Stresslet whose action on
a vector $v$ is given by

.. math::
Expand All @@ -218,19 +218,59 @@ a vector $v$ is given by
(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{bmatrix} \, ,
let $\mathcal{R}^{\textrm{stok}}(x,y)$ denote the Rotlet whose action on
a vector $v$ is given by

.. math::
v\cdot \mathcal{R}^{\textrm{stok}}(x,y) =
\frac{v \cdot (x-y)}{4\pi\|x-y \|^3}
\begin{bmatrix}
1 & 0 & 0 \\
0 & 1 & 0 \\
0 & 0 & 1
\end{bmatrix} \, ,
and $\mathcal{D}^{\textrm{stok}}(x,y)$ denote the symmetric part of Doublet whose action on
a vector $v$ is given by

.. math::
v\cdot \mathcal{D}^{\textrm{stok}}(x,y) =
\frac{3 v \cdot (x-y)}{4\pi\|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} - \\
\frac{1}{4\pi\|x-y \|^3}
\begin{bmatrix}
v_1(x_{1}-y_{1}) & v_2(x_{1}-y_{1}) & v_3(x_{1}-y_{1}) \\
v_2(x_{2}-y_{2}) & v_2(x_{2}-y_{2}) & v_3(x_{2}-y_{2}) \\
v_3(x_{3}-y_{3}) & v_3(x_{3}-y_{3}) & v_3(x_{3}-y_{3})
\end{bmatrix} \, .
This subroutine computes the N-body Stokes
interactions, its gradients and the corresponding pressure
in three dimensions given by
interactions, its gradients and the corresponding pressure
in three dimensions given by

.. math::
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}
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} +
\nu^{r}_{j} \cdot \mathcal{R}^{\textrm{stok}}(x,x_{j}) \cdot \mu^{r}_{j} - \mu^{r}_{j} \cdot \mathcal{R}^{\textrm{stok}}(x,x_{j}) \cdot \nu^{r}_{j} + \\
\nu^{d}_{j} \cdot \mathcal{R}^{\textrm{stok}}(x,x_{j}) \cdot \mu^{d}_{j} - \mu^{d}_{j} \cdot \mathcal{R}^{\textrm{stok}}(x,x_{j}) \cdot \nu^{d}_{j} +
\nu^{d}_{j} \cdot \mathcal{D}^{\textrm{stok}}(x,x_{j}) \cdot \mu^{d}_{j} \, .
where $\sigma_{j}$ are the Stokeslet densities,
$\nu_{j}$ are the stresslet orientation vectors, $\mu_{j}$
are the stresslet densities, and
are the stresslet densities,
$\nu^{r}_{j}$ are the rotlet orientation vectors, $\mu^{r}_{j}$
are the rotlet densities,
$\nu^{d}_{j}$ are the doublet orientation vectors, $\mu^{d}_{j}$
are the doublet densities, and
$x_{j}$ are the source locations.
When $x=x_{j}$, the term corresponding to $x_{j}$ is dropped
from the sum.
Expand Down Expand Up @@ -263,6 +303,19 @@ Args:
* srcinfo.strsvec: double(nd,3,n)
Stresslet orientiation vectors, $\nu_{j}$ (optional
default - term corresponding to stresslet dropped)
* srcinfo.rotlet: double(nd,3,n)
Rotlet densities, $\mu^{r}_{j}$ (optional
default - term corresponding to rotlet dropped)
* srcinfo.rotvec: double(nd,3,n)
Rotlet orientiation vectors, $\nu^{r}_{j}$ (optional
default - term corresponding to rotlet dropped)
* srcinfo.doublet: double(nd,3,n)
Doublet densities, $\mu^{d}_{j}$ (optional
default - term corresponding to doublet dropped)
* srcinfo.doubvec: double(nd,3,n)
Doublet orientiation vectors, $\nu^{d}_{j}$ (optional
default - term corresponding to doublet dropped)


- ifppreg: integer
| source eval flag
Expand Down
Loading

0 comments on commit c0dbe29

Please sign in to comment.