fdfd_tools/meanas/fdmath/__init__.py

503 lines
24 KiB
Python
Raw Normal View History

2019-11-26 01:47:52 -08:00
"""
2019-11-30 01:24:16 -08:00
2019-11-26 01:47:52 -08:00
Basic discrete calculus for finite difference (fd) simulations.
2019-12-08 01:46:47 -08:00
TODO: short description of functional vs operator form
2019-11-30 01:24:16 -08:00
Discrete calculus
=================
2019-11-26 01:47:52 -08:00
This documentation and approach is roughly based on W.C. Chew's excellent
"Electromagnetic Theory on a Lattice" (doi:10.1063/1.355770),
which covers a superset of this material with similar notation and more detail.
2019-12-09 21:28:26 -08:00
Scalar derivatives and cell shifts
----------------------------------
2019-11-30 01:24:16 -08:00
2019-11-26 01:47:52 -08:00
Define the discrete forward derivative as
2019-12-01 02:32:31 -08:00
$$ [\\tilde{\\partial}_x f ]_{m + \\frac{1}{2}} = \\frac{1}{\\Delta_{x, m}} (f_{m + 1} - f_m) $$
2019-12-08 01:46:47 -08:00
where \\( f \\) is a function defined at discrete locations on the x-axis (labeled using \\( m \\)).
The value at \\( m \\) occupies a length \\( \\Delta_{x, m} \\) along the x-axis. Note that \\( m \\)
is an index along the x-axis, _not_ necessarily an x-coordinate, since each length
\\( \\Delta_{x, m}, \\Delta_{x, m+1}, ...\\) is independently chosen.
If we treat `f` as a 1D array of values, with the `i`-th value `f[i]` taking up a length `dx[i]`
along the x-axis, the forward derivative is
deriv_forward(f)[i] = (f[i + 1] - f[i]) / dx[i]
2019-11-26 01:47:52 -08:00
Likewise, discrete reverse derivative is
2019-12-01 02:32:31 -08:00
$$ [\\hat{\\partial}_x f ]_{m - \\frac{1}{2}} = \\frac{1}{\\Delta_{x, m}} (f_{m} - f_{m - 1}) $$
2019-11-30 01:24:16 -08:00
or
2019-11-26 01:47:52 -08:00
2019-12-08 01:46:47 -08:00
deriv_back(f)[i] = (f[i] - f[i - 1]) / dx[i]
The derivatives' values are shifted by a half-cell relative to the original function, and
will have different cell widths if all the `dx[i]` ( \\( \\Delta_{x, m} \\) ) are not
identical:
[figure: derivatives and cell sizes]
dx0 dx1 dx2 dx3 cell sizes for function
----- ----- ----------- -----
______________________________
| | | |
f0 | f1 | f2 | f3 | function
_____|_____|___________|_____|
| | | |
| Df0 | Df1 | Df2 | Df3 forward derivative (periodic boundary)
__|_____|________|________|___
dx'3] dx'0 dx'1 dx'2 [dx'3 cell sizes for forward derivative
-- ----- -------- -------- ---
dx'0] dx'1 dx'2 dx'3 [dx'0 cell sizes for reverse derivative
______________________________
| | | |
| df1 | df2 | df3 | df0 reverse derivative (periodic boundary)
__|_____|________|________|___
Periodic boundaries are used here and elsewhere unless otherwise noted.
In the above figure,
`f0 =` \\(f_0\\), `f1 =` \\(f_1\\)
`Df0 =` \\([\\tilde{\\partial}f]_{0 + \\frac{1}{2}}\\)
`Df1 =` \\([\\tilde{\\partial}f]_{1 + \\frac{1}{2}}\\)
`df0 =` \\([\\hat{\\partial}f]_{0 - \\frac{1}{2}}\\)
etc.
The fractional subscript \\( m + \\frac{1}{2} \\) is used to indicate values defined
at shifted locations relative to the original \\( m \\), with corresponding lengths
$$ \\Delta_{x, m + \\frac{1}{2}} = \\frac{1}{2} * (\\Delta_{x, m} + \\Delta_{x, m + 1}) $$
Just as \\( m \\) is not itself an x-coordinate, neither is \\( m + \\frac{1}{2} \\);
carefully note the positions of the various cells in the above figure vs their labels.
2019-12-09 21:28:26 -08:00
If the positions labeled with \\( m \\) are considered the "base" or "original" grid,
the positions labeled with \\( m + \\frac{1}{2} \\) are said to lie on a "dual" or
"derived" grid.
2019-12-08 01:46:47 -08:00
For the remainder of the `Discrete calculus` section, all figures will show
constant-length cells in order to focus on the vector derivatives themselves.
2019-12-09 21:28:26 -08:00
See the `Grid description` section below for additional information on this topic
and generalization to three dimensions.
2019-11-26 01:47:52 -08:00
2019-11-30 01:24:16 -08:00
Gradients and fore-vectors
--------------------------
2019-11-26 01:47:52 -08:00
Expanding to three dimensions, we can define two gradients
2019-12-01 02:32:31 -08:00
$$ [\\tilde{\\nabla} f]_{m,n,p} = \\vec{x} [\\tilde{\\partial}_x f]_{m + \\frac{1}{2},n,p} +
2019-11-26 01:47:52 -08:00
\\vec{y} [\\tilde{\\partial}_y f]_{m,n + \\frac{1}{2},p} +
\\vec{z} [\\tilde{\\partial}_z f]_{m,n,p + \\frac{1}{2}} $$
$$ [\\hat{\\nabla} f]_{m,n,p} = \\vec{x} [\\hat{\\partial}_x f]_{m + \\frac{1}{2},n,p} +
\\vec{y} [\\hat{\\partial}_y f]_{m,n + \\frac{1}{2},p} +
\\vec{z} [\\hat{\\partial}_z f]_{m,n,p + \\frac{1}{2}} $$
2019-12-01 02:32:31 -08:00
or
[code: gradients]
grad_forward(f)[i,j,k] = [Dx_forward(f)[i, j, k],
Dy_forward(f)[i, j, k],
Dz_forward(f)[i, j, k]]
= [(f[i + 1, j, k] - f[i, j, k]) / dx[i],
(f[i, j + 1, k] - f[i, j, k]) / dy[i],
(f[i, j, k + 1] - f[i, j, k]) / dz[i]]
grad_back(f)[i,j,k] = [Dx_back(f)[i, j, k],
Dy_back(f)[i, j, k],
Dz_back(f)[i, j, k]]
= [(f[i, j, k] - f[i - 1, j, k]) / dx[i],
(f[i, j, k] - f[i, j - 1, k]) / dy[i],
(f[i, j, k] - f[i, j, k - 1]) / dz[i]]
2019-11-26 01:47:52 -08:00
The three derivatives in the gradient cause shifts in different
directions, so the x/y/z components of the resulting "vector" are defined
at different points: the x-component is shifted in the x-direction,
y in y, and z in z.
We call the resulting object a "fore-vector" or "back-vector", depending
on the direction of the shift. We write it as
$$ \\tilde{g}_{m,n,p} = \\vec{x} g^x_{m + \\frac{1}{2},n,p} +
\\vec{y} g^y_{m,n + \\frac{1}{2},p} +
\\vec{z} g^z_{m,n,p + \\frac{1}{2}} $$
$$ \\hat{g}_{m,n,p} = \\vec{x} g^x_{m - \\frac{1}{2},n,p} +
\\vec{y} g^y_{m,n - \\frac{1}{2},p} +
\\vec{z} g^z_{m,n,p - \\frac{1}{2}} $$
2019-12-01 02:32:31 -08:00
[figure: gradient / fore-vector]
(m, n+1, p+1) ______________ (m+1, n+1, p+1)
/: /|
/ : / |
/ : / |
(m, n, p+1)/_____________/ | The forward derivatives are defined
| : | | at the Dx, Dy, Dz points,
| :.........|...| but the forward-gradient fore-vector
2019-12-09 21:28:26 -08:00
z y Dz / | / is the set of all three
|/_x | Dy | / and is said to be "located" at (m,n,p)
|/ |/
(m, n, p)|_____Dx______| (m+1, n, p)
2019-11-30 01:24:16 -08:00
Divergences
-----------
2019-11-26 01:47:52 -08:00
There are also two divergences,
$$ d_{n,m,p} = [\\tilde{\\nabla} \\cdot \\hat{g}]_{n,m,p}
= [\\tilde{\\partial}_x g^x]_{m,n,p} +
[\\tilde{\\partial}_y g^y]_{m,n,p} +
[\\tilde{\\partial}_z g^z]_{m,n,p} $$
$$ d_{n,m,p} = [\\hat{\\nabla} \\cdot \\tilde{g}]_{n,m,p}
= [\\hat{\\partial}_x g^x]_{m,n,p} +
[\\hat{\\partial}_y g^y]_{m,n,p} +
[\\hat{\\partial}_z g^z]_{m,n,p} $$
2019-12-01 02:32:31 -08:00
or
[code: divergences]
div_forward(g)[i,j,k] = Dx_forward(gx)[i, j, k] +
Dy_forward(gy)[i, j, k] +
Dz_forward(gz)[i, j, k]
= (gx[i + 1, j, k] - gx[i, j, k]) / dx[i] +
(gy[i, j + 1, k] - gy[i, j, k]) / dy[i] +
(gz[i, j, k + 1] - gz[i, j, k]) / dz[i]
div_back(g)[i,j,k] = Dx_back(gx)[i, j, k] +
Dy_back(gy)[i, j, k] +
Dz_back(gz)[i, j, k]
= (gx[i, j, k] - gx[i - 1, j, k]) / dx[i] +
(gy[i, j, k] - gy[i, j - 1, k]) / dy[i] +
(gz[i, j, k] - gz[i, j, k - 1]) / dz[i]
where `g = [gx, gy, gz]` is a fore- or back-vector field.
2019-11-26 01:47:52 -08:00
Since we applied the forward divergence to the back-vector (and vice-versa), the resulting scalar value
is defined at the back-vector's (fore-vectors) location \\( (m,n,p) \\) and not at the locations of its components
\\( (m \\pm \\frac{1}{2},n,p) \\) etc.
2019-12-01 02:32:31 -08:00
[figure: divergence]
^^
(m-1/2, n+1/2, p+1/2) _____||_______ (m+1/2, n+1/2, p+1/2)
/: || ,, /|
/ : || // / | The divergence at (m, n, p) (the center
/ : // / | of this cube) of a fore-vector field
(m-1/2, n-1/2, p+1/2)/_____________/ | is the sum of the outward-pointing
| : | | fore-vector components, which are
2019-12-09 21:28:26 -08:00
z y <==|== :.........|.====> located at the face centers.
|/_x | / | /
| / // | / Note that in a nonuniform grid, each
|/ // || |/ dimension is normalized by the cell width.
(m-1/2, n-1/2, p-1/2)|____//_______| (m+1/2, n-1/2, p-1/2)
2019-12-01 02:32:31 -08:00
'' ||
VV
2019-11-26 01:47:52 -08:00
2019-11-30 01:24:16 -08:00
Curls
-----
2019-11-26 01:47:52 -08:00
The two curls are then
2019-11-30 01:24:16 -08:00
2019-12-01 02:32:31 -08:00
$$ \\begin{align*}
2019-11-26 01:47:52 -08:00
\\hat{h}_{m + \\frac{1}{2}, n + \\frac{1}{2}, p + \\frac{1}{2}} &= \\\\
[\\tilde{\\nabla} \\times \\tilde{g}]_{m + \\frac{1}{2}, n + \\frac{1}{2}, p + \\frac{1}{2}} &=
\\vec{x} (\\tilde{\\partial}_y g^z_{m,n,p + \\frac{1}{2}} - \\tilde{\\partial}_z g^y_{m,n + \\frac{1}{2},p}) \\\\
&+ \\vec{y} (\\tilde{\\partial}_z g^x_{m + \\frac{1}{2},n,p} - \\tilde{\\partial}_x g^z_{m,n,p + \\frac{1}{2}}) \\\\
2019-12-01 02:32:31 -08:00
&+ \\vec{z} (\\tilde{\\partial}_x g^y_{m,n + \\frac{1}{2},p} - \\tilde{\\partial}_y g^z_{m + \\frac{1}{2},n,p})
\\end{align*} $$
2019-11-30 01:24:16 -08:00
and
2019-11-26 01:47:52 -08:00
$$ \\tilde{h}_{m - \\frac{1}{2}, n - \\frac{1}{2}, p - \\frac{1}{2}} =
[\\hat{\\nabla} \\times \\hat{g}]_{m - \\frac{1}{2}, n - \\frac{1}{2}, p - \\frac{1}{2}} $$
where \\( \\hat{g} \\) and \\( \\tilde{g} \\) are located at \\((m,n,p)\\)
with components at \\( (m \\pm \\frac{1}{2},n,p) \\) etc.,
2019-11-28 01:27:10 -08:00
while \\( \\hat{h} \\) and \\( \\tilde{h} \\) are located at \\((m \\pm \\frac{1}{2}, n \\pm \\frac{1}{2}, p \\pm \\frac{1}{2})\\)
2019-11-26 01:47:52 -08:00
with components at \\((m, n \\pm \\frac{1}{2}, p \\pm \\frac{1}{2})\\) etc.
2019-12-01 02:32:31 -08:00
[code: curls]
curl_forward(g)[i,j,k] = [Dy_forward(gz)[i, j, k] - Dz_forward(gy)[i, j, k],
Dz_forward(gx)[i, j, k] - Dx_forward(gz)[i, j, k],
Dx_forward(gy)[i, j, k] - Dy_forward(gx)[i, j, k]]
curl_back(g)[i,j,k] = [Dy_back(gz)[i, j, k] - Dz_back(gy)[i, j, k],
Dz_back(gx)[i, j, k] - Dx_back(gz)[i, j, k],
Dx_back(gy)[i, j, k] - Dy_back(gx)[i, j, k]]
For example, consider the forward curl, at (m, n, p), of a back-vector field `g`, defined
on a grid containing (m + 1/2, n + 1/2, p + 1/2).
The curl will be a fore-vector, so its z-component will be defined at (m, n, p + 1/2).
Take the nearest x- and y-components of `g` in the xy plane where the curl's z-component
is located; these are
[curl components]
(m, n + 1/2, p + 1/2) : x-component of back-vector at (m + 1/2, n + 1/2, p + 1/2)
(m + 1, n + 1/2, p + 1/2) : x-component of back-vector at (m + 3/2, n + 1/2, p + 1/2)
(m + 1/2, n , p + 1/2) : y-component of back-vector at (m + 1/2, n + 1/2, p + 1/2)
(m + 1/2, n + 1 , p + 1/2) : y-component of back-vector at (m + 1/2, n + 3/2, p + 1/2)
These four xy-components can be used to form a loop around the curl's z-component; its magnitude and sign
is set by their loop-oriented sum (i.e. two have their signs flipped to complete the loop).
[figure: z-component of curl]
2019-12-09 21:28:26 -08:00
: |
z y : ^^ |
|/_x :....||.<.....| (m, n+1, p+1/2)
/ || /
| v || | ^
|/ |/
(m, n, p+1/2) |_____>______| (m+1, n, p+1/2)
2019-12-01 02:32:31 -08:00
Maxwell's Equations
===================
If we discretize both space (m,n,p) and time (l), Maxwell's equations become
$$ \\begin{align*}
2019-12-08 01:46:47 -08:00
\\tilde{\\nabla} \\times \\tilde{E}_{l,\\vec{r}} &= -\\tilde{\\partial}_t \\hat{B}_{l-\\frac{1}{2}, \\vec{r} + \\frac{1}{2}}
+ \\hat{M}_{l-1, \\vec{r} + \\frac{1}{2}} \\\\
\\hat{\\nabla} \\times \\hat{H}_{l,\\vec{r} + \\frac{1}{2}} &= \\hat{\\partial}_t \\tilde{D}_{l, \\vec{r}}
+ \\tilde{J}_{l-\\frac{1}{2},\\vec{r}} \\\\
2019-12-01 02:32:31 -08:00
\\tilde{\\nabla} \\cdot \\hat{B}_{l-\\frac{1}{2}, \\vec{r} + \\frac{1}{2}} &= 0 \\\\
\\hat{\\nabla} \\cdot \\tilde{D}_{l,\\vec{r}} &= \\rho_{l,\\vec{r}}
\\end{align*} $$
with
$$ \\begin{align*}
\\hat{B}_\\vec{r} &= \\mu_{\\vec{r} + \\frac{1}{2}} \\cdot \\hat{H}_{\\vec{r} + \\frac{1}{2}} \\\\
\\tilde{D}_\\vec{r} &= \\epsilon_\\vec{r} \\cdot \\tilde{E}_\\vec{r}
\\end{align*} $$
where the spatial subscripts are abbreviated as \\( \\vec{r} = (m, n, p) \\) and
2019-12-08 01:46:47 -08:00
\\( \\vec{r} + \\frac{1}{2} = (m + \\frac{1}{2}, n + \\frac{1}{2}, p + \\frac{1}{2}) \\),
\\( \\tilde{E} \\) and \\( \\hat{H} \\) are the electric and magnetic fields,
\\( \\tilde{J} \\) and \\( \\hat{M} \\) are the electric and magnetic current distributions,
and \\( \\epsilon \\) and \\( \\mu \\) are the dielectric permittivity and magnetic permeability.
The above is Yee's algorithm, written in a form analogous to Maxwell's equations.
The time derivatives can be expanded to form the update equations:
[code: Maxwell's equations]
H[i, j, k] -= (curl_forward(E[t])[i, j, k] - M[t, i, j, k]) / mu[i, j, k]
E[i, j, k] += (curl_back( H[t])[i, j, k] + J[t, i, j, k]) / epsilon[i, j, k]
Note that the E-field fore-vector and H-field back-vector are offset by a half-cell, resulting
in distinct locations for all six E- and H-field components:
[figure: Yee cell]
(m, n+1, p+1) _________________________ (m+1, n+1, p+1)
/: /|
2019-12-09 21:28:26 -08:00
z y / : / |
|/_x / : / | Locations of the
2019-12-08 01:46:47 -08:00
/ : / | E- and H-field components
/ : / | for the E fore-vector at
/ : / | r = (m, n, p) and its associated
(m, n, p+1)/________________________/ | H back-vector at r + 1/2 =
| : | | (m + 1/2, n + 1/2, p + 1/2)
| : | | (the large cube's center)
| Hx : | |
| /: :.................|......| (m+1, n+1, p)
2019-12-09 21:28:26 -08:00
|/ : / | /
Ez..........Hy | /
| Ey.......:..Hz | / This is the Yee discretization
| / : / | / scheme ("Yee cell").
| / : / | /
|/ :/ |/
r=(m, n, p)|___________Ex___________| (m+1, n, p)
2019-12-08 01:46:47 -08:00
Each component forms its own grid, offset from the others:
[figure: E-fields for adjacent cells]
________Ex(p+1, m+1)_____
2019-12-09 21:28:26 -08:00
z y /: /|
|/_x / : / |
2019-12-08 01:46:47 -08:00
/ : / |
Ey(p+1) Ey(m+1, p+1)
/ : / |
/ Ez(n+1) / Ez(m+1, n+1)
/__________Ex(p+1)_______/ |
| : | |
| : | | This figure shows which fore-vector
| : | | each e-field component belongs to.
| :.........Ex(n+1).|......| Indices are shortened; e.g. Ex(p+1)
2019-12-09 21:28:26 -08:00
| / | / means "Ex for the fore-vector located
Ez / Ez(m+1) at (m, n, p+1)".
| Ey | /
| / | Ey(m+1)
| / | /
|/ |/
r=(m, n, p)|___________Ex___________|
2019-12-08 01:46:47 -08:00
2019-12-01 02:32:31 -08:00
The divergence equations can be derived by taking the divergence of the curl equations
and combining them with charge continuity,
$$ \\hat{\\nabla} \\cdot \\tilde{J} + \\hat{\\partial}_t \\rho = 0 $$
implying that the discrete Maxwell's equations do not produce spurious charges.
Wave equation
-------------
2019-12-08 01:46:47 -08:00
Taking the backward curl of the \\( \\tilde{\\nabla} \\times \\tilde{E} \\) equation and
replacing the resulting \\( \\hat{\\nabla} \\times \\hat{H} \\) term using its respective equation,
and setting \\( \\hat{M} \\) to zero, we can form the discrete wave equation:
2019-12-01 02:32:31 -08:00
$$
2019-12-08 01:46:47 -08:00
\\begin{align*}
\\tilde{\\nabla} \\times \\tilde{E}_{l,\\vec{r}} &=
-\\tilde{\\partial}_t \\hat{B}_{l-\\frac{1}{2}, \\vec{r} + \\frac{1}{2}}
+ \\hat{M}_{l-1, \\vec{r} + \\frac{1}{2}} \\\\
\\mu^{-1}_{\\vec{r} + \\frac{1}{2}} \\cdot \\tilde{\\nabla} \\times \\tilde{E}_{l,\\vec{r}} &=
-\\tilde{\\partial}_t \\hat{H}_{l-\\frac{1}{2}, \\vec{r} + \\frac{1}{2}} \\\\
\\hat{\\nabla} \\times (\\mu^{-1}_{\\vec{r} + \\frac{1}{2}} \\cdot \\tilde{\\nabla} \\times \\tilde{E}_{l,\\vec{r}}) &=
\\hat{\\nabla} \\times (-\\tilde{\\partial}_t \\hat{H}_{l-\\frac{1}{2}, \\vec{r} + \\frac{1}{2}}) \\\\
\\hat{\\nabla} \\times (\\mu^{-1}_{\\vec{r} + \\frac{1}{2}} \\cdot \\tilde{\\nabla} \\times \\tilde{E}_{l,\\vec{r}}) &=
-\\tilde{\\partial}_t \\hat{\\nabla} \\times \\hat{H}_{l-\\frac{1}{2}, \\vec{r} + \\frac{1}{2}} \\\\
\\hat{\\nabla} \\times (\\mu^{-1}_{\\vec{r} + \\frac{1}{2}} \\cdot \\tilde{\\nabla} \\times \\tilde{E}_{l,\\vec{r}}) &=
-\\tilde{\\partial}_t \\hat{\\partial}_t \\epsilon_\\vec{r} \\tilde{E}_{l, \\vec{r}} + \\hat{\\partial}_t \\tilde{J}_{l-\\frac{1}{2},\\vec{r}} \\\\
\\hat{\\nabla} \\times (\\mu^{-1}_{\\vec{r} + \\frac{1}{2}} \\cdot \\tilde{\\nabla} \\times \\tilde{E}_{l, \\vec{r}})
2019-12-01 02:32:31 -08:00
+ \\tilde{\\partial}_t \\hat{\\partial}_t \\epsilon_\\vec{r} \\cdot \\tilde{E}_{l, \\vec{r}}
2019-12-08 01:46:47 -08:00
&= \\tilde{\\partial}_t \\tilde{J}_{l - \\frac{1}{2}, \\vec{r}}
\\end{align*}
$$
2019-12-01 02:32:31 -08:00
Grid description
================
2019-12-08 01:46:47 -08:00
2019-12-09 21:28:26 -08:00
As described in the section on scalar discrete derivatives above, cell widths along
each axis can be arbitrary and independently defined. Moreover, all field components
are defined at "derived" or "dual" positions, in between the "base" grid points on
one or more axes.
[figure: 3D base and derived grids]
_____________________________ _____________________________
z y /: /: /: /| z y /: /: /:
|/_x / : / : / : / | |/_x / : / : / :
/ : / : / : / | / : / : / :
/___________________________/ | dz[1] ________________________/____
/ : / : / : /| | /: : / : /: : dz[1]
/: : / : / : / | | / : : / : / : :
/ : :..../......:/......:/..|...| / .:...:../......:/..:...:.....
/___________/_______/_______/ | /| ______/_________/_______/___: :
| : / : | | | | / | | : : | | : :
| : / : | | | |/ | | : : | | : :
| :/ : | | | | | dz[0] | : : | | : : dz[0]
| / : | | | /| | | : : | | : :
| /: :...|.......|.......|./ |...| | ..:...:.|.......|...:...:.....
|/ : / | /| /|/ | / | : / | /| : /
|___________|_______|_______| | / dy[1] ______|_________|_______|___: / dy[1]
| : / | / | / | |/ | :/ | / | :/
| :/.......|.../...|.../...|...| ..|...:.....|.../...|...:...
| / | / | / | / | / | / | / dy[0]
| / | / | / | / dy[0] | / | / | /
|/ |/ |/ |/ |/ |/ |/
|___________|_______|_______| ______|_________|_______|___
dx[0] dx[1] dx[2] dx'[0] dx'[1] dx'[2]
Base grid Shifted one half-cell right
(e.g. for 1D forward x derivative of all components)
Some lines are omitted for clarity.
z y : / : / :dz'[1]
|/_x :/ :/ :/
.......:..........:.......:...
| /: | /: | /:
| / : | / : | / :
|/ : |/ : |/ :dz'[0]
______________________________
/| :/ /| :/ /| :/dy'[1]
/.|...:..../.|...:./.|.. :....
| /: | /: | /:
| / : | / : | /dy'[0]
|/ : |/ : |/ :
_______________________________
/| /| /|
/ | / | / |
| | |
dx'[0] dx'[1] dx'[2]
All three dimensions shifted by one half-
cell. This is quite hard to visualize
(and probably not entirely to scale).
Nevertheless, while the spacing
[figure: Component centers]
___________________________________________
z y /: /: /|
|/_x / : / : / |
/ : / : / |
Ey...........Hz Ey.....Hz / |
/ : / / : / / |
/ : / / : / / |
/ : / / :/ / |
/___________Ex____________/______Ex________/ |
| : | : | |
| : | : | |
| Hx : | Hx : | Hx |
| /: :.................|../:...:........|../:...|
| / : / | / : / | / : /
|/ : / |/ : / |/ : /
Ez...........Hy Ez......Hy Ez :/
| Ey........:..Hz | Ey...:..Hz | Ey
| / : / | / : / | /
| / : / | / : / | /
|/ :/ |/ :/ |/
|___________Ex____________|_______Ex_______|
Part of a nonuniform "base grid", with labels specifying
positions of the various field components.
z y mmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmm
|/_x m: m:
m : m :
Ey...........m..:.........Ey......m..:.....Ey
m : m :
m : m :
_____m_____:______________m_____:________
mmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmm
| / : | / :
| / : | / :
| / : | / :
wwww|w/wwww:wwwwwwwwwwwww|w/wwww:wwwwwwwwww
|/ w |/ w
____________|____________________|__________
Ey.......|...w............Ey..|...w........Ey
| w | w
| w | w
|w |w
wwwwwwwwwwww|wwwwwwwwwwwwwwwwwwww|wwwwwwwwww
The Ey values are positioned on the y-edges of the base
grid, but they represent the Ey field in a volume that
contains (but isn't necessarily centered on) the points
at which they are defined.
Here, the 'Ey' labels represent the same points as before;
the grid lines _|:/ are edges of the area represented
by each Ey value, and the lines drawn using m.w represent
areas where a cell's faces extend beyond the drawn area
(i.e. where the drawing is truncated in the z-direction).
2019-12-08 01:46:47 -08:00
2019-11-27 22:59:52 -08:00
TODO: explain dxes
2019-11-26 01:47:52 -08:00
"""
2019-11-27 22:59:52 -08:00
from .types import fdfield_t, vfdfield_t, dx_lists_t, fdfield_updater_t
from .vectorization import vec, unvec
from . import operators, functional, types, vectorization