blob: f8a17df99affcab14ee2eecab42f62402664c2ea (
plain)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
|
#include <cctk.h>
#include "fasterp.hh"
namespace CarpetInterp2 {
extern "C"
CCTK_INT
CarpetInterp2_InterpGridArrays (CCTK_POINTER_TO_CONST const cctkGH_,
CCTK_INT const N_dims,
CCTK_INT const order,
CCTK_INT const N_interp_points,
CCTK_POINTER_TO_CONST const interp_coords_,
CCTK_INT const N_input_arrays,
CCTK_INT const * const input_array_indices,
CCTK_INT const N_output_arrays,
CCTK_POINTER const output_arrays_)
{
// Check input values and convert types
cGH const * const cctkGH = static_cast<cGH const *> (cctkGH_);
assert (cctkGH);
assert (N_dims == dim);
assert (N_interp_points >= 0);
CCTK_REAL const * const * const interp_coords =
static_cast<CCTK_REAL const * const *> (interp_coords_);
assert (interp_coords);
for (int d=0; d<dim; ++d) {
assert (N_interp_points==0 or interp_coords[d]);
}
assert (N_input_arrays >= 0);
assert (input_array_indices);
for (int n=0; n<N_input_arrays; ++n) {
assert (input_array_indices[n]>=0 and
input_array_indices[n]<CCTK_NumVars());
}
assert (N_output_arrays >= 0);
assert (N_output_arrays == N_input_arrays);
CCTK_REAL * const * const output_arrays =
static_cast<CCTK_REAL * const *> (output_arrays_);
assert (output_arrays);
for (int n=0; n<N_output_arrays; ++n) {
assert (N_output_arrays==0 or output_arrays[n]);
}
// Set up interpolation
fasterp_glocs_t locations (N_interp_points);
for (int d=0; d<dim; ++d) {
for (int i=0; i<N_interp_points; ++i) {
locations.coords[d].AT(i) = interp_coords[d][i];
}
}
fasterp_setup_t const setup (cctkGH, locations, order);
// Interpolate
vector<int> varinds (N_input_arrays);
for (int n=0; n<N_input_arrays; ++n) {
varinds.AT(n) = input_array_indices[n];
}
vector<CCTK_REAL *> values (N_input_arrays);
for (int n=0; n<N_input_arrays; ++n) {
values.AT(n) = output_arrays[n];
}
setup.interpolate (cctkGH, varinds, values);
// Done
return 0;
}
} // namespace CarpetInterp2
|
