From 9e4252037171f7d6c3bbc1059435590a82deec42 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Marek=20Ne=C4=8Dada?= Date: Mon, 20 Jul 2020 12:18:01 +0300 Subject: [PATCH] Docstrings, array dtypes, fix _sswk "constructor" --- qpms/qpms_c.pyx | 62 ++++++++++++++++++++++++++++++++++++++++++++----- 1 file changed, 56 insertions(+), 6 deletions(-) diff --git a/qpms/qpms_c.pyx b/qpms/qpms_c.pyx index 7ea3b9c..0b5c6be 100644 --- a/qpms/qpms_c.pyx +++ b/qpms/qpms_c.pyx @@ -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