WarpX
sample.H
Go to the documentation of this file.
1 /* Copyright 2022 Edoardo Zoni, Remi Lehe, David Grote, Axel Huebl
2  *
3  * This file is part of ABLASTR.
4  *
5  * License: BSD-3-Clause-LBNL
6  */
7 #ifndef ABLASTR_COARSEN_SAMPLE_H_
8 #define ABLASTR_COARSEN_SAMPLE_H_
9 
10 
11 #include <AMReX_Array.H>
12 #include <AMReX_Array4.H>
13 #include <AMReX_Extension.H>
14 #include <AMReX_GpuQualifiers.H>
15 #include <AMReX_IntVect.H>
16 #include <AMReX_Math.H>
17 #include <AMReX_REAL.H>
18 
19 #include <AMReX_BaseFwd.H>
20 
21 #include <cstdlib>
22 
23 
29 {
49  amrex::Real Interp (
50  amrex::Array4<amrex::Real const> const& arr_src,
51  amrex::GpuArray<int,3> const& sf,
52  amrex::GpuArray<int,3> const& sc,
54  const int i,
55  const int j,
56  const int k,
57  const int comp
58  )
59  {
60  using namespace amrex::literals;
61 
62  // Indices of destination array (coarse)
63  const int ic[3] = { i, j, k };
64 
65  // Number of points and starting indices of source array (fine)
66  int np[3], idx_min[3];
67 
68  // Compute number of points
69  for ( int l = 0; l < 3; ++l ) {
70  if ( cr[l] == 1 ) { np[l] = 1+amrex::Math::abs(sf[l]-sc[l]); // no coarsening
71  } else { np[l] = 2-sf[l]; }
72  }
73 
74  // Compute starting indices of source array (fine)
75  for ( int l = 0; l < 3; ++l ) {
76  if ( cr[l] == 1 ) { idx_min[l] = ic[l]-sc[l]*(1-sf[l]); // no coarsening
77  } else { idx_min[l] = ic[l]*cr[l]+static_cast<int>(cr[l]/2)*(1-sc[l])-(1-sf[l]); }
78  }
79 
80  // Auxiliary integer variables
81  const int numx = np[0];
82  const int numy = np[1];
83  const int numz = np[2];
84  const int imin = idx_min[0];
85  const int jmin = idx_min[1];
86  const int kmin = idx_min[2];
87  amrex::Real const wx = 1.0_rt / static_cast<amrex::Real>(numx);
88  amrex::Real const wy = 1.0_rt / static_cast<amrex::Real>(numy);
89  amrex::Real const wz = 1.0_rt / static_cast<amrex::Real>(numz);
90 
91  // Interpolate over points computed above
92  amrex::Real c = 0.0_rt;
93  for (int kref = 0; kref < numz; ++kref) {
94  for (int jref = 0; jref < numy; ++jref) {
95  for (int iref = 0; iref < numx; ++iref) {
96  const int ii = imin+iref;
97  const int jj = jmin+jref;
98  const int kk = kmin+kref;
99  c += wx*wy*wz*arr_src(ii,jj,kk,comp);
100  }
101  }
102  }
103  return c;
104  }
105 
125  void Loop ( amrex::MultiFab& mf_dst,
126  const amrex::MultiFab& mf_src,
127  int dcomp,
128  int scomp,
129  int ncomp,
130  amrex::IntVect ngrow,
131  amrex::IntVect crse_ratio=amrex::IntVect(1) );
132 
152  void Coarsen ( amrex::MultiFab& mf_dst,
153  const amrex::MultiFab& mf_src,
154  int dcomp,
155  int scomp,
156  int ncomp,
157  int ngrow,
158  amrex::IntVect crse_ratio=amrex::IntVect(1) );
159  void Coarsen ( amrex::MultiFab& mf_dst,
160  const amrex::MultiFab& mf_src,
161  int dcomp,
162  int scomp,
163  int ncomp,
164  amrex::IntVect ngrowvect,
165  amrex::IntVect crse_ratio=amrex::IntVect(1) );
166 
167 } // namespace ablastr::coarsen::sample
168 
169 #endif // ABLASTR_COARSEN_SAMPLE_H_
#define AMREX_FORCE_INLINE
#define AMREX_GPU_DEVICE
Definition: sample.cpp:28
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real Interp(amrex::Array4< amrex::Real const > const &arr_src, amrex::GpuArray< int, 3 > const &sf, amrex::GpuArray< int, 3 > const &sc, amrex::GpuArray< int, 3 > const &cr, const int i, const int j, const int k, const int comp)
Interpolates the floating point data contained in the source Array4 arr_src, extracted from a fine Mu...
Definition: sample.H:49
void Coarsen(amrex::MultiFab &mf_dst, const amrex::MultiFab &mf_src, const int dcomp, const int scomp, const int ncomp, const int ngrow, const amrex::IntVect crse_ratio)
Stores in the coarsened MultiFab mf_dst the values obtained by interpolating the data contained in th...
Definition: sample.cpp:83
void Loop(amrex::MultiFab &mf_dst, const amrex::MultiFab &mf_src, const int dcomp, const int scomp, const int ncomp, const amrex::IntVect ngrowvect, const amrex::IntVect crse_ratio)
Loops over the boxes of the coarsened MultiFab mf_dst and fills them by interpolating the data contai...
Definition: sample.cpp:30
static constexpr auto c
vacuum speed of light [m/s]
Definition: constant.H:44
i
Definition: check_interp_points_and_weights.py:174
ii
Definition: check_interp_points_and_weights.py:148
jj
Definition: check_interp_points_and_weights.py:160
imin
Definition: check_interp_points_and_weights.py:124
cr
TODO Coarsening for IO: interpolation points and weights def coarsening_points_and_weights_for_IO( i,...
Definition: check_interp_points_and_weights.py:101