QPMS
/
qpms
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity

Docstrings, array dtypes, fix _sswk "constructor"

This commit is contained in:
Marek Nečada 2020-07-20 12:18:01 +03:00
parent a68e1d8f8c
commit 9e42520371
1 changed files with 56 additions and 6 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

62
qpms/qpms_c.pyx
View File

@ -947,13 +947,16 @@ cdef class ScatteringSystem:
return retdict
@boundscheck(False)
def scattered_E(self, cdouble wavenumber, scatcoeffvector_full, evalpos, bint alt=False, btyp=BesselType.HANKEL_PLUS):
def scattered_E(self, cdouble wavenumber, scatcoeffvector_full, evalpos, btyp=BesselType.HANKEL_PLUS, bint alt=False):
# TODO DOC, periodic systems
if self.lattice_dimension != 0: # TODO enable also periodic systems
raise NotImplementedError("For periodic system, use _ScatteringSystemAtOmegaK.scattered_E")
cdef qpms_bessel_t btyp_c = BesselType(btyp)
evalpos = np.array(evalpos, dtype=float, copy=False)
if evalpos.shape[-1] != 3:
raise ValueError("Last dimension of evalpos has to be 3")
cdef np.ndarray[double,ndim=2] evalpos_a = evalpos.reshape(-1,3)
cdef np.ndarray[dtype=complex, ndim=1] scv = np.array(scatcoeffvector_full, copy=False)
cdef np.ndarray[dtype=complex, ndim=1] scv = np.array(scatcoeffvector_full, copy=False, dtype=complex)
cdef cdouble[::1] scv_view = scv
cdef np.ndarray[complex, ndim=2] results = np.empty((evalpos_a.shape[0],3), dtype=complex)
cdef ccart3_t res
@ -1003,6 +1006,8 @@ def empty_lattice_modes_xy(EpsMu epsmu, reciprocal_basis, wavevector, double max
cdef class _ScatteringSystemAtOmegaK:
'''
Wrapper over the C qpms_scatsys_at_omega_k_t structure
This represents an infinite periodic system with a given (fixed) frequency and Bloch vector.
'''
cdef qpms_scatsys_at_omega_k_t sswk
cdef _ScatteringSystemAtOmega ssw_pyref
@ -1018,7 +1023,22 @@ cdef class _ScatteringSystemAtOmegaK:
self.sswk.eta = eta
@boundscheck(False)
def scattered_E(self, scatcoeffvector_full, evalpos, btyp=QPMS_HANKEL_PLUS): # TODO DOC!!!
def scattered_E(self, scatcoeffvector_full, evalpos, btyp=QPMS_HANKEL_PLUS):
"""Evaluate electric field for a given excitation coefficient vector (periodic system)
Parameters
----------
scatcoeffvector_full: array_like of shape (self.fecv_size,)
Onedimensional excitation coefficient vector, describing SVWFs leaving
the scatterers inside the representative unit cell.
evalpos: array_like of floats and shape (..., 3)
Evaluation points in cartesian coordinates.
Returns
-------
array_like of complex, with the same shape as `evalpos`
Electric field at the positions given in `evalpos`.
"""
if(btyp != QPMS_HANKEL_PLUS):
raise NotImplementedError("Only first kind Bessel function-based fields are supported")
cdef qpms_bessel_t btyp_c = BesselType(btyp)
@ -1026,7 +1046,7 @@ cdef class _ScatteringSystemAtOmegaK:
if evalpos.shape[-1] != 3:
raise ValueError("Last dimension of evalpos has to be 3")
cdef np.ndarray[double,ndim=2] evalpos_a = evalpos.reshape(-1,3)
cdef np.ndarray[dtype=complex, ndim=1] scv = np.array(scatcoeffvector_full, copy=False)
cdef np.ndarray[dtype=complex, ndim=1] scv = np.array(scatcoeffvector_full, dtype=complex, copy=False)
cdef cdouble[::1] scv_view = scv
cdef np.ndarray[complex, ndim=2] results = np.empty((evalpos_a.shape[0],3), dtype=complex)
cdef ccart3_t res
@ -1042,6 +1062,9 @@ cdef class _ScatteringSystemAtOmegaK:
results[i,1] = res.y
results[i,2] = res.z
return results.reshape(evalpos.shape)
property fecv_size:
def __get__(self): return self.ssw_pyref.ss_pyref.fecv_size
cdef class _ScatteringSystemAtOmega:
'''
@ -1097,6 +1120,7 @@ cdef class _ScatteringSystemAtOmega:
def _sswk(self, k):
cdef _ScatteringSystemAtOmegaK sswk = _ScatteringSystemAtOmegaK()
sswk.sswk.ssw = self.ssw
sswk.ssw_pyref=self
sswk.sswk.k[0] = k[0]
sswk.sswk.k[1] = k[1]
sswk.sswk.k[2] = k[2]
@ -1164,13 +1188,39 @@ cdef class _ScatteringSystemAtOmega:
return self.ss_pyref.translation_matrix_packed(wavenumber=self.wavenumber, iri=iri, J=J)
@boundscheck(False)
def scattered_E(self, scatcoeffvector_full, evalpos, bint alt=False, btyp=QPMS_HANKEL_PLUS):
def scattered_E(self, scatcoeffvector_full, evalpos, blochvector=None, btyp=QPMS_HANKEL_PLUS, bint alt=False):
# FIXME TODO this obviously does not work for periodic systems
"""Evaluate electric field for a given excitation coefficient vector
Parameters
----------
scatcoeffvector_full: array_like of shape (self.fecv_size,)
Onedimensional excitation coefficient vector, describing SVWFs leaving
the scatterers.
evalpos: array_like of floats and shape (..., 3)
Evaluation points in cartesian coordinates.
blochvector: array_like or None
Bloch vector, must be supplied (non-None) for periodic systems, else None.
Returns
-------
array_like of complex, with the same shape as `evalpos`
Electric field at the positions given in `evalpos`.
See Also
--------
_ScatteringSystemAtOmegaK.scattered_E: Evaluate scattered field in an infinite periodic system.
"""
if(self.ssw[0].ss[0].lattice_dimension != 0):
if blochvector is None:
raise ValueError("For a periodic system, blochvector must be specified.")
return self._sswk(blochvector).scattered_E(scatcoeffvector_full, evalpos=evalpos, btyp=btyp)
cdef qpms_bessel_t btyp_c = BesselType(btyp)
evalpos = np.array(evalpos, dtype=float, copy=False)
if evalpos.shape[-1] != 3:
raise ValueError("Last dimension of evalpos has to be 3")
cdef np.ndarray[double,ndim=2] evalpos_a = evalpos.reshape(-1,3)
cdef np.ndarray[dtype=complex, ndim=1] scv = np.array(scatcoeffvector_full, copy=False)
cdef np.ndarray[dtype=complex, ndim=1] scv = np.array(scatcoeffvector_full, dtype=complex, copy=False)
cdef cdouble[::1] scv_view = scv
cdef np.ndarray[complex, ndim=2] results = np.empty((evalpos_a.shape[0],3), dtype=complex)
cdef ccart3_t res