gauss_seidel_x¶
-
emg3d.core.
gauss_seidel_x
(ex, ey, ez, sx, sy, sz, eta_x, eta_y, eta_z, zeta, hx, hy, hz, nu)[source]¶ Gauss-Seidel method with line relaxation in x-direction.
This is the equivalent to
gauss_seidel()
, but with line relaxation in the x-direction. Seegauss_seidel()
for more details.The resulting system A x = b to solve consists of n unknowns (x-vector), and the corresponding matrix A is a banded matrix with the main diagonal and five upper and lower diagonals:
.-0 |X|\ 0 0-.-0 left: middle: right: \|X|\ (not used) 0-.-0 0- .- 0 \|X|\ \ |X |\ 0-.-0 0 \|X| 0-. . 1*1, - 4*1, | 1*4, X 4*4, \ 4*4 upper or lower
The matrix A is complex and symmetric (A = A^T), and therefore only the main diagonal and the lower five off-diagonals are required.
- The right-hand-side b has length 5*nCx-4 (nCx even).
- The matrix A has length of b and 1+2*5 diagonals; we use for it an array of length 6*len(b).
The values are computed in rows of 5 lines, with the indicated middle and left matrices as indicated in the above scheme. These blocks are filled into the main matrix A and vector b, and subsequently solved with a non-standard Cholesky factorisation,
solve()
.Tangential components at the boundaries are assumed to be 0 (PEC boundaries).
The result is stored in the provided electric fields ex, ey, and ez.
Parameters: - ex, ey, ez : ndarray
Electric fields in x-, y-, and z-directions, as obtained from
emg3d.fields.Field
.- sx, sy, sz :
Source fields in x-, y-, and z-directions, as obtained from
emg3d.fields.Field
.- eta_x, eta_y, eta_z, zeta :
VolumeModel parameters (multiplied by volumes) as obtained from
emg3d.models.VolumeModel()
.- hx, hy, hz : ndarray
Cell widths in x-, y-, and z-directions.
- nu : int
Number of Gauss-Seidel iterations.