Block-Structured AMR Software Framework
AMReX_MLPoisson_1D_K.H
Go to the documentation of this file.
1 #ifndef AMREX_MLPOISSON_1D_K_H_
2 #define AMREX_MLPOISSON_1D_K_H_
3 #include <AMReX_Config.H>
4 
5 namespace amrex {
6 
7 template <typename T>
9 void mlpoisson_adotx (int i, Array4<T> const& y,
10  Array4<T const> const& x,
11  T dhx) noexcept
12 {
13  y(i,0,0) = dhx * (x(i-1,0,0) - T(2.0)*x(i,0,0) + x(i+1,0,0));
14 }
15 
16 template <typename T>
18 void mlpoisson_adotx_os (int i, Array4<T> const& y,
19  Array4<T const> const& x,
20  Array4<int const> const& osm,
21  T dhx) noexcept
22 {
23  if (osm(i,0,0) == 0) {
24  y(i,0,0) = T(0.0);
25  } else {
26  y(i,0,0) = dhx * (x(i-1,0,0) - T(2.0)*x(i,0,0) + x(i+1,0,0));
27  }
28 }
29 
30 template <typename T>
32 void mlpoisson_adotx_m (int i, Array4<T> const& y,
33  Array4<T const> const& x,
34  T dhx, T dx, T probxlo) noexcept
35 {
36  T rel = (probxlo + i *dx) * (probxlo + i *dx);
37  T rer = (probxlo +(i+1)*dx) * (probxlo +(i+1)*dx);
38  y(i,0,0) = dhx * (rel*x(i-1,0,0) - (rel+rer)*x(i,0,0) + rer*x(i+1,0,0));
39 }
40 
41 template <typename T>
43 void mlpoisson_flux_x (Box const& box, Array4<T> const& fx,
44  Array4<T const> const& sol, T dxinv) noexcept
45 {
46  const auto lo = amrex::lbound(box);
47  const auto hi = amrex::ubound(box);
48 
50  for (int i = lo.x; i <= hi.x; ++i) {
51  fx(i,0,0) = dxinv*(sol(i,0,0)-sol(i-1,0,0));
52  }
53 }
54 
55 template <typename T>
57 void mlpoisson_flux_x_m (Box const& box, Array4<T> const& fx,
58  Array4<T const> const& sol, T dxinv,
59  T dx, T probxlo) noexcept
60 {
61  const auto lo = amrex::lbound(box);
62  const auto hi = amrex::ubound(box);
63 
65  for (int i = lo.x; i <= hi.x; ++i) {
66  T re = (probxlo + i*dx) * (probxlo + i*dx);
67  fx(i,0,0) = dxinv*re*(sol(i,0,0)-sol(i-1,0,0));
68  }
69 }
70 
71 template <typename T>
73 void mlpoisson_flux_xface (Box const& box, Array4<T> const& fx,
74  Array4<T const> const& sol, T dxinv, int xlen) noexcept
75 {
76  const auto lo = amrex::lbound(box);
77 
78  int i = lo.x;
79  fx(i,0,0) = dxinv*(sol(i,0,0)-sol(i-1,0,0));
80  i += xlen;
81  fx(i,0,0) = dxinv*(sol(i,0,0)-sol(i-1,0,0));
82 }
83 
84 template <typename T>
86 void mlpoisson_flux_xface_m (Box const& box, Array4<T> const& fx,
87  Array4<T const> const& sol, T dxinv, int xlen,
88  T dx, T probxlo) noexcept
89 {
90  const auto lo = amrex::lbound(box);
91 
92  int i = lo.x;
93  T re = (probxlo + i*dx) * (probxlo + i*dx);
94  fx(i,0,0) = dxinv*re*(sol(i,0,0)-sol(i-1,0,0));
95  i += xlen;
96  re = (probxlo + i*dx) * (probxlo + i*dx);
97  fx(i,0,0) = dxinv*re*(sol(i,0,0)-sol(i-1,0,0));
98 }
99 
100 template <typename T>
102 void mlpoisson_gsrb (int i, int, int, Array4<T> const& phi, Array4<T const> const& rhs,
103  T dhx,
104  Array4<T const> const& f0, Array4<int const> const& m0,
105  Array4<T const> const& f1, Array4<int const> const& m1,
106  Box const& vbox, int redblack) noexcept
107 {
108  const auto vlo = amrex::lbound(vbox);
109  const auto vhi = amrex::ubound(vbox);
110 
111  T gamma = -dhx*T(2.0);
112 
113  if ((i+redblack)%2 == 0) {
114  T cf0 = (i == vlo.x && m0(vlo.x-1,0,0) > 0)
115  ? f0(vlo.x,0,0) : T(0.0);
116  T cf1 = (i == vhi.x && m1(vhi.x+1,0,0) > 0)
117  ? f1(vhi.x,0,0) : T(0.0);
118 
119  T g_m_d = gamma + dhx*(cf0+cf1);
120 
121  T res = rhs(i,0,0) - gamma*phi(i,0,0)
122  - dhx*(phi(i-1,0,0) + phi(i+1,0,0));
123 
124  phi(i,0,0) = phi(i,0,0) + res /g_m_d;
125  }
126 }
127 
128 template <typename T>
130 void mlpoisson_gsrb_os (int i, int, int, Array4<T> const& phi, Array4<T const> const& rhs,
131  Array4<int const> const& osm, T dhx,
132  Array4<T const> const& f0, Array4<int const> const& m0,
133  Array4<T const> const& f1, Array4<int const> const& m1,
134  Box const& vbox, int redblack) noexcept
135 {
136  const auto vlo = amrex::lbound(vbox);
137  const auto vhi = amrex::ubound(vbox);
138 
139  T gamma = -dhx*T(2.0);
140 
141  if ((i+redblack)%2 == 0) {
142  if (osm(i,0,0) == 0) {
143  phi(i,0,0) = T(0.0);
144  } else {
145  T cf0 = (i == vlo.x && m0(vlo.x-1,0,0) > 0)
146  ? f0(vlo.x,0,0) : T(0.0);
147  T cf1 = (i == vhi.x && m1(vhi.x+1,0,0) > 0)
148  ? f1(vhi.x,0,0) : T(0.0);
149 
150  T g_m_d = gamma + dhx*(cf0+cf1);
151 
152  T res = rhs(i,0,0) - gamma*phi(i,0,0)
153  - dhx*(phi(i-1,0,0) + phi(i+1,0,0));
154 
155  phi(i,0,0) = phi(i,0,0) + res /g_m_d;
156  }
157  }
158 }
159 
160 template <typename T>
162 void mlpoisson_gsrb_m (int i, int, int, Array4<T> const& phi, Array4<T const> const& rhs,
163  T dhx,
164  Array4<T const> const& f0, Array4<int const> const& m0,
165  Array4<T const> const& f1, Array4<int const> const& m1,
166  Box const& vbox, int redblack, T dx, T probxlo) noexcept
167 {
168  const auto vlo = amrex::lbound(vbox);
169  const auto vhi = amrex::ubound(vbox);
170 
171  if ((i+redblack)%2 == 0) {
172  T cf0 = (i == vlo.x && m0(vlo.x-1,0,0) > 0)
173  ? f0(vlo.x,0,0) : T(0.0);
174  T cf1 = (i == vhi.x && m1(vhi.x+1,0,0) > 0)
175  ? f1(vhi.x,0,0) : T(0.0);
176 
177  T rel = (probxlo + i *dx) * (probxlo + i *dx);
178  T rer = (probxlo +(i+1)*dx) * (probxlo +(i+1)*dx);
179 
180  T gamma = -dhx*(rel+rer);
181 
182  T g_m_d = gamma + dhx*(rel*cf0+rer*cf1);
183 
184  T res = rhs(i,0,0) - gamma*phi(i,0,0)
185  - dhx*(rel*phi(i-1,0,0) + rer*phi(i+1,0,0));
186 
187  phi(i,0,0) = phi(i,0,0) + res /g_m_d;
188  }
189 }
190 
191 template <typename T>
193 void mlpoisson_jacobi (int i, int, int, Array4<T> const& phi, Array4<T const> const& rhs,
194  Array4<T const> const& Ax, T dhx,
195  Array4<T const> const& f0, Array4<int const> const& m0,
196  Array4<T const> const& f1, Array4<int const> const& m1,
197  Box const& vbox) noexcept
198 {
199  const auto vlo = amrex::lbound(vbox);
200  const auto vhi = amrex::ubound(vbox);
201 
202  T gamma = -dhx*T(2.0);
203 
204  T cf0 = (i == vlo.x && m0(vlo.x-1,0,0) > 0)
205  ? f0(vlo.x,0,0) : T(0.0);
206  T cf1 = (i == vhi.x && m1(vhi.x+1,0,0) > 0)
207  ? f1(vhi.x,0,0) : T(0.0);
208 
209  T g_m_d = gamma + dhx*(cf0+cf1);
210 
211  phi(i,0,0) += T(2.0/3.0) * (rhs(i,0,0) - Ax(i,0,0)) / g_m_d;
212 }
213 
214 template <typename T>
216 void mlpoisson_jacobi_os (int i, int, int, Array4<T> const& phi, Array4<T const> const& rhs,
217  Array4<T const> const& Ax, Array4<int const> const& osm, T dhx,
218  Array4<T const> const& f0, Array4<int const> const& m0,
219  Array4<T const> const& f1, Array4<int const> const& m1,
220  Box const& vbox) noexcept
221 {
222  const auto vlo = amrex::lbound(vbox);
223  const auto vhi = amrex::ubound(vbox);
224 
225  if (osm(i,0,0) == 0) {
226  phi(i,0,0) = T(0.0);
227  } else {
228  T cf0 = (i == vlo.x && m0(vlo.x-1,0,0) > 0)
229  ? f0(vlo.x,0,0) : T(0.0);
230  T cf1 = (i == vhi.x && m1(vhi.x+1,0,0) > 0)
231  ? f1(vhi.x,0,0) : T(0.0);
232 
233  T gamma = -dhx*T(2.0);
234  T g_m_d = gamma + dhx*(cf0+cf1);
235 
236  phi(i,0,0) += T(2.0/3.0) * (rhs(i,0,0) - Ax(i,0,0)) / g_m_d;
237  }
238 }
239 
240 template <typename T>
242 void mlpoisson_jacobi_m (int i, int, int, Array4<T> const& phi, Array4<T const> const& rhs,
243  Array4<T const> const& Ax, T dhx,
244  Array4<T const> const& f0, Array4<int const> const& m0,
245  Array4<T const> const& f1, Array4<int const> const& m1,
246  Box const& vbox, T dx, T probxlo) noexcept
247 {
248  const auto vlo = amrex::lbound(vbox);
249  const auto vhi = amrex::ubound(vbox);
250 
251  T cf0 = (i == vlo.x && m0(vlo.x-1,0,0) > 0)
252  ? f0(vlo.x,0,0) : T(0.0);
253  T cf1 = (i == vhi.x && m1(vhi.x+1,0,0) > 0)
254  ? f1(vhi.x,0,0) : T(0.0);
255 
256  T rel = (probxlo + i *dx) * (probxlo + i *dx);
257  T rer = (probxlo +(i+1)*dx) * (probxlo +(i+1)*dx);
258 
259  T gamma = -dhx*(rel+rer);
260 
261  T g_m_d = gamma + dhx*(rel*cf0+rer*cf1);
262 
263  phi(i,0,0) += T(2.0/3.0) * (rhs(i,0,0) - Ax(i,0,0)) / g_m_d;
264 }
265 
266 template <typename T>
268 void mlpoisson_normalize (int i, int, int, Array4<T> const& x,
269  T dhx, T dx, T probxlo) noexcept
270 {
271  T rel = (probxlo + i *dx) * (probxlo + i *dx);
272  T rer = (probxlo +(i+1)*dx) * (probxlo +(i+1)*dx);
273  x(i,0,0) /= (-dhx*(rel+rer));
274 }
275 
276 }
277 
278 #endif
#define AMREX_PRAGMA_SIMD
Definition: AMReX_Extension.H:80
#define AMREX_FORCE_INLINE
Definition: AMReX_Extension.H:119
#define AMREX_GPU_HOST_DEVICE
Definition: AMReX_GpuQualifiers.H:20
Definition: AMReX_Amr.cpp:49
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_jacobi(int i, int, int, Array4< T > const &phi, Array4< T const > const &rhs, Array4< T const > const &Ax, T dhx, Array4< T const > const &f0, Array4< int const > const &m0, Array4< T const > const &f1, Array4< int const > const &m1, Box const &vbox) noexcept
Definition: AMReX_MLPoisson_1D_K.H:193
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_adotx_m(int i, Array4< T > const &y, Array4< T const > const &x, T dhx, T dx, T probxlo) noexcept
Definition: AMReX_MLPoisson_1D_K.H:32
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_gsrb_m(int i, int, int, Array4< T > const &phi, Array4< T const > const &rhs, T dhx, Array4< T const > const &f0, Array4< int const > const &m0, Array4< T const > const &f1, Array4< int const > const &m1, Box const &vbox, int redblack, T dx, T probxlo) noexcept
Definition: AMReX_MLPoisson_1D_K.H:162
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_flux_x(Box const &box, Array4< T > const &fx, Array4< T const > const &sol, T dxinv) noexcept
Definition: AMReX_MLPoisson_1D_K.H:43
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_adotx(int i, Array4< T > const &y, Array4< T const > const &x, T dhx) noexcept
Definition: AMReX_MLPoisson_1D_K.H:9
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_flux_xface_m(Box const &box, Array4< T > const &fx, Array4< T const > const &sol, T dxinv, int xlen, T dx, T probxlo) noexcept
Definition: AMReX_MLPoisson_1D_K.H:86
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_jacobi_os(int i, int, int, Array4< T > const &phi, Array4< T const > const &rhs, Array4< T const > const &Ax, Array4< int const > const &osm, T dhx, Array4< T const > const &f0, Array4< int const > const &m0, Array4< T const > const &f1, Array4< int const > const &m1, Box const &vbox) noexcept
Definition: AMReX_MLPoisson_1D_K.H:216
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_flux_x_m(Box const &box, Array4< T > const &fx, Array4< T const > const &sol, T dxinv, T dx, T probxlo) noexcept
Definition: AMReX_MLPoisson_1D_K.H:57
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_gsrb(int i, int, int, Array4< T > const &phi, Array4< T const > const &rhs, T dhx, Array4< T const > const &f0, Array4< int const > const &m0, Array4< T const > const &f1, Array4< int const > const &m1, Box const &vbox, int redblack) noexcept
Definition: AMReX_MLPoisson_1D_K.H:102
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 ubound(Array4< T > const &a) noexcept
Definition: AMReX_Array4.H:315
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 lbound(Array4< T > const &a) noexcept
Definition: AMReX_Array4.H:308
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_flux_xface(Box const &box, Array4< T > const &fx, Array4< T const > const &sol, T dxinv, int xlen) noexcept
Definition: AMReX_MLPoisson_1D_K.H:73
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_normalize(int i, int, int, Array4< T > const &x, T dhx, T dx, T probxlo) noexcept
Definition: AMReX_MLPoisson_1D_K.H:268
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_gsrb_os(int i, int, int, Array4< T > const &phi, Array4< T const > const &rhs, Array4< int const > const &osm, T dhx, Array4< T const > const &f0, Array4< int const > const &m0, Array4< T const > const &f1, Array4< int const > const &m1, Box const &vbox, int redblack) noexcept
Definition: AMReX_MLPoisson_1D_K.H:130
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_adotx_os(int i, Array4< T > const &y, Array4< T const > const &x, Array4< int const > const &osm, T dhx) noexcept
Definition: AMReX_MLPoisson_1D_K.H:18
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_jacobi_m(int i, int, int, Array4< T > const &phi, Array4< T const > const &rhs, Array4< T const > const &Ax, T dhx, Array4< T const > const &f0, Array4< int const > const &m0, Array4< T const > const &f1, Array4< int const > const &m1, Box const &vbox, T dx, T probxlo) noexcept
Definition: AMReX_MLPoisson_1D_K.H:242
Definition: AMReX_Array4.H:61