XGCa
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
velocity_interpolate.hpp
Go to the documentation of this file.
1 
15 #ifndef VELOCITY_INTERPOLATE_HPP
16 #define VELOCITY_INTERPOLATE_HPP
17 #include <Kokkos_Core.hpp>
18 #include "space_settings.hpp"
19 #include "sml.hpp"
20 #include "plasma.hpp"
21 #include "DM_wrapper.hpp"
22 
23 #ifndef NO_PETSC
24 #include <petscversion.h>
25 #include <petscsys.h>
26 #include <petscdmplex.h>
27 #if PETSC_VERSION_GE(3,14,0)
28 #include <petscds.h>
29 #include "petscdmswarm.h"
30 #include "petscksp.h"
31 #endif
32 #endif
33 
34 #ifndef NO_PETSC
35 // For some reason Kokkos cannot handle Views of PETSc objects, but it works when putting them in a wrapper
36 // Pseudo-inverse PETSc objects
38  DM swarm;
39  Vec rho;
40  Vec rhs;
41  Mat M_p;
42 };
43 #endif
44 
45 template<class Device>
47 
48  public:
49 
50 #ifndef NO_PETSC
51  // All PETSc objects on host
52  Kokkos::View<PseudoInversePetscObjects*, Kokkos::LayoutRight, HostType> obj;
53 #endif
54 
55  // Arrays on Device and Host:
56  Kokkos::View<double***, Kokkos::LayoutRight, Device> ptl_coords;
57  Kokkos::View<double***, Kokkos::LayoutRight, Kokkos::HostSpace> ptl_coords_h;
58 
59  Kokkos::View<double***, Kokkos::LayoutRight, Device> ptl_weights;
60  Kokkos::View<double***, Kokkos::LayoutRight, Kokkos::HostSpace> ptl_weights_h;
61 
62  Kokkos::View<double***, Kokkos::LayoutRight, Device> ptl_alpha_weights;
63  Kokkos::View<double***, Kokkos::LayoutRight, Kokkos::HostSpace> ptl_alpha_weights_h;
64 
65  Kokkos::View<int***, Kokkos::LayoutRight, Device> ptl_indices;
66  Kokkos::View<int***, Kokkos::LayoutRight, Kokkos::HostSpace> ptl_indices_h;
67 
68  Kokkos::View<int**, Kokkos::LayoutRight, Device> n_ptl;
69  Kokkos::View<int**, Kokkos::LayoutRight, Kokkos::HostSpace> n_ptl_h;
70 
71  Kokkos::View<double***, Kokkos::LayoutRight, Device> grid_coords;
72  Kokkos::View<double***, Kokkos::LayoutRight, Kokkos::HostSpace> grid_coords_h;
73 
74  Kokkos::View<double***, Kokkos::LayoutRight, Device> grid_weights;
75  Kokkos::View<double***, Kokkos::LayoutRight, Kokkos::HostSpace> grid_weights_h;
76 
77  Kokkos::View<double***, Kokkos::LayoutRight, Device> grid_alpha_weights;
78  Kokkos::View<double***, Kokkos::LayoutRight, Kokkos::HostSpace> grid_alpha_weights_h;
79 
80  Kokkos::View<bool**, Kokkos::LayoutRight, Device> vgrid_filled;
81  Kokkos::View<bool**, Kokkos::LayoutRight, Kokkos::HostSpace> vgrid_filled_h;
82 
84  bool force;
87 
88  // Initialize these Views but don't allocate any memory yet; Instead, resize when needed
90  : ptl_coords("ptl_coords", 0, 0, 0),
91  ptl_coords_h("ptl_coords_h", 0, 0, 0),
92  ptl_weights("ptl_weights", 0, 0, 0),
93  ptl_weights_h("ptl_weights_h", 0, 0, 0),
94  ptl_alpha_weights("ptl_alpha_weights", 0, 0, 0),
95  ptl_alpha_weights_h("ptl_alpha_weights_h", 0, 0, 0),
96  ptl_indices("ptl_indices", 0, 0, 0),
97  ptl_indices_h("ptl_indices_h", 0, 0, 0),
98  n_ptl("n_ptl", 0, 0),
99  n_ptl_h("n_ptl_h", 0, 0),
100  grid_coords("grid_coords", 0, 0, 0),
101  grid_coords_h("grid_coords_h", 0, 0, 0),
102  grid_weights("grid_weights", 0, 0, 0),
103  grid_weights_h("grid_weights_h", 0, 0, 0),
104  grid_alpha_weights("grid_alpha_weights", 0, 0, 0),
105  grid_alpha_weights_h("grid_alpha_weights_h", 0, 0, 0),
106  vgrid_filled("vgrid_filled", 0, 0),
107  vgrid_filled_h("vgrid_filled_h", 0, 0),
109  force(false),
111  use_delta_grid(false){
112 #ifndef NO_PETSC
113  Kokkos::resize(obj, 0);
114 #endif
115  }
116 
119  // Allocate memory for the arrays by resizing them
120  void resize_on_device_and_host(int nthreads, int n_non_ad, int nnode, int max_n_ptl_on_node, int n_grid_points){
121 #ifndef NO_PETSC
122  Kokkos::resize(obj, nthreads);
123 #endif
124  Kokkos::resize(ptl_coords, n_non_ad, nnode, 2*max_n_ptl_on_node);
125  Kokkos::resize(ptl_coords_h, n_non_ad, nnode, 2*max_n_ptl_on_node);
126  Kokkos::resize(ptl_weights, n_non_ad, nnode, max_n_ptl_on_node);
127  Kokkos::resize(ptl_weights_h, n_non_ad, nnode, max_n_ptl_on_node);
128  Kokkos::resize(ptl_alpha_weights, n_non_ad, nnode, max_n_ptl_on_node);
129  Kokkos::resize(ptl_alpha_weights_h, n_non_ad, nnode, max_n_ptl_on_node);
130  Kokkos::resize(ptl_indices, n_non_ad, nnode, max_n_ptl_on_node);
131  Kokkos::resize(ptl_indices_h, n_non_ad, nnode, max_n_ptl_on_node);
132  Kokkos::resize(n_ptl, n_non_ad, nnode);
133  Kokkos::resize(n_ptl_h, n_non_ad, nnode);
134  Kokkos::resize(grid_coords, n_non_ad, nnode, 2*n_grid_points);
135  Kokkos::resize(grid_coords_h, n_non_ad, nnode, 2*n_grid_points);
136  Kokkos::resize(grid_weights, n_non_ad, nnode, n_grid_points);
137  Kokkos::resize(grid_weights_h, n_non_ad, nnode, n_grid_points);
138  Kokkos::resize(grid_alpha_weights, n_non_ad, nnode, n_grid_points);
139  Kokkos::resize(grid_alpha_weights_h, n_non_ad, nnode, n_grid_points);
140  Kokkos::resize(vgrid_filled, n_non_ad, nnode);
141  Kokkos::resize(vgrid_filled_h, n_non_ad, nnode);
142  }
143 
144  // Deallocate memory of the device arrays by resizing back to 0
146 #ifndef NO_PETSC
147  Kokkos::resize(obj, 0);
148 #endif
149  Kokkos::resize(ptl_coords, 0, 0, 0);
150  Kokkos::resize(ptl_coords_h, 0, 0, 0);
151  Kokkos::resize(ptl_weights, 0, 0, 0);
152  Kokkos::resize(ptl_weights_h, 0, 0, 0);
153  Kokkos::resize(ptl_alpha_weights, 0, 0, 0);
154  Kokkos::resize(ptl_alpha_weights_h, 0, 0, 0);
155  Kokkos::resize(ptl_indices, 0, 0, 0);
156  Kokkos::resize(ptl_indices_h, 0, 0, 0);
157  Kokkos::resize(n_ptl, 0, 0);
158  Kokkos::resize(n_ptl_h, 0, 0);
159  Kokkos::resize(grid_coords, 0, 0, 0);
160  Kokkos::resize(grid_coords_h, 0, 0, 0);
161  Kokkos::resize(grid_weights, 0, 0, 0);
162  Kokkos::resize(grid_weights_h, 0, 0, 0);
163  Kokkos::resize(grid_alpha_weights, 0, 0, 0);
164  Kokkos::resize(grid_alpha_weights_h, 0, 0, 0);
165  Kokkos::resize(vgrid_filled, 0, 0);
166  Kokkos::resize(vgrid_filled_h, 0, 0);
167  }
168 
169  // Deallocate memory of the alpha weights device arrays by resizing back to 0
171  Kokkos::resize(ptl_alpha_weights, 0, 0, 0);
172  Kokkos::resize(ptl_alpha_weights_h, 0, 0, 0);
173  Kokkos::resize(grid_alpha_weights, 0, 0, 0);
174  Kokkos::resize(grid_alpha_weights_h, 0, 0, 0);
175  }
176 
177  // Send arrays between device and host
179  Kokkos::deep_copy(ptl_coords_h, ptl_coords);
180  Kokkos::deep_copy(ptl_weights_h, ptl_weights);
181  Kokkos::deep_copy(ptl_alpha_weights_h, ptl_alpha_weights);
182  Kokkos::deep_copy(ptl_indices_h, ptl_indices);
183  Kokkos::deep_copy(n_ptl_h, n_ptl);
184  Kokkos::deep_copy(grid_coords_h, grid_coords);
185  Kokkos::deep_copy(grid_weights_h, grid_weights);
186  Kokkos::deep_copy(grid_alpha_weights_h, grid_alpha_weights);
187  Kokkos::deep_copy(vgrid_filled_h, vgrid_filled);
188  }
189 
191  Kokkos::deep_copy(ptl_coords, ptl_coords_h);
192  Kokkos::deep_copy(ptl_weights, ptl_weights_h);
193  Kokkos::deep_copy(ptl_alpha_weights, ptl_alpha_weights_h);
194  Kokkos::deep_copy(ptl_indices, ptl_indices_h);
195  Kokkos::deep_copy(n_ptl, n_ptl_h);
196  Kokkos::deep_copy(grid_coords, grid_coords_h);
197  Kokkos::deep_copy(grid_weights, grid_weights_h);
198  Kokkos::deep_copy(grid_alpha_weights, grid_alpha_weights_h);
199  Kokkos::deep_copy(vgrid_filled, vgrid_filled_h);
200  }
201 
203  Kokkos::deep_copy(ptl_coords_h, ptl_coords);
204  Kokkos::deep_copy(ptl_weights_h, ptl_weights);
205  Kokkos::deep_copy(ptl_alpha_weights_h, ptl_alpha_weights);
206  Kokkos::deep_copy(ptl_indices_h, ptl_indices);
207  Kokkos::deep_copy(n_ptl_h, n_ptl);
208  //Kokkos::deep_copy(vgrid_filled_h, vgrid_filled);
209  }
210 
212  Kokkos::deep_copy(grid_coords, grid_coords_h);
213  Kokkos::deep_copy(grid_weights, grid_weights_h);
214  Kokkos::deep_copy(grid_alpha_weights, grid_alpha_weights_h);
215  }
216 
217  void set_vgrid_filled(Kokkos::View<bool**,Kokkos::LayoutRight,DeviceType> vgrid_filled_in){
218  Kokkos::deep_copy(vgrid_filled, vgrid_filled_in);
219  Kokkos::deep_copy(vgrid_filled_h, vgrid_filled);
220  }
221 
223  Kokkos::deep_copy(grid_weights, 0.0);
224  Kokkos::deep_copy(grid_weights_h, 0.0);
225  }
226 };
227 
228 #ifndef NO_PETSC
231  Kokkos::View<PetscReal*, Kokkos::LayoutRight, Kokkos::HostSpace> particle_coords,
232  Kokkos::View<PetscReal*, Kokkos::LayoutRight, Kokkos::HostSpace> particle_weights,
233  const PetscInt Nparticles, PetscInt order, const PetscInt Nspecies,
234  Kokkos::View<PetscReal*, Kokkos::LayoutRight, Kokkos::HostSpace> grid_coords,
235  Kokkos::View<PetscReal*, Kokkos::LayoutRight, Kokkos::HostSpace> grid_weights);
236 
238 PetscErrorCode pseudo_inv_create_Swarm(const DM dm, DM *sw);
239 
241 PetscErrorCode pseudo_inv_create_projection(const DM dm, DM sw,
242  const PetscInt thread_id, const PetscInt target_id,
243  Kokkos::View<PetscReal*, Kokkos::LayoutRight, Kokkos::HostSpace> particle_coords,
244  Kokkos::View<PetscReal*, Kokkos::LayoutRight, Kokkos::HostSpace> particle_weights,
245  const PetscInt Nparticles,
246  Vec rho, Mat *Mp_out);
247 
249 PetscErrorCode pseudo_inv_set_rhs(const DM dm, Mat MM,
250  Kokkos::View<PetscReal*, Kokkos::LayoutRight, Kokkos::HostSpace> grid_weights,
251  const PetscInt Ngrid, Vec rho, Vec rhs);
252 
254 PetscErrorCode pseudo_inv_grid_to_particles(const DM dm, DM sw,
255  Vec rhs, Mat M_p, Mat PM_p,
256  const PetscInt Ngrid,
257  const PetscInt Nparticles,
258  Kokkos::View<PetscReal*, Kokkos::LayoutRight, Kokkos::HostSpace> particle_weights,
259  bool &ksp_converged);
260 
263  const Grid<DeviceType>& grid,
265  Plasma& plasma,
266  const VelocityGrid& vgrid,
267  const DomainDecomposition<DeviceType>& pol_decomp,
268  Pseudo_inverse<DeviceType>& pseudo_inv);
269 
272  const Grid<DeviceType>& grid,
274  Plasma& plasma,
275  const VelocityGrid& vgrid,
276  const DomainDecomposition<DeviceType>& pol_decomp,
277  DMWrapper& pseudo_inv_dm,
278  Pseudo_inverse<DeviceType>& pseudo_inv);
279 #endif
280 //NO_PETSC
281 
284  const Grid<DeviceType>& grid,
286  Plasma& plasma,
287  const VelocityGrid& vgrid,
288  const DomainDecomposition<DeviceType>& pol_decomp,
289  int n_non_adi_species,
290  Kokkos::View<bool**,Kokkos::LayoutRight,DeviceType> vgrid_filled);
291 
292 template<class Device> KOKKOS_INLINE_FUNCTION bool pseudo_inv_do_node(const Pseudo_inverse<Device>& pseudo_inv, int isp_non_ad, int node){
293  return (pseudo_inv.n_ptl(isp_non_ad, node) > 0 && (pseudo_inv.vgrid_filled(isp_non_ad, node) || pseudo_inv.force));
294 }
295 
296 inline bool pseudo_inv_do_node_h(Pseudo_inverse<DeviceType>& pseudo_inv, int isp_non_ad, int node){
297  return (pseudo_inv.n_ptl_h(isp_non_ad, node) > 0 && (pseudo_inv.vgrid_filled_h(isp_non_ad, node) || pseudo_inv.force));
298 }
299 #endif
void send_host_to_device()
Definition: velocity_interpolate.hpp:190
Kokkos::View< double ***, Kokkos::LayoutRight, Device > grid_weights
Stores velocity grid weights (non-adiabatic species, node, weight)
Definition: velocity_interpolate.hpp:74
Kokkos::View< double ***, Kokkos::LayoutRight, Kokkos::HostSpace > ptl_weights_h
host version of ptl_weights
Definition: velocity_interpolate.hpp:60
PetscInt pseudo_inverse_interpolate_particles_to_velocity_grid(const DM &dm, Kokkos::View< PetscReal *, Kokkos::LayoutRight, Kokkos::HostSpace > particle_coords, Kokkos::View< PetscReal *, Kokkos::LayoutRight, Kokkos::HostSpace > particle_weights, const PetscInt Nparticles, PetscInt order, const PetscInt Nspecies, Kokkos::View< PetscReal *, Kokkos::LayoutRight, Kokkos::HostSpace > grid_coords, Kokkos::View< PetscReal *, Kokkos::LayoutRight, Kokkos::HostSpace > grid_weights)
For description see source file.
Definition: velocity_interpolate.cpp:162
For description see source file.
Definition: velocity_interpolate.hpp:46
Kokkos::View< int **, Kokkos::LayoutRight, Device > n_ptl
Stores number of particles (non-adiabatic species, node)
Definition: velocity_interpolate.hpp:68
Definition: velocity_grid.hpp:5
Definition: sml.hpp:8
Kokkos::View< int ***, Kokkos::LayoutRight, Kokkos::HostSpace > ptl_indices_h
host version of ptl_indices
Definition: velocity_interpolate.hpp:66
bool use_delta_grid
Whether to only map delta-grid weights.
Definition: velocity_interpolate.hpp:86
PetscErrorCode pseudo_inv_create_projection(const DM dm, DM sw, const PetscInt thread_id, const PetscInt target_id, Kokkos::View< PetscReal *, Kokkos::LayoutRight, Kokkos::HostSpace > particle_coords, Kokkos::View< PetscReal *, Kokkos::LayoutRight, Kokkos::HostSpace > particle_weights, const PetscInt Nparticles, Vec rho, Mat *Mp_out)
For description see source file.
Definition: velocity_interpolate.cpp:347
Definition: DM_wrapper.hpp:33
Kokkos::View< int **, Kokkos::LayoutRight, Kokkos::HostSpace > n_ptl_h
host version of n_ptl
Definition: velocity_interpolate.hpp:69
Definition: magnetic_field.hpp:9
void pseudo_inv_allocate(const Simulation< DeviceType > &sml, const Grid< DeviceType > &grid, const MagneticField< DeviceType > &magnetic_field, Plasma &plasma, const VelocityGrid &vgrid, const DomainDecomposition< DeviceType > &pol_decomp, Pseudo_inverse< DeviceType > &pseudo_inv)
For description see source file.
Definition: velocity_interpolate.cpp:711
Vec rho
PETSc helper vector.
Definition: velocity_interpolate.hpp:39
void resize_to_zero()
Definition: velocity_interpolate.hpp:145
Definition: grid.hpp:10
Mat M_p
PETSc projection matrix.
Definition: velocity_interpolate.hpp:41
Kokkos::View< double ***, Kokkos::LayoutRight, Device > grid_alpha_weights
Stores velocity grid alpha (see sml_f0_grid_alpha) weights (non-adiabatic species, node, weight)
Definition: velocity_interpolate.hpp:77
Kokkos::View< bool **, Kokkos::LayoutRight, Device > vgrid_filled
Stores whether a node has no empty velocity cells (non-adiabatic species, node)
Definition: velocity_interpolate.hpp:80
Kokkos::View< double ***, Kokkos::LayoutRight, Device > grid_coords
Stores velocity grid coordinates (non-adiabatic species, node, v_para-v_perp coordinates) ...
Definition: velocity_interpolate.hpp:71
Kokkos::View< double ***, Kokkos::LayoutRight, Device > ptl_alpha_weights
Stores particle alpha (see sml_f0_grid_alpha) weights (non-adiabatic species, node, weight)
Definition: velocity_interpolate.hpp:62
void resize_on_device_and_host(int nthreads, int n_non_ad, int nnode, int max_n_ptl_on_node, int n_grid_points)
Definition: velocity_interpolate.hpp:120
Kokkos::View< double ***, Kokkos::LayoutRight, Device > ptl_weights
Stores particle weights (non-adiabatic species, node, weight)
Definition: velocity_interpolate.hpp:59
Kokkos::View< int ***, Kokkos::LayoutRight, Device > ptl_indices
Stores particle indices (non-adiabatic species, node, index)
Definition: velocity_interpolate.hpp:65
bool force
Whether to force pseudo-inverse even if a node has empty velocity cells.
Definition: velocity_interpolate.hpp:84
PetscErrorCode pseudo_inv_set_rhs(const DM dm, Mat MM, Kokkos::View< PetscReal *, Kokkos::LayoutRight, Kokkos::HostSpace > grid_weights, const PetscInt Ngrid, Vec rho, Vec rhs)
For description see source file.
Definition: velocity_interpolate.cpp:435
bool do_particles_to_grid_petsc
Whether the particles -&gt; velocity grid interpolation should be done without or with (probably slower)...
Definition: velocity_interpolate.hpp:85
void resize_alpha_weights_to_zero()
Definition: velocity_interpolate.hpp:170
bool pseudo_inv_do_node_h(Pseudo_inverse< DeviceType > &pseudo_inv, int isp_non_ad, int node)
Definition: velocity_interpolate.hpp:296
int max_n_ptl_on_node_and_vgrid_filled(const Simulation< DeviceType > &sml, const Grid< DeviceType > &grid, const MagneticField< DeviceType > &magnetic_field, Plasma &plasma, const VelocityGrid &vgrid, const DomainDecomposition< DeviceType > &pol_decomp, int n_non_adi_species, Kokkos::View< bool **, Kokkos::LayoutRight, DeviceType > vgrid_filled)
&lt; For description see source file
Definition: velocity_interpolate.cpp:979
Kokkos::View< double ***, Kokkos::LayoutRight, Kokkos::HostSpace > ptl_coords_h
host version of ptl_coords
Definition: velocity_interpolate.hpp:57
DM swarm
PETSc particle swarm object.
Definition: velocity_interpolate.hpp:38
Kokkos::View< double ***, Kokkos::LayoutRight, Kokkos::HostSpace > grid_alpha_weights_h
host version of grid_alpha_weights
Definition: velocity_interpolate.hpp:78
Pseudo_inverse()
Definition: velocity_interpolate.hpp:89
Kokkos::View< bool **, Kokkos::LayoutRight, Kokkos::HostSpace > vgrid_filled_h
host version of vgrid_filled
Definition: velocity_interpolate.hpp:81
void send_grid_host_to_device()
Definition: velocity_interpolate.hpp:211
Definition: velocity_interpolate.hpp:37
void zero_out_grid_weights()
Definition: velocity_interpolate.hpp:222
subroutine plasma(grid, itr, p, dene_out, deni_out, Te_out, Ti_out, Vparai_out)
Calculate the plasma density, temperature, and parallel velocity for a point in triangle itr using pl...
Definition: neutral_totalf.F90:1302
Kokkos::View< double ***, Kokkos::LayoutRight, Device > ptl_coords
Stores particle velocity coordinates (non-adiabatic species, node, v_para-v_perp coordinates) ...
Definition: velocity_interpolate.hpp:56
Definition: magnetic_field.F90:1
PetscErrorCode pseudo_inv_create_Swarm(const DM dm, DM *sw)
For description see source file.
Definition: velocity_interpolate.cpp:304
Definition: domain_decomposition.hpp:11
Definition: plasma.hpp:8
Vec rhs
PETSc right-hand side vector.
Definition: velocity_interpolate.hpp:40
int expected_n_grid_points
Expected number of grid points returned from PETSc.
Definition: velocity_interpolate.hpp:83
PetscErrorCode pseudo_inv_grid_to_particles(const DM dm, DM sw, Vec rhs, Mat M_p, Mat PM_p, const PetscInt Ngrid, const PetscInt Nparticles, Kokkos::View< PetscReal *, Kokkos::LayoutRight, Kokkos::HostSpace > particle_weights, bool &ksp_converged)
For description see source file.
Definition: velocity_interpolate.cpp:506
void all_species_update_pseudo_inv_weights(const Simulation< DeviceType > &sml, const Grid< DeviceType > &grid, const MagneticField< DeviceType > &magnetic_field, Plasma &plasma, const VelocityGrid &vgrid, const DomainDecomposition< DeviceType > &pol_decomp, DMWrapper &pseudo_inv_dm, Pseudo_inverse< DeviceType > &pseudo_inv)
Definition: velocity_interpolate.cpp:754
Kokkos::View< double ***, Kokkos::LayoutRight, Kokkos::HostSpace > grid_weights_h
host version of grid_weights
Definition: velocity_interpolate.hpp:75
void set_vgrid_filled(Kokkos::View< bool **, Kokkos::LayoutRight, DeviceType > vgrid_filled_in)
Definition: velocity_interpolate.hpp:217
void send_device_to_host()
Definition: velocity_interpolate.hpp:178
KOKKOS_INLINE_FUNCTION bool pseudo_inv_do_node(const Pseudo_inverse< Device > &pseudo_inv, int isp_non_ad, int node)
Definition: velocity_interpolate.hpp:292
void send_ptl_device_to_host()
Definition: velocity_interpolate.hpp:202
Kokkos::View< PseudoInversePetscObjects *, Kokkos::LayoutRight, HostType > obj
Stores #threads PseudoInversePetscObjects.
Definition: velocity_interpolate.hpp:52
Kokkos::View< double ***, Kokkos::LayoutRight, Kokkos::HostSpace > ptl_alpha_weights_h
host version of ptl_alpha_weights
Definition: velocity_interpolate.hpp:63
Kokkos::View< double ***, Kokkos::LayoutRight, Kokkos::HostSpace > grid_coords_h
host version of grid_coords
Definition: velocity_interpolate.hpp:72