diff options
Diffstat (limited to 'teukolsky_data.py')
-rw-r--r-- | teukolsky_data.py | 21 |
1 files changed, 21 insertions, 0 deletions
diff --git a/teukolsky_data.py b/teukolsky_data.py index 00ec5c4..a9bd490 100644 --- a/teukolsky_data.py +++ b/teukolsky_data.py @@ -170,6 +170,27 @@ class TeukolskyData(object): return metric, curv + def eval_residual(self, r, theta): + if r.ndim != 1 or theta.ndim != 1: + raise TypeError('r and theta must both be 1-dimensional NumPy arrays') + + out = np.empty((3, r.shape[0] * theta.shape[0])) + + c_out = (ctypes.POINTER(ctypes.c_double) * 3)() + for i in range(3): + c_out[i] = ctypes.cast(np.ctypeslib.as_ctypes(out[i]), ctypes.POINTER(ctypes.c_double)) + + c_r = np.ctypeslib.as_ctypes(r) + c_theta = np.ctypeslib.as_ctypes(theta) + + ret = self._libtd.td_eval_residual(self._tdctx, r.shape[0], c_r, + theta.shape[0], c_theta, c_out) + if ret < 0: + raise RuntimeError('Error evaluating the variable: %d' % ret) + + out.shape = (3, theta.shape[0], r.shape[0]) + return out + @property def amplitude(self): return self._tdctx.contents.__getattribute__('amplitude') |