Block-Structured AMR Software Framework
AMReX_MLABecLap_1D_K.H
Go to the documentation of this file.
1 #ifndef AMREX_MLABECLAP_1D_K_H_
2 #define AMREX_MLABECLAP_1D_K_H_
3 #include <AMReX_Config.H>
4 
5 namespace amrex {
6 
7 template <typename T>
8 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
9 void mlabeclap_adotx (int i, int, int, int n, Array4<T> const& y,
10  Array4<T const> const& x,
11  Array4<T const> const& a,
12  Array4<T const> const& bX,
13  GpuArray<T,AMREX_SPACEDIM> const& dxinv,
14  T alpha, T beta) noexcept
15 {
16  const T dhx = beta*dxinv[0]*dxinv[0];
17  y(i,0,0,n) = alpha*a(i,0,0)*x(i,0,0,n)
18  - dhx * (bX(i+1,0,0,n)*(x(i+1,0,0,n) - x(i ,0,0,n))
19  - bX(i ,0,0,n)*(x(i ,0,0,n) - x(i-1,0,0,n)));
20 }
21 
22 template <typename T>
23 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
24 void mlabeclap_adotx_os (int i, int, int, int n, Array4<T> const& y,
25  Array4<T const> const& x,
26  Array4<T const> const& a,
27  Array4<T const> const& bX,
28  Array4<int const> const& osm,
29  GpuArray<T,AMREX_SPACEDIM> const& dxinv,
30  T alpha, T beta) noexcept
31 {
32  const T dhx = beta*dxinv[0]*dxinv[0];
33  if (osm(i,0,0) == 0) {
34  y(i,0,0,n) = T(0.0);
35  } else {
36  y(i,0,0,n) = alpha*a(i,0,0)*x(i,0,0,n)
37  - dhx * (bX(i+1,0,0,n)*(x(i+1,0,0,n) - x(i ,0,0,n))
38  - bX(i ,0,0,n)*(x(i ,0,0,n) - x(i-1,0,0,n)));
39  }
40 }
41 
42 template <typename T>
43 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
44 void mlabeclap_normalize (int i, int, int, int n, Array4<T> const& x,
45  Array4<T const> const& a,
46  Array4<T const> const& bX,
47  GpuArray<T,AMREX_SPACEDIM> const& dxinv,
48  T alpha, T beta) noexcept
49 {
50  const T dhx = beta*dxinv[0]*dxinv[0];
51  x(i,0,0,n) /= alpha*a(i,0,0) + dhx*(bX(i,0,0,n)+bX(i+1,0,0,n));
52 }
53 
54 template <typename T>
55 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
56 void mlabeclap_flux_x (Box const& box, Array4<T> const& fx, Array4<T const> const& sol,
57  Array4<T const> const& bx, T fac, int ncomp) noexcept
58 {
59  const auto lo = amrex::lbound(box);
60  const auto hi = amrex::ubound(box);
61 
62  for (int n = 0; n < ncomp; ++n) {
63  AMREX_PRAGMA_SIMD
64  for (int i = lo.x; i <= hi.x; ++i) {
65  fx(i,0,0,n) = -fac*bx(i,0,0,n)*(sol(i,0,0,n)-sol(i-1,0,0,n));
66  }
67  }
68 }
69 
70 template <typename T>
71 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
72 void mlabeclap_flux_xface (Box const& box, Array4<T> const& fx, Array4<T const> const& sol,
73  Array4<T const> const& bx, T fac, int xlen, int ncomp) noexcept
74 {
75  const auto lo = amrex::lbound(box);
76 
77  for (int n = 0; n < ncomp; ++n) {
78  int i = lo.x;
79  fx(i,0,0,n) = -fac*bx(i,0,0,n)*(sol(i,0,0,n)-sol(i-1,0,0,n));
80  i += xlen;
81  fx(i,0,0,n) = -fac*bx(i,0,0,n)*(sol(i,0,0,n)-sol(i-1,0,0,n));
82  }
83 }
84 
85 template <typename T>
86 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
87 void abec_gsrb (int i, int, int, int n, Array4<T> const& phi, Array4<T const> const& rhs,
88  T alpha, Array4<T const> const& a,
89  T dhx,
90  Array4<T const> const& bX,
91  Array4<int const> const& m0,
92  Array4<int const> const& m1,
93  Array4<T const> const& f0,
94  Array4<T const> const& f1,
95  Box const& vbox, int redblack) noexcept
96 {
97  if ((i+redblack)%2 == 0) {
98  const auto vlo = amrex::lbound(vbox);
99  const auto vhi = amrex::ubound(vbox);
100 
101  T cf0 = (i == vlo.x && m0(vlo.x-1,0,0) > 0)
102  ? f0(vlo.x,0,0,n) : T(0.0);
103  T cf1 = (i == vhi.x && m1(vhi.x+1,0,0) > 0)
104  ? f1(vhi.x,0,0,n) : T(0.0);
105 
106  T delta = dhx*(bX(i,0,0,n)*cf0 + bX(i+1,0,0,n)*cf1);
107 
108  T gamma = alpha*a(i,0,0)
109  + dhx*( bX(i,0,0,n) + bX(i+1,0,0,n) );
110 
111  T rho = dhx*(bX(i ,0 ,0,n)*phi(i-1,0 ,0,n)
112  + bX(i+1,0 ,0,n)*phi(i+1,0 ,0,n));
113 
114  phi(i,0,0,n) = (rhs(i,0,0,n) + rho - phi(i,0,0,n)*delta)
115  / (gamma - delta);
116  }
117 }
118 
119 template <typename T>
120 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
121 void abec_gsrb_os (int i, int, int, int n, Array4<T> const& phi, Array4<T const> const& rhs,
122  T alpha, Array4<T const> const& a,
123  T dhx,
124  Array4<T const> const& bX,
125  Array4<int const> const& m0,
126  Array4<int const> const& m1,
127  Array4<T const> const& f0,
128  Array4<T const> const& f1,
129  Array4<int const> const& osm,
130  Box const& vbox, int redblack) noexcept
131 {
132  if ((i+redblack)%2 == 0) {
133  if (osm(i,0,0) == 0) {
134  phi(i,0,0) = T(0.0);
135  } else {
136  const auto vlo = amrex::lbound(vbox);
137  const auto vhi = amrex::ubound(vbox);
138 
139  T cf0 = (i == vlo.x && m0(vlo.x-1,0,0) > 0)
140  ? f0(vlo.x,0,0,n) : T(0.0);
141  T cf1 = (i == vhi.x && m1(vhi.x+1,0,0) > 0)
142  ? f1(vhi.x,0,0,n) : T(0.0);
143 
144  T delta = dhx*(bX(i,0,0,n)*cf0 + bX(i+1,0,0,n)*cf1);
145 
146  T gamma = alpha*a(i,0,0)
147  + dhx*( bX(i,0,0,n) + bX(i+1,0,0,n) );
148 
149  T rho = dhx*(bX(i ,0 ,0,n)*phi(i-1,0 ,0,n)
150  + bX(i+1,0 ,0,n)*phi(i+1,0 ,0,n));
151 
152  phi(i,0,0,n) = (rhs(i,0,0,n) + rho - phi(i,0,0,n)*delta)
153  / (gamma - delta);
154  }
155  }
156 }
157 
158 template <typename T>
159 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
160 void abec_jacobi (int i, int, int, int n, Array4<T> const& phi,
161  Array4<T const> const& rhs, Array4<T const> const& Ax,
162  T alpha, Array4<T const> const& a,
163  T dhx,
164  Array4<T const> const& bX,
165  Array4<int const> const& m0,
166  Array4<int const> const& m1,
167  Array4<T const> const& f0,
168  Array4<T const> const& f1,
169  Box const& vbox) noexcept
170 {
171  const auto vlo = amrex::lbound(vbox);
172  const auto vhi = amrex::ubound(vbox);
173 
174  T cf0 = (i == vlo.x && m0(vlo.x-1,0,0) > 0)
175  ? f0(vlo.x,0,0,n) : T(0.0);
176  T cf1 = (i == vhi.x && m1(vhi.x+1,0,0) > 0)
177  ? f1(vhi.x,0,0,n) : T(0.0);
178 
179  T delta = dhx*(bX(i,0,0,n)*cf0 + bX(i+1,0,0,n)*cf1);
180 
181  T gamma = alpha*a(i,0,0)
182  + dhx*( bX(i,0,0,n) + bX(i+1,0,0,n) );
183 
184  phi(i,0,0,n) += T(2.0/3.0) * (rhs(i,0,0,n) - Ax(i,0,0,n)) / (gamma - delta);
185 }
186 
187 template <typename T>
188 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
189 void abec_jacobi_os (int i, int, int, int n, Array4<T> const& phi,
190  Array4<T const> const& rhs, Array4<T const> const& Ax,
191  T alpha, Array4<T const> const& a,
192  T dhx,
193  Array4<T const> const& bX,
194  Array4<int const> const& m0,
195  Array4<int const> const& m1,
196  Array4<T const> const& f0,
197  Array4<T const> const& f1,
198  Array4<int const> const& osm,
199  Box const& vbox) noexcept
200 {
201  if (osm(i,0,0) == 0) {
202  phi(i,0,0) = T(0.0);
203  } else {
204  const auto vlo = amrex::lbound(vbox);
205  const auto vhi = amrex::ubound(vbox);
206 
207  T cf0 = (i == vlo.x && m0(vlo.x-1,0,0) > 0)
208  ? f0(vlo.x,0,0,n) : T(0.0);
209  T cf1 = (i == vhi.x && m1(vhi.x+1,0,0) > 0)
210  ? f1(vhi.x,0,0,n) : T(0.0);
211 
212  T delta = dhx*(bX(i,0,0,n)*cf0 + bX(i+1,0,0,n)*cf1);
213 
214  T gamma = alpha*a(i,0,0)
215  + dhx*( bX(i,0,0,n) + bX(i+1,0,0,n) );
216 
217  phi(i,0,0,n) += T(2.0/3.0) * (rhs(i,0,0,n) - Ax(i,0,0,n)) / (gamma - delta);
218  }
219 }
220 
221 template <typename T>
222 AMREX_FORCE_INLINE
223 void abec_gsrb_with_line_solve (
224  Box const& /*box*/, Array4<T> const& /*phi*/, Array4<T const> const& /*rhs*/,
225  T /*alpha*/, Array4<T const> const& /*a*/,
226  T /*dhx*/,
227  Array4<T const> const& /*bX*/,
228  Array4<int const> const& /*m0*/,
229  Array4<int const> const& /*m1*/,
230  Array4<T const> const& /*f0*/,
231  Array4<T const> const& /*f1*/,
232  Box const& /*vbox*/, int /*redblack*/, int /*nc*/) noexcept
233 {
234  amrex::Abort("abec_gsrb_with_line_solve not implemented in 1D");
235 }
236 
237 template <typename T>
238 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
239 void overset_rescale_bcoef_x (Box const& box, Array4<T> const& bX, Array4<int const> const& osm,
240  int ncomp, T osfac) noexcept
241 {
242  const auto lo = amrex::lbound(box);
243  const auto hi = amrex::ubound(box);
244  for (int n = 0; n < ncomp; ++n) {
245  for (int i = lo.x; i <= hi.x; ++i) {
246  if ((osm(i-1,0,0)+osm(i,0,0)) == 1) {
247  bX(i,0,0,n) *= osfac;
248  }
249  }
250  }
251 }
252 
253 }
254 #endif
AMREX_PRAGMA_SIMD
#define AMREX_PRAGMA_SIMD
Definition: AMReX_Extension.H:80
AMREX_FORCE_INLINE
#define AMREX_FORCE_INLINE
Definition: AMReX_Extension.H:119
AMREX_GPU_HOST_DEVICE
#define AMREX_GPU_HOST_DEVICE
Definition: AMReX_GpuQualifiers.H:20
amrex
Definition: AMReX_Amr.cpp:49
amrex::overset_rescale_bcoef_x
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void overset_rescale_bcoef_x(Box const &box, Array4< T > const &bX, Array4< int const > const &osm, int ncomp, T osfac) noexcept
Definition: AMReX_MLABecLap_1D_K.H:239
amrex::mlabeclap_adotx_os
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlabeclap_adotx_os(int i, int, int, int n, Array4< T > const &y, Array4< T const > const &x, Array4< T const > const &a, Array4< T const > const &bX, Array4< int const > const &osm, GpuArray< T, AMREX_SPACEDIM > const &dxinv, T alpha, T beta) noexcept
Definition: AMReX_MLABecLap_1D_K.H:24
amrex::mlabeclap_flux_x
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlabeclap_flux_x(Box const &box, Array4< T > const &fx, Array4< T const > const &sol, Array4< T const > const &bx, T fac, int ncomp) noexcept
Definition: AMReX_MLABecLap_1D_K.H:56
amrex::abec_gsrb_os
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void abec_gsrb_os(int i, int, int, int n, Array4< T > const &phi, Array4< T const > const &rhs, T alpha, Array4< T const > const &a, T dhx, Array4< T const > const &bX, Array4< int const > const &m0, Array4< int const > const &m1, Array4< T const > const &f0, Array4< T const > const &f1, Array4< int const > const &osm, Box const &vbox, int redblack) noexcept
Definition: AMReX_MLABecLap_1D_K.H:121
amrex::ubound
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 ubound(Array4< T > const &a) noexcept
Definition: AMReX_Array4.H:315
amrex::lbound
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 lbound(Array4< T > const &a) noexcept
Definition: AMReX_Array4.H:308
amrex::abec_jacobi
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void abec_jacobi(int i, int, int, int n, Array4< T > const &phi, Array4< T const > const &rhs, Array4< T const > const &Ax, T alpha, Array4< T const > const &a, T dhx, Array4< T const > const &bX, Array4< int const > const &m0, Array4< int const > const &m1, Array4< T const > const &f0, Array4< T const > const &f1, Box const &vbox) noexcept
Definition: AMReX_MLABecLap_1D_K.H:160
amrex::abec_gsrb
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void abec_gsrb(int i, int, int, int n, Array4< T > const &phi, Array4< T const > const &rhs, T alpha, Array4< T const > const &a, T dhx, Array4< T const > const &bX, Array4< int const > const &m0, Array4< int const > const &m1, Array4< T const > const &f0, Array4< T const > const &f1, Box const &vbox, int redblack) noexcept
Definition: AMReX_MLABecLap_1D_K.H:87
amrex::abec_gsrb_with_line_solve
AMREX_FORCE_INLINE void abec_gsrb_with_line_solve(Box const &, Array4< T > const &, Array4< T const > const &, T, Array4< T const > const &, T, Array4< T const > const &, Array4< int const > const &, Array4< int const > const &, Array4< T const > const &, Array4< T const > const &, Box const &, int, int) noexcept
Definition: AMReX_MLABecLap_1D_K.H:223
amrex::mlabeclap_normalize
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlabeclap_normalize(int i, int, int, int n, Array4< T > const &x, Array4< T const > const &a, Array4< T const > const &bX, GpuArray< T, AMREX_SPACEDIM > const &dxinv, T alpha, T beta) noexcept
Definition: AMReX_MLABecLap_1D_K.H:44
amrex::Abort
void Abort(const std::string &msg)
Print out message to cerr and exit via abort().
Definition: AMReX.cpp:221
amrex::abec_jacobi_os
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void abec_jacobi_os(int i, int, int, int n, Array4< T > const &phi, Array4< T const > const &rhs, Array4< T const > const &Ax, T alpha, Array4< T const > const &a, T dhx, Array4< T const > const &bX, Array4< int const > const &m0, Array4< int const > const &m1, Array4< T const > const &f0, Array4< T const > const &f1, Array4< int const > const &osm, Box const &vbox) noexcept
Definition: AMReX_MLABecLap_1D_K.H:189
amrex::mlabeclap_flux_xface
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlabeclap_flux_xface(Box const &box, Array4< T > const &fx, Array4< T const > const &sol, Array4< T const > const &bx, T fac, int xlen, int ncomp) noexcept
Definition: AMReX_MLABecLap_1D_K.H:72
amrex::mlabeclap_adotx
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlabeclap_adotx(int i, int, int, int n, Array4< T > const &y, Array4< T const > const &x, Array4< T const > const &a, Array4< T const > const &bX, GpuArray< T, AMREX_SPACEDIM > const &dxinv, T alpha, T beta) noexcept
Definition: AMReX_MLABecLap_1D_K.H:9
amrex::Array4
Definition: AMReX_Array4.H:61
amrex::GpuArray
Definition: AMReX_Array.H:33