Commit e2f3056d authored by Leon Merten Lohse's avatar Leon Merten Lohse
Browse files

Merge branch 'master' of gitlab.gwdg.de:irp/fresnel

parents c0aa297c 2676e61d
import numpy as np
import _solver
_k = 2 * np.pi # all lengths are in units of the wavelength
class Solver2d():
......@@ -9,157 +9,157 @@ class Solver2d():
Solver for 2d PDEs of the form
u_z = A * u_xx + F * u
"""
dtype = np.complex128
def __init__(self, A0, F0, u0, dz, dx):
self._nx = u0.shape[-1]
self._ones = np.ones((self._nx,), dtype = self.dtype)
self._boundary_slice = np.zeros_like(self._ones)
self.u = u0 * self._ones
self.dx = dx
self.dz = dz
self._compute_coefficients(A0, F0)
def _compute_coefficients(self, A, F):
self.rx = A * self.dz / ( 2. * self.dx**2 ) * self._ones
self.rc = F * self.dz / 2. * self._ones
def step(self, A, F, boundary):
up = self.u
rxp = self.rx
rcp = self.rc
self._compute_coefficients(A, F)
u = self._boundary_slice
rx = self.rx
rc = self.rc
# set boundary values
u[0] = boundary[0]
u[-1] = boundary[1]
self.u = _solver.step_AF(up, rxp, rcp, u, rx, rc)
return self.u
class Propagator2d(Solver2d):
def __init__(self, n0, u0, dz, dx):
A, F = self._compute_helmholtz_coefficients(n0)
super().__init__(A, F, u0, dz, dx)
def _compute_helmholtz_coefficients(self,n):
A = -1j / (2. * _k)
F = -1j * _k * (n*n - 1) / 2
return A, F
def step(self, n, boundary):
A, F = self._compute_helmholtz_coefficients(n)
return super().step(A, F, boundary)
class Solver3d():
"""
Solver for equations of the form
u_z = A * u_xx + C * u_yy + F * u
Abbreviations:
'rx' : A * dz / dx**2
'rc' : C * dz / dy**2
'rf : F * dz / 2
"""
dtype = np.complex128
def __init__(self, A0, F0, u0, dz, dy, dx):
self._nx = u0.shape[-1]
self._ny = u0.shape[-2]
self._ones = np.ones((self._ny, self._nx,), dtype = self.dtype)
self._boundary_slice = np.zeros_like(self._ones)
self.u = u0 * self._ones
self.dx = dx
self.dy = dy
self.dz = dz
self._compute_coefficients(A0, F0)
def _compute_coefficients(self, A, F):
self.rx = A * self.dz / ( 2. * self.dx * self.dx) * self._ones
self.ry = A * self.dz / ( 2. * self.dy * self.dy) * self._ones
self.rc = F * self.dz / 4. * self._ones
def step(self, A, F, boundary):
up = self.u
rxp = self.rx
ryp = self.ry
rcp = self.rc
self._compute_coefficients(A, F)
u = self._boundary_slice
rx = self.rx
ry = self.ry
rc = self.rc
# set boundary values
u[ 0,:] = boundary[0]
u[-1,:] = boundary[1]
u[:, 0] = boundary[2]
u[:,-1] = boundary[3]
self.u = _solver.step_ACF(up, rxp, ryp, rcp, u, rx, ry, rc)
return self.u
class Propagator3d(Solver3d):
def __init__(self, n0, u0, dz, dy, dx):
A, F = self._compute_helmholtz_coefficients(n0)
super().__init__(A, F, u0, dz, dy, dx)
def _compute_helmholtz_coefficients(self,n):
A = -1j / (2. * _k)
F = -1j * _k * (n*n - 1) / 2
return A, F
def step(self, n, boundary):
A, F = self._compute_helmholtz_coefficients(n)
return super().step(A, F, boundary)
\ No newline at end of file
return super().step(A, F, boundary)
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment