XGCa
pseudo_inverse.hpp
Go to the documentation of this file.
1 
15 #ifndef PSEUDO_INVERSE_HPP
16 #define PSEUDO_INVERSE_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 #include <petscds.h>
28 #include "petscdmswarm.h"
29 #include "petscksp.h"
30 #endif
31 
32 #ifndef NO_PETSC
33 // For some reason Kokkos cannot handle Views of PETSc objects, but it works when putting them in a wrapper
34 // Pseudo-inverse PETSc objects
36  DM swarm;
37  Vec rho;
38  Vec rhs;
39  Mat M_p;
40 };
41 #endif
42 
43 template<class Device>
45 
46  public:
47 
48 #ifndef NO_PETSC
49  // All PETSc objects on host
50  Kokkos::View<PseudoInversePetscObjects*, Kokkos::LayoutRight, HostType> obj;
51 #endif
52 
53  // Arrays on Device and Host:
54  Kokkos::View<double***, Kokkos::LayoutRight, Device> ptl_coords;
55  Kokkos::View<double***, Kokkos::LayoutRight, Kokkos::HostSpace> ptl_coords_h;
56 
57  Kokkos::View<double***, Kokkos::LayoutRight, Device> ptl_weights;
58  Kokkos::View<double***, Kokkos::LayoutRight, Kokkos::HostSpace> ptl_weights_h;
59 
60  Kokkos::View<double***, Kokkos::LayoutRight, Device> ptl_alpha_weights;
61  Kokkos::View<double***, Kokkos::LayoutRight, Kokkos::HostSpace> ptl_alpha_weights_h;
62 
63  Kokkos::View<int***, Kokkos::LayoutRight, Device> ptl_indices;
64  Kokkos::View<int***, Kokkos::LayoutRight, Kokkos::HostSpace> ptl_indices_h;
65 
66  Kokkos::View<int**, Kokkos::LayoutRight, Device> n_ptl;
67  Kokkos::View<int**, Kokkos::LayoutRight, Kokkos::HostSpace> n_ptl_h;
68 
69  Kokkos::View<double***, Kokkos::LayoutRight, Device> grid_coords;
70  Kokkos::View<double***, Kokkos::LayoutRight, Kokkos::HostSpace> grid_coords_h;
71 
72  Kokkos::View<double***, Kokkos::LayoutRight, Device> grid_weights;
73  Kokkos::View<double***, Kokkos::LayoutRight, Kokkos::HostSpace> grid_weights_h;
74 
75  Kokkos::View<double***, Kokkos::LayoutRight, Device> grid_alpha_weights;
76  Kokkos::View<double***, Kokkos::LayoutRight, Kokkos::HostSpace> grid_alpha_weights_h;
77 
78  Kokkos::View<bool**, Kokkos::LayoutRight, Device> vgrid_filled;
79  Kokkos::View<bool**, Kokkos::LayoutRight, Kokkos::HostSpace> vgrid_filled_h;
80 
82  bool force;
86 
87  // Initialize these Views but don't allocate any memory yet; Instead, resize when needed
89  : ptl_coords("ptl_coords", 0, 0, 0),
90  ptl_coords_h("ptl_coords_h", 0, 0, 0),
91  ptl_weights("ptl_weights", 0, 0, 0),
92  ptl_weights_h("ptl_weights_h", 0, 0, 0),
93  ptl_alpha_weights("ptl_alpha_weights", 0, 0, 0),
94  ptl_alpha_weights_h("ptl_alpha_weights_h", 0, 0, 0),
95  ptl_indices("ptl_indices", 0, 0, 0),
96  ptl_indices_h("ptl_indices_h", 0, 0, 0),
97  n_ptl("n_ptl", 0, 0),
98  n_ptl_h("n_ptl_h", 0, 0),
99  grid_coords("grid_coords", 0, 0, 0),
100  grid_coords_h("grid_coords_h", 0, 0, 0),
101  grid_weights("grid_weights", 0, 0, 0),
102  grid_weights_h("grid_weights_h", 0, 0, 0),
103  grid_alpha_weights("grid_alpha_weights", 0, 0, 0),
104  grid_alpha_weights_h("grid_alpha_weights_h", 0, 0, 0),
105  vgrid_filled("vgrid_filled", 0, 0),
106  vgrid_filled_h("vgrid_filled_h", 0, 0),
108  force(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  n_nonadiabatic_species = n_non_ad;
143  }
144 
145  // Deallocate memory of the device arrays by resizing back to 0
147 #ifndef NO_PETSC
148  Kokkos::resize(obj, 0);
149 #endif
150  Kokkos::resize(ptl_coords, 0, 0, 0);
151  Kokkos::resize(ptl_coords_h, 0, 0, 0);
152  Kokkos::resize(ptl_weights, 0, 0, 0);
153  Kokkos::resize(ptl_weights_h, 0, 0, 0);
154  Kokkos::resize(ptl_alpha_weights, 0, 0, 0);
155  Kokkos::resize(ptl_alpha_weights_h, 0, 0, 0);
156  Kokkos::resize(ptl_indices, 0, 0, 0);
157  Kokkos::resize(ptl_indices_h, 0, 0, 0);
158  Kokkos::resize(n_ptl, 0, 0);
159  Kokkos::resize(n_ptl_h, 0, 0);
160  Kokkos::resize(grid_coords, 0, 0, 0);
161  Kokkos::resize(grid_coords_h, 0, 0, 0);
162  Kokkos::resize(grid_weights, 0, 0, 0);
163  Kokkos::resize(grid_weights_h, 0, 0, 0);
164  Kokkos::resize(grid_alpha_weights, 0, 0, 0);
165  Kokkos::resize(grid_alpha_weights_h, 0, 0, 0);
166  Kokkos::resize(vgrid_filled, 0, 0);
167  Kokkos::resize(vgrid_filled_h, 0, 0);
169  }
170 
171  // Deallocate memory of the alpha weights device arrays by resizing back to 0
173  Kokkos::resize(ptl_alpha_weights, 0, 0, 0);
174  Kokkos::resize(ptl_alpha_weights_h, 0, 0, 0);
175  Kokkos::resize(grid_alpha_weights, 0, 0, 0);
176  Kokkos::resize(grid_alpha_weights_h, 0, 0, 0);
177  }
178 
179  // Send arrays between device and host
181  Kokkos::deep_copy(ptl_coords_h, ptl_coords);
182  Kokkos::deep_copy(ptl_weights_h, ptl_weights);
183  Kokkos::deep_copy(ptl_alpha_weights_h, ptl_alpha_weights);
184  Kokkos::deep_copy(ptl_indices_h, ptl_indices);
185  Kokkos::deep_copy(n_ptl_h, n_ptl);
186  Kokkos::deep_copy(grid_coords_h, grid_coords);
187  Kokkos::deep_copy(grid_weights_h, grid_weights);
188  Kokkos::deep_copy(grid_alpha_weights_h, grid_alpha_weights);
189  Kokkos::deep_copy(vgrid_filled_h, vgrid_filled);
190  }
191 
193  Kokkos::deep_copy(ptl_coords, ptl_coords_h);
194  Kokkos::deep_copy(ptl_weights, ptl_weights_h);
195  Kokkos::deep_copy(ptl_alpha_weights, ptl_alpha_weights_h);
196  Kokkos::deep_copy(ptl_indices, ptl_indices_h);
197  Kokkos::deep_copy(n_ptl, n_ptl_h);
198  Kokkos::deep_copy(grid_coords, grid_coords_h);
199  Kokkos::deep_copy(grid_weights, grid_weights_h);
200  Kokkos::deep_copy(grid_alpha_weights, grid_alpha_weights_h);
201  Kokkos::deep_copy(vgrid_filled, vgrid_filled_h);
202  }
203 
205  Kokkos::deep_copy(ptl_coords_h, ptl_coords);
206  Kokkos::deep_copy(ptl_weights_h, ptl_weights);
207  Kokkos::deep_copy(ptl_alpha_weights_h, ptl_alpha_weights);
208  Kokkos::deep_copy(ptl_indices_h, ptl_indices);
209  Kokkos::deep_copy(n_ptl_h, n_ptl);
210  //Kokkos::deep_copy(vgrid_filled_h, vgrid_filled);
211  }
212 
214  Kokkos::deep_copy(grid_coords, grid_coords_h);
215  Kokkos::deep_copy(grid_weights, grid_weights_h);
216  Kokkos::deep_copy(grid_alpha_weights, grid_alpha_weights_h);
217  }
218 
219  void set_vgrid_filled(Kokkos::View<bool**,Kokkos::LayoutRight,DeviceType> vgrid_filled_in){
220  Kokkos::deep_copy(vgrid_filled, vgrid_filled_in);
221  Kokkos::deep_copy(vgrid_filled_h, vgrid_filled);
222  }
223 
225  Kokkos::deep_copy(grid_weights, 0.0);
226  Kokkos::deep_copy(grid_weights_h, 0.0);
227  }
228 };
229 
230 #ifndef NO_PETSC
233  Kokkos::View<PetscReal*, Kokkos::LayoutRight, Kokkos::HostSpace> particle_coords,
234  Kokkos::View<PetscReal*, Kokkos::LayoutRight, Kokkos::HostSpace> particle_weights,
235  const PetscInt Nparticles, const PetscInt Nspecies,
236  Kokkos::View<PetscReal*, Kokkos::LayoutRight, Kokkos::HostSpace> grid_coords,
237  Kokkos::View<PetscReal*, Kokkos::LayoutRight, Kokkos::HostSpace> grid_weights);
238 
240 PetscErrorCode pseudo_inv_create_Swarm(const DM dm, DM *sw);
241 
243 PetscErrorCode pseudo_inv_create_projection(const DM dm, DM sw,
244  const PetscInt thread_id, const PetscInt target_id,
245  Kokkos::View<PetscReal*, Kokkos::LayoutRight, Kokkos::HostSpace> particle_coords,
246  Kokkos::View<PetscReal*, Kokkos::LayoutRight, Kokkos::HostSpace> particle_weights,
247  const PetscInt Nparticles,
248  Vec rho, Mat *Mp_out);
249 
251 PetscErrorCode pseudo_inv_set_rhs(const DM dm, Mat MM,
252  Kokkos::View<PetscReal*, Kokkos::LayoutRight, Kokkos::HostSpace> grid_weights,
253  const PetscInt Ngrid, Vec rho, Vec rhs);
254 
256 PetscErrorCode pseudo_inv_grid_to_particles(const DM dm, DM sw,
257  Vec rhs, Mat M_p, Mat PM_p,
258  const PetscInt Ngrid,
259  const PetscInt Nparticles,
260  const PetscInt order,
261  Kokkos::View<PetscReal*, Kokkos::LayoutRight, Kokkos::HostSpace> particle_weights,
262  int &ksp_num_its,
263  bool &ksp_converged);
264 
265 #endif
266 //NO_PETSC
267 
270  const Grid<DeviceType>& grid,
272  Plasma& plasma,
273  const VelocityGrid& vgrid,
274  const DomainDecomposition<DeviceType>& pol_decomp,
275  Pseudo_inverse<DeviceType>& pseudo_inv);
276 
279  const Grid<DeviceType>& grid,
281  Plasma& plasma,
282  const VelocityGrid& vgrid,
283  const DomainDecomposition<DeviceType>& pol_decomp,
284  DMWrapper& pseudo_inv_dm,
285  Pseudo_inverse<DeviceType>& pseudo_inv);
286 
289  const Grid<DeviceType>& grid,
291  Plasma& plasma,
292  const VelocityGrid& vgrid,
293  const DomainDecomposition<DeviceType>& pol_decomp,
294  const std::string &msg,
295  Kokkos::View<bool**,Kokkos::LayoutRight,DeviceType> vgrid_filled);
296 
297 template<class Device> KOKKOS_INLINE_FUNCTION bool pseudo_inv_do_node(const Pseudo_inverse<Device>& pseudo_inv, int isp_non_ad, int node){
298  return (pseudo_inv.n_ptl(isp_non_ad, node) > 0 && (pseudo_inv.vgrid_filled(isp_non_ad, node) || pseudo_inv.force));
299 }
300 
301 inline bool pseudo_inv_do_node_h(Pseudo_inverse<DeviceType>& pseudo_inv, int isp_non_ad, int node){
302  return (pseudo_inv.n_ptl_h(isp_non_ad, node) > 0 && (pseudo_inv.vgrid_filled_h(isp_non_ad, node) || pseudo_inv.force));
303 }
304 #endif
Definition: magnetic_field.hpp:12
Definition: plasma.hpp:13
For description see source file.
Definition: pseudo_inverse.hpp:44
void resize_alpha_weights_to_zero()
Definition: pseudo_inverse.hpp:172
bool use_delta_grid
Whether to only map delta-grid weights.
Definition: pseudo_inverse.hpp:84
Kokkos::View< bool **, Kokkos::LayoutRight, Kokkos::HostSpace > vgrid_filled_h
host version of vgrid_filled
Definition: pseudo_inverse.hpp:79
Kokkos::View< int **, Kokkos::LayoutRight, Kokkos::HostSpace > n_ptl_h
host version of n_ptl
Definition: pseudo_inverse.hpp:67
void send_device_to_host()
Definition: pseudo_inverse.hpp:180
void resize_to_zero()
Definition: pseudo_inverse.hpp:146
Kokkos::View< double ***, Kokkos::LayoutRight, Kokkos::HostSpace > ptl_coords_h
host version of ptl_coords
Definition: pseudo_inverse.hpp:55
bool do_particles_to_grid_petsc
Whether the particles -> velocity grid interpolation should be done without or with (probably slower)...
Definition: pseudo_inverse.hpp:83
Kokkos::View< double ***, Kokkos::LayoutRight, Kokkos::HostSpace > grid_alpha_weights_h
host version of grid_alpha_weights
Definition: pseudo_inverse.hpp:76
void send_host_to_device()
Definition: pseudo_inverse.hpp:192
Kokkos::View< double ***, Kokkos::LayoutRight, Device > ptl_alpha_weights
Stores particle alpha (see sml_f0_grid_alpha) weights: weight(non-adiabatic species,...
Definition: pseudo_inverse.hpp:60
Kokkos::View< double ***, Kokkos::LayoutRight, Device > grid_weights
Stores velocity grid weights: weight(non-adiabatic species, node, grid point index)
Definition: pseudo_inverse.hpp:72
void zero_out_grid_weights()
Definition: pseudo_inverse.hpp:224
int expected_n_grid_points
Expected number of grid points returned from PETSc.
Definition: pseudo_inverse.hpp:81
Kokkos::View< double ***, Kokkos::LayoutRight, Device > ptl_coords
Stores particle velocity coordinates: v_para(non-adiabatic species, node, 2*(particle index on node))...
Definition: pseudo_inverse.hpp:54
Kokkos::View< int ***, Kokkos::LayoutRight, Device > ptl_indices
Stores particle indices: index(non-adiabatic species, node, particle index on node)
Definition: pseudo_inverse.hpp:63
Kokkos::View< double ***, Kokkos::LayoutRight, Kokkos::HostSpace > grid_weights_h
host version of grid_weights
Definition: pseudo_inverse.hpp:73
void send_ptl_device_to_host()
Definition: pseudo_inverse.hpp:204
Kokkos::View< double ***, Kokkos::LayoutRight, Device > grid_alpha_weights
Stores velocity grid alpha (see sml_f0_grid_alpha) weights: weight(non-adiabatic species,...
Definition: pseudo_inverse.hpp:75
Kokkos::View< double ***, Kokkos::LayoutRight, Device > grid_coords
Stores velocity grid coordinates: v_para(non-adiabatic species, node, 2*(grid point index)); v_perp(n...
Definition: pseudo_inverse.hpp:69
Kokkos::View< int ***, Kokkos::LayoutRight, Kokkos::HostSpace > ptl_indices_h
host version of ptl_indices
Definition: pseudo_inverse.hpp:64
Kokkos::View< double ***, Kokkos::LayoutRight, Kokkos::HostSpace > ptl_weights_h
host version of ptl_weights
Definition: pseudo_inverse.hpp:58
void set_vgrid_filled(Kokkos::View< bool **, Kokkos::LayoutRight, DeviceType > vgrid_filled_in)
Definition: pseudo_inverse.hpp:219
bool force
Whether to force pseudo-inverse even if a node has empty velocity cells.
Definition: pseudo_inverse.hpp:82
Kokkos::View< double ***, Kokkos::LayoutRight, Kokkos::HostSpace > ptl_alpha_weights_h
host version of ptl_alpha_weights
Definition: pseudo_inverse.hpp:61
Kokkos::View< bool **, Kokkos::LayoutRight, Device > vgrid_filled
Stores whether a node has no empty velocity cells (non-adiabatic species, node)
Definition: pseudo_inverse.hpp:78
Kokkos::View< double ***, Kokkos::LayoutRight, Kokkos::HostSpace > grid_coords_h
host version of grid_coords
Definition: pseudo_inverse.hpp:70
Kokkos::View< int **, Kokkos::LayoutRight, Device > n_ptl
Stores number of particles (non-adiabatic species, node)
Definition: pseudo_inverse.hpp:66
Kokkos::View< PseudoInversePetscObjects *, Kokkos::LayoutRight, HostType > obj
Stores number of threads PseudoInversePetscObjects.
Definition: pseudo_inverse.hpp:50
Pseudo_inverse()
Definition: pseudo_inverse.hpp:88
int n_nonadiabatic_species
Number of nonadiabatic species.
Definition: pseudo_inverse.hpp:85
void send_grid_host_to_device()
Definition: pseudo_inverse.hpp:213
Kokkos::View< double ***, Kokkos::LayoutRight, Device > ptl_weights
Stores particle weights: weight(non-adiabatic species, node, particle index on node)
Definition: pseudo_inverse.hpp:57
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: pseudo_inverse.hpp:120
Definition: sml.hpp:9
Definition: magnetic_field.F90:1
logical false
Definition: module.F90:102
subroutine plasma(grid, itr, p, dene_out, deni_out, Te_out, Ti_out, Vparai_out, ignore_vacuum)
Calculate the plasma density, temperature, and parallel velocity for a point in triangle itr using pl...
Definition: neutral_totalf.F90:1548
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: pseudo_inverse.cpp:652
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, const PetscInt order, Kokkos::View< PetscReal *, Kokkos::LayoutRight, Kokkos::HostSpace > particle_weights, int &ksp_num_its, bool &ksp_converged)
Definition: pseudo_inverse.cpp:463
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, 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: pseudo_inverse.cpp:147
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, const std::string &msg, Kokkos::View< bool **, Kokkos::LayoutRight, DeviceType > vgrid_filled)
Definition: pseudo_inverse.cpp:900
bool pseudo_inv_do_node_h(Pseudo_inverse< DeviceType > &pseudo_inv, int isp_non_ad, int node)
Definition: pseudo_inverse.hpp:301
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: pseudo_inverse.cpp:316
PetscErrorCode pseudo_inv_create_Swarm(const DM dm, DM *sw)
For description see source file.
Definition: pseudo_inverse.cpp:280
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)
For description see source file.
Definition: pseudo_inverse.cpp:698
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: pseudo_inverse.cpp:397
KOKKOS_INLINE_FUNCTION bool pseudo_inv_do_node(const Pseudo_inverse< Device > &pseudo_inv, int isp_non_ad, int node)
Definition: pseudo_inverse.hpp:297
Definition: DM_wrapper.hpp:34
Definition: pseudo_inverse.hpp:35
Mat M_p
PETSc projection matrix.
Definition: pseudo_inverse.hpp:39
DM swarm
PETSc particle swarm object.
Definition: pseudo_inverse.hpp:36
Vec rhs
PETSc right-hand side vector.
Definition: pseudo_inverse.hpp:38
Vec rho
PETSc helper vector.
Definition: pseudo_inverse.hpp:37
Definition: velocity_grid.hpp:8