1 #ifndef AMREX_TRACERPARTICLE_MOD_K_H
2 #define AMREX_TRACERPARTICLE_MOD_K_H
4 #include <AMReX_Config.H>
31 amrex::ParticleReal * val,
int M = AMREX_SPACEDIM)
42 amrex::ParticleReal * val,
int M = AMREX_SPACEDIM)
45 int ncomp_per_array =
M;
60 amrex::ParticleReal * val,
int M = AMREX_SPACEDIM)
63 int ncomp_per_array =
M;
79 amrex::ParticleReal * val)
82 int ncomp_per_array = 1;
83 int num_arrays = AMREX_SPACEDIM;
84 IntVect is_nodal[AMREX_SPACEDIM];
85 for (
int d=0; d < AMREX_SPACEDIM; ++d) {
104 amrex::ParticleReal * val,
106 int start_comp,
int ncomp,
int num_arrays)
112 for (
int d = 0; d < num_arrays; d++)
114 AMREX_D_TERM(amrex::Real lx = (Real(p.pos(0))-plo[0])*dxi[0] -
static_cast<Real
>(!is_nodal[d][0])*Real(0.5);,
115 amrex::Real ly = (Real(p.pos(1))-plo[1])*dxi[1] -
static_cast<Real
>(!is_nodal[d][1])*Real(0.5);,
116 amrex::Real lz = (Real(p.pos(2))-plo[2])*dxi[2] -
static_cast<Real
>(!is_nodal[d][2])*Real(0.5));
121 AMREX_D_TERM(
int const i0 =
static_cast<int>(amrex::Math::floor(lx));,
122 int const j0 =
static_cast<int>(amrex::Math::floor(ly));,
123 int const k0 =
static_cast<int>(amrex::Math::floor(lz)));
125 AMREX_D_TERM(amrex::Real
const xint = lx -
static_cast<Real
>(i0);,
126 amrex::Real
const yint = ly -
static_cast<Real
>(j0);,
127 amrex::Real
const zint = lz -
static_cast<Real
>(k0));
129 amrex::Real sx[] = {amrex::Real(1.0) - xint, xint};
130 #if (AMREX_SPACEDIM > 1)
131 amrex::Real sy[] = {amrex::Real(1.0) - yint, yint};
133 #if (AMREX_SPACEDIM > 2)
134 amrex::Real sz[] = {amrex::Real(1.0) - zint, zint};
137 for (
int comp = start_comp; comp < start_comp + ncomp; ++comp) {
138 val[ctr] = ParticleReal(0.0);
139 #if (AMREX_SPACEDIM > 2)
140 for (
int kk = 0; kk <=1; ++kk) {
143 #if (AMREX_SPACEDIM > 1)
144 for (
int jj = 0; jj <= 1; ++jj) {
146 for (
int ii = 0; ii <= 1; ++ii) {
147 val[ctr] +=
static_cast<ParticleReal
>((data_arr[d])(
IntVect(
AMREX_D_DECL(i0+ii, j0+jj, k0+kk)), comp) *
165 template <
typename P>
172 amrex::ParticleReal * val,
int M = AMREX_SPACEDIM)
181 template <
typename P>
188 amrex::ParticleReal * val,
int M = AMREX_SPACEDIM)
191 int ncomp_per_array =
M;
195 val, &is_nodal, icomp, ncomp_per_array, num_arrays);
202 template <
typename P>
209 amrex::ParticleReal * val,
int M = AMREX_SPACEDIM)
212 int ncomp_per_array =
M;
216 val, &is_nodal, icomp, ncomp_per_array, num_arrays);
224 template <
typename P>
231 amrex::ParticleReal * val)
234 int ncomp_per_array = 1;
235 int num_arrays = AMREX_SPACEDIM;
236 IntVect is_nodal[AMREX_SPACEDIM];
237 for (
int d=0; d < AMREX_SPACEDIM; ++d) {
242 val, &is_nodal[0], icomp, ncomp_per_array, num_arrays);
250 template <
typename P>
255 const Array4<amrex::Real const>* data_arr,
257 amrex::ParticleReal * val,
259 int start_comp,
int ncomp,
int num_arrays)
261 #if (AMREX_SPACEDIM == 1)
262 amrex::ignore_unused(p, plo, dxi, data_arr, height_arr, val, is_nodal, start_comp, ncomp, num_arrays);
263 amrex::Abort(
" Terrain fitted grid interpolation is not supported in 1D\n");
269 for (
int d = 0; d < num_arrays; d++)
271 AMREX_D_TERM(amrex::Real lx = (Real(p.pos(0))-plo[0])*dxi[0] -
static_cast<Real
>(!is_nodal[d][0])*Real(0.5);,
272 amrex::Real ly = (Real(p.pos(1))-plo[1])*dxi[1] -
static_cast<Real
>(!is_nodal[d][1])*Real(0.5);,);
274 int const i0 =
static_cast<int>(amrex::Math::floor(lx));
277 #if (AMREX_SPACEDIM == 2)
279 int const j = p.idata(0);
281 amrex::Real hlo_xlo = amrex::Real(0.25)
282 * ( height_arr(i0 , j , k)
283 + height_arr(i0 + (!(is_nodal[d][0])) , j , k)
284 + height_arr(i0 , j + (!is_nodal[d][1]) , k)
285 + height_arr(i0 + (!(is_nodal[d][0])) , j + (!is_nodal[d][1]) , k) );
287 amrex::Real hlo_xhi = amrex::Real(0.25)
288 * ( height_arr(i0 + 1 , j , k )
289 + height_arr(i0 + 1 + (!(is_nodal[d][0])) , j , k )
290 + height_arr(i0 + 1 , j + (!is_nodal[d][1]), k )
291 + height_arr(i0 + 1 + (!(is_nodal[d][0])) , j + (!is_nodal[d][1]), k ) );
294 amrex::Real
const xint = lx -
static_cast<Real
>(i0);
295 amrex::Real sx[] = { amrex::Real(1.) - xint, xint};
296 amrex::Real height_at_px = sx[0] * hlo_xlo + sx[1] * hlo_xhi;
298 int const j0 = (amrex::Real(p.pos(1)) >= height_at_px) ? j : j-1;
302 for (
int ii=0; ii < 2; ++ii) {
303 for (
int jj=0; jj < 2; ++jj) {
304 ht[yctr] = amrex::Real(0.25)
305 * ( height_arr(i0 + ii , j0 + jj , k )
306 + height_arr(i0 + ii + (!(is_nodal[d][0])) , j0 + jj , k )
307 + height_arr(i0 + ii , j0 + jj + (!is_nodal[d][1]), k )
308 + height_arr(i0 + ii + (!(is_nodal[d][0])) , j0 + jj + (!is_nodal[d][1]), k ) );
312 amrex::Real hint_ilo = (p.pos(1) - ht[0]) / (ht[1] - ht[0]);
313 amrex::Real hint_ihi = (p.pos(1) - ht[2]) / (ht[3] - ht[2]);
315 amrex::Real sy[] = { amrex::Real(1.) - hint_ilo, amrex::Real(1.) - hint_ihi,
318 #elif (AMREX_SPACEDIM == 3)
320 int const j0 =
static_cast<int>(amrex::Math::floor(ly));
322 amrex::Real
const xint = lx -
static_cast<Real
>(i0);
323 amrex::Real
const yint = ly -
static_cast<Real
>(j0);
324 amrex::Real sx[] = { amrex::Real(1.) - xint, xint};
325 amrex::Real sy[] = { amrex::Real(1.) - yint, yint};
329 amrex::Real height_at_pxy = 0.;
330 for (
int ii = 0; ii < 2; ++ii) {
331 for (
int jj = 0; jj < 2; ++jj) {
332 hlo[ilo] = amrex::Real(0.125)
333 * ( height_arr(i0 + ii , j0 + jj , k )
334 + height_arr(i0 + ii + (!is_nodal[d][0]), j0 + jj , k )
335 + height_arr(i0 + ii , j0 + jj + (!is_nodal[d][1]), k )
336 + height_arr(i0 + ii + (!is_nodal[d][0]), j0 + jj + (!is_nodal[d][1]), k )
337 + height_arr(i0 + ii , j0 + jj , k + (!is_nodal[d][2]))
338 + height_arr(i0 + ii + (!is_nodal[d][0]), j0 + jj , k + (!is_nodal[d][2]))
339 + height_arr(i0 + ii , j0 + jj + (!is_nodal[d][1]), k + (!is_nodal[d][2]))
340 + height_arr(i0 + ii + (!is_nodal[d][0]), j0 + jj + (!is_nodal[d][1]), k + (!is_nodal[d][2]))
342 height_at_pxy += hlo[ilo] * sx[ii] * sy[jj];
347 int const k0 = (amrex::Real(p.pos(2)) >= height_at_pxy ) ? k : k-1;
351 for (
int ii = 0; ii < 2; ++ii) {
352 for (
int jj = 0; jj < 2; ++jj) {
353 for (
int kk = 0; kk < 2; ++kk) {
354 ht[zctr] = amrex::Real(0.125) *
355 ( height_arr(i0 + ii , j0 + jj , k0 + kk )
356 + height_arr(i0 + ii , j0 + jj , k0 + kk + (!is_nodal[d][2]))
357 + height_arr(i0 + ii , j0 + jj + (!is_nodal[d][1]), k0 + kk )
358 + height_arr(i0 + ii , j0 + jj + (!is_nodal[d][1]), k0 + kk + (!is_nodal[d][2]))
359 + height_arr(i0 + ii + (!is_nodal[d][0]), j0 + jj , k0 + kk )
360 + height_arr(i0 + ii + (!is_nodal[d][0]), j0 + jj , k0 + kk + (!is_nodal[d][2]))
361 + height_arr(i0 + ii + (!is_nodal[d][0]), j0 + jj + (!is_nodal[d][1]), k0 + kk )
362 + height_arr(i0 + ii + (!is_nodal[d][0]), j0 + jj + (!is_nodal[d][1]), k0 + kk + (!is_nodal[d][2]))
367 amrex::Real hint_ilojlo = ( p.pos(2) - ht[0] ) / (ht[1] - ht[0]);
368 amrex::Real hint_ilojhi = ( p.pos(2) - ht[2] ) / (ht[3] - ht[2]);
369 amrex::Real hint_ihijlo = ( p.pos(2) - ht[4] ) / (ht[5] - ht[4]);
370 amrex::Real hint_ihijhi = ( p.pos(2) - ht[6] ) / (ht[7] - ht[6]);
372 amrex::Real sz[] = { amrex::Real(1.) - hint_ilojlo,
373 amrex::Real(1.) - hint_ihijlo,
374 amrex::Real(1.) - hint_ilojhi,
375 amrex::Real(1.) - hint_ihijhi,
381 for (
int comp = start_comp; comp < start_comp + ncomp; ++comp) {
382 val[ctr] = amrex::ParticleReal(0.);
383 #if (AMREX_SPACEDIM == 2)
386 for (
int jj = 0; jj <= 1; ++jj) {
387 for (
int ii = 0; ii <=1; ++ii) {
388 val[ctr] +=
static_cast<ParticleReal
>( (data_arr[d])(i0+ii, j0+jj, k0 ,comp)*sx[ii]*sy[sy_ctr] );
392 #elif (AMREX_SPACEDIM == 3)
394 for (
int kk = 0; kk <= 1; ++kk) {
395 for (
int jj = 0; jj <= 1; ++jj) {
396 for (
int ii = 0; ii <= 1; ++ii) {
397 val[ctr] +=
static_cast<ParticleReal
>(
398 (data_arr[d])(i0+ii, j0+jj, k0 + kk, comp)*sx[ii]*sy[jj]*sz[sz_ctr]);
424 amrex::Real amult,apm,apn;
425 for (
int n=1; n<neq+1; n++)
430 for (
int i=1; i<neq; i++)
437 for (
int j=ip1; j<neq+1; j++)
450 for (
int l=ip1; l<neq+1; l++)
455 for (
int kk = ip1; kk < neq+1; kk++)
457 a(n,kk) -= amult*a(j,kk);
464 void particle_interp_solve (Array2D<Real,1,7,1,7>& a, Array1D<Real,1,7>& b, Array1D<int,1,7>& ip,
int neq)
467 for (
int l=2; l<neq+1; l++)
471 for (
int k=1; k<l; k++){
472 b(n) -= a(n,k) * b(ip(k));
475 b(ip(neq))=b(ip(neq))/a(ip(neq),neq);
477 for (
int l=1; l<neq; l++)
482 for (
int k=jp1; k<neq+1; k++){
483 b(n) -= a(n,k ) * b(ip(k));
487 for (
int n=1; n<neq+1; n++)
504 template <
typename P>
509 amrex::ParticleReal * val,
int M = AMREX_SPACEDIM)
512 int ncomp_per_array =
M;
516 val, &is_nodal, icomp, ncomp_per_array, num_arrays);
524 template <
typename P>
527 const Array4<amrex::Real const>* data_arr,
529 amrex::ParticleReal * val,
531 int start_comp,
int ncomp,
int num_arrays)
533 #if (AMREX_SPACEDIM == 1)
535 amrex::Abort(
" General mapped interpolation is not supported in 1D\n");
540 for (
int d = 0; d < num_arrays; d++)
545 #if (AMREX_SPACEDIM == 3)
553 for (
int comp = start_comp; comp < start_comp + ncomp; ++comp) {
555 #if (AMREX_SPACEDIM == 2)
557 amrex::Real f1 = (data_arr[n_of_data])(i ,j ,0,comp);
558 amrex::Real f2 = (data_arr[n_of_data])(i ,j+1,0,comp) - f1;
559 amrex::Real f3 = (data_arr[n_of_data])(i+1,j ,0,comp) - f1;
560 amrex::Real f4 = (data_arr[n_of_data])(i+1,j+1,0,comp) - f1;
563 amrex::Real x1 = loc_arr(i ,j ,0,0);
564 amrex::Real x2 = loc_arr(i ,j+1,0,0) - x1;
565 amrex::Real x3 = loc_arr(i+1,j ,0,0) - x1;
566 amrex::Real x4 = loc_arr(i+1,j+1,0,0) - x1;
568 amrex::Real y1 = loc_arr(i ,j ,0,1);
569 amrex::Real y2 = loc_arr(i ,j+1,0,1) - y1;
570 amrex::Real y3 = loc_arr(i+1,j ,0,1) - y1;
571 amrex::Real y4 = loc_arr(i+1,j+1,0,1) - y1;
574 amrex::Real x = p.pos(0) - x1;
575 amrex::Real y = p.pos(1) - y1;
577 amrex::Real det = x2*x4*y3*(y2 - y4) - x3*x4*y2*(y3 - y4) - x2*x3*(y2 - y3)*y4;
579 amrex::Real b = ( f4*(x2*y2*y3 - x3*y2*y3) + (f2*(x3 - x4)*y3)*y4 + f3*(-(x2*y2*y4) + x4*y2*y4) ) / det;
580 amrex::Real c = ( -f2*x3*x4*y3 + f4*x2*x3*(-y2 + y3) + f3*x2*x4*(y2 - y4) + f2*x3*x4*y4 ) / det;
581 amrex::Real d = ( f2*x4*y3 + f4*(x3*y2 - x2*y3) - f2*x3*y4 + f3*(-(x4*y2) + x2*y4) ) / det;
583 amrex::Real
f = b*x + c*y + d*x*y;
584 val[comp-start_comp] = f1 +
f;
586 #elif (AMREX_SPACEDIM == 3)
594 Array1D<Real,0,7> x_vals;
595 Array1D<Real,0,7> y_vals;
596 Array1D<Real,0,7> z_vals;
599 amrex::Real x1, y1, z1, f1;
601 x1 = loc_arr(i ,j ,k ,0);
602 x_vals(1) = loc_arr(i+1,j ,k ,0) - x1;
603 x_vals(2) = loc_arr(i ,j+1,k ,0) - x1;
604 x_vals(3) = loc_arr(i+1,j+1,k ,0) - x1;
605 x_vals(4) = loc_arr(i ,j ,k+1,0) - x1;
606 x_vals(5) = loc_arr(i+1,j ,k+1,0) - x1;
607 x_vals(6) = loc_arr(i ,j+1,k+1,0) - x1;
608 x_vals(7) = loc_arr(i+1,j+1,k+1,0) - x1;
610 y1 = loc_arr(i ,j ,k ,1);
611 y_vals(1) = loc_arr(i+1,j ,k ,1) - y1;
612 y_vals(2) = loc_arr(i ,j+1,k ,1) - y1;
613 y_vals(3) = loc_arr(i+1,j+1,k ,1) - y1;
614 y_vals(4) = loc_arr(i ,j ,k+1,1) - y1;
615 y_vals(5) = loc_arr(i+1,j ,k+1,1) - y1;
616 y_vals(6) = loc_arr(i ,j+1,k+1,1) - y1;
617 y_vals(7) = loc_arr(i+1,j+1,k+1,1) - y1;
619 z1 = loc_arr(i ,j ,k ,2);
620 z_vals(1) = loc_arr(i+1,j ,k ,2) - z1;
621 z_vals(2) = loc_arr(i ,j+1,k ,2) - z1;
622 z_vals(3) = loc_arr(i+1,j+1,k ,2) - z1;
623 z_vals(4) = loc_arr(i ,j ,k+1,2) - z1;
624 z_vals(5) = loc_arr(i+1,j ,k+1,2) - z1;
625 z_vals(6) = loc_arr(i ,j+1,k+1,2) - z1;
626 z_vals(7) = loc_arr(i+1,j+1,k+1,2) - z1;
628 f1 = (data_arr[n_of_data])(i ,j ,k ,comp);
629 b(1) = (data_arr[n_of_data])(i+1,j ,k ,comp) - f1;
630 b(2) = (data_arr[n_of_data])(i ,j+1,k ,comp) - f1;
631 b(3) = (data_arr[n_of_data])(i+1,j+1,k ,comp) - f1;
632 b(4) = (data_arr[n_of_data])(i ,j ,k+1,comp) - f1;
633 b(5) = (data_arr[n_of_data])(i+1,j ,k+1,comp) - f1;
634 b(6) = (data_arr[n_of_data])(i ,j+1,k+1,comp) - f1;
635 b(7) = (data_arr[n_of_data])(i+1,j+1,k+1,comp) - f1;
637 Array2D<Real,1,7,1,7> a;
642 for (
int n=1; n<neq+1; n++) {
646 a(n,4) = x_vals(n)*y_vals(n);
647 a(n,5) = x_vals(n)*z_vals(n);
648 a(n,6) = z_vals(n)*y_vals(n);
649 a(n,7) = x_vals(n)*y_vals(n)*z_vals(n);
655 amrex::Real px, py, pz;
660 amrex::Real
f = b(1)*px + b(2)*py + b(3)*pz + b(4)*px*py + b(5)*px*pz + b(6)*py*pz + b(7)*px*py*pz;
661 val[comp-start_comp] =
static_cast<ParticleReal
>(
f + f1);
#define AMREX_ALWAYS_ASSERT_WITH_MESSAGE(EX, MSG)
Definition: AMReX_BLassert.H:49
#define AMREX_ASSERT(EX)
Definition: AMReX_BLassert.H:38
#define AMREX_FORCE_INLINE
Definition: AMReX_Extension.H:119
#define AMREX_GPU_HOST_DEVICE
Definition: AMReX_GpuQualifiers.H:20
#define AMREX_D_TERM(a, b, c)
Definition: AMReX_SPACE.H:129
#define AMREX_D_DECL(a, b, c)
Definition: AMReX_SPACE.H:104
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void cic_interpolate_mapped_z(const P &p, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &plo, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &dxi, const amrex::Array4< amrex::Real const > &data_arr, const amrex::Array4< amrex::Real const > &height_arr, amrex::ParticleReal *val, int M=AMREX_SPACEDIM)
Linearly interpolates the mesh data to the particle position from cell-centered data on a terrain-fit...
Definition: AMReX_TracerParticle_mod_K.H:167
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void linear_interpolate_to_particle_z(const P &p, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &plo, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &dxi, const Array4< amrex::Real const > *data_arr, const amrex::Array4< amrex::Real const > &height_arr, amrex::ParticleReal *val, const IntVect *is_nodal, int start_comp, int ncomp, int num_arrays)
Linearly interpolates the mesh data to the particle position from mesh data. This general form can ha...
Definition: AMReX_TracerParticle_mod_K.H:252
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void cic_interpolate_nd_mapped_z(const P &p, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &plo, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &dxi, const amrex::Array4< amrex::Real const > &data_arr, const amrex::Array4< amrex::Real const > &height_arr, amrex::ParticleReal *val, int M=AMREX_SPACEDIM)
Linearly interpolates the mesh data to the particle position from node-centered data....
Definition: AMReX_TracerParticle_mod_K.H:204
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mac_interpolate_mapped_z(const P &p, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &plo, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &dxi, amrex::GpuArray< amrex::Array4< amrex::Real const >, AMREX_SPACEDIM > const &data_arr, const amrex::Array4< amrex::Real const > &height_arr, amrex::ParticleReal *val)
Linearly interpolates the mesh data to the particle position from face-centered data on a terrain-fit...
Definition: AMReX_TracerParticle_mod_K.H:226
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void cic_interpolate_nd_mapped(const P &p, const amrex::Array4< amrex::Real const > &data_arr, const amrex::Array4< amrex::Real const > &loc_arr, amrex::ParticleReal *val, int M=AMREX_SPACEDIM)
Linearly interpolates the mesh data to the particle position from node-centered data on a general map...
Definition: AMReX_TracerParticle_mod_K.H:506
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void particle_interp_decomp(Array2D< Real, 1, 7, 1, 7 > &a, Array1D< int, 1, 7 > &ip, int neq)
Helper functions for 3D interpolation from values on nodes/vertices of arbitrary hexahedra to particl...
Definition: AMReX_TracerParticle_mod_K.H:422
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void particle_interp_solve(Array2D< Real, 1, 7, 1, 7 > &a, Array1D< Real, 1, 7 > &b, Array1D< int, 1, 7 > &ip, int neq)
Definition: AMReX_TracerParticle_mod_K.H:464
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void cic_interpolate_cc_mapped_z(const P &p, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &plo, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &dxi, const amrex::Array4< amrex::Real const > &data_arr, const amrex::Array4< amrex::Real const > &height_arr, amrex::ParticleReal *val, int M=AMREX_SPACEDIM)
Linearly interpolates the mesh data to the particle position from cell-centered data on a terrain-fit...
Definition: AMReX_TracerParticle_mod_K.H:183
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void linear_interpolate_to_particle_mapped(const P &p, const Array4< amrex::Real const > *data_arr, const amrex::Array4< amrex::Real const > &loc_arr, amrex::ParticleReal *val, const IntVect *is_nodal, int start_comp, int ncomp, int num_arrays)
Linearly interpolates the mesh data to the particle position on a general mapped grid....
Definition: AMReX_TracerParticle_mod_K.H:526
AMREX_GPU_HOST_DEVICE static constexpr AMREX_FORCE_INLINE IntVectND< dim > TheUnitVector() noexcept
This static member function returns a reference to a constant IntVectND object, all of whose dim argu...
Definition: AMReX_IntVect.H:682
AMREX_GPU_HOST_DEVICE static constexpr AMREX_FORCE_INLINE IntVectND< dim > TheZeroVector() noexcept
This static member function returns a reference to a constant IntVectND object, all of whose dim argu...
Definition: AMReX_IntVect.H:672
static int f(amrex::Real t, N_Vector y_data, N_Vector y_rhs, void *user_data)
Definition: AMReX_SundialsIntegrator.H:44
static constexpr int M
Definition: AMReX_OpenBC.H:13
static constexpr int P
Definition: AMReX_OpenBC.H:14
Definition: AMReX_Amr.cpp:49
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void cic_interpolate_nd(const P &p, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &plo, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &dxi, const amrex::Array4< amrex::Real const > &data_arr, amrex::ParticleReal *val, int M=AMREX_SPACEDIM)
Linearly interpolates the mesh data to the particle position from node-centered data.
Definition: AMReX_TracerParticle_mod_K.H:56
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void linear_interpolate_to_particle(const P &p, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &plo, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &dxi, const Array4< amrex::Real const > *data_arr, amrex::ParticleReal *val, const IntVect *is_nodal, int start_comp, int ncomp, int num_arrays)
Linearly interpolates the mesh data to the particle position from mesh data. This general form can ha...
Definition: AMReX_TracerParticle_mod_K.H:100
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void cic_interpolate(const P &p, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &plo, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &dxi, const amrex::Array4< amrex::Real const > &data_arr, amrex::ParticleReal *val, int M=AMREX_SPACEDIM)
Linearly interpolates the mesh data to the particle position from cell-centered data.
Definition: AMReX_TracerParticle_mod_K.H:27
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE T abs(const GpuComplex< T > &a_z) noexcept
Return the absolute value of a complex number.
Definition: AMReX_GpuComplex.H:356
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mac_interpolate(const P &p, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &plo, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &dxi, amrex::GpuArray< amrex::Array4< amrex::Real const >, AMREX_SPACEDIM > const &data_arr, amrex::ParticleReal *val)
Linearly interpolates the mesh data to the particle position from face-centered data....
Definition: AMReX_TracerParticle_mod_K.H:75
IntVectND< AMREX_SPACEDIM > IntVect
Definition: AMReX_BaseFwd.H:30
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void ignore_unused(const Ts &...)
This shuts up the compiler about unused variables.
Definition: AMReX.H:111
void Abort(const std::string &msg)
Print out message to cerr and exit via abort().
Definition: AMReX.cpp:225
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void cic_interpolate_cc(const P &p, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &plo, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &dxi, const amrex::Array4< amrex::Real const > &data_arr, amrex::ParticleReal *val, int M=AMREX_SPACEDIM)
Definition: AMReX_TracerParticle_mod_K.H:38
Definition: AMReX_Array4.H:61
Definition: AMReX_Array.H:34
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE const T * data() const noexcept
Definition: AMReX_Array.H:56