Block-Structured AMR Software Framework
AMReX_Habec_2D_K.H
Go to the documentation of this file.
1 #ifndef AMREX_HABEC_2D_H_
2 #define AMREX_HABEC_2D_H_
3 #include <AMReX_Config.H>
4 
5 #ifdef AMREX_USE_EB
6 #include <AMReX_EBMultiFabUtil.H>
7 #include <AMReX_MultiCutFab.H>
8 #include <AMReX_EBFabFactory.H>
9 #include <AMReX_MLLinOp_K.H>
10 #endif
11 
12 namespace amrex {
13 
14 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
15 void habec_mat (GpuArray<Real,2*AMREX_SPACEDIM+1>& sten, int i, int j, int k,
16  Dim3 const& boxlo, Dim3 const& boxhi,
17  Real sa, Array4<Real const> const& a,
18  Real sb, GpuArray<Real,AMREX_SPACEDIM> const& dx,
19  GpuArray<Array4<Real const>, AMREX_SPACEDIM> const& b,
20  GpuArray<int,AMREX_SPACEDIM*2> const& bctype,
21  GpuArray<Real,AMREX_SPACEDIM*2> const& bcl, int bho,
22  GpuArray<Array4<int const>, AMREX_SPACEDIM*2> const& msk,
23  Array4<Real> const& diaginv)
24 {
25  sten[1] = -(sb / (dx[0]*dx[0])) * b[0](i,j,k);
26  sten[2] = -(sb / (dx[0]*dx[0])) * b[0](i+1,j,k);
27  sten[3] = -(sb / (dx[1]*dx[1])) * b[1](i,j,k);
28  sten[4] = -(sb / (dx[1]*dx[1])) * b[1](i,j+1,k);
29  sten[0] = -(sten[1] + sten[2] + sten[3] + sten[4]);
30  if (sa != Real(0.0)) {
31  sten[0] += sa * a(i,j,k);
32  }
33 
34  // xlo
35  if (i == boxlo.x) {
36  int cdir = Orientation(Direction::x, Orientation::low);
37  if (msk[cdir](i-1,j,k) > 0) {
38  Real bf1, bf2;
39  detail::comp_bf(bf1, bf2, sb, dx[0], bctype[cdir], bcl[cdir], bho);
40  sten[0] += bf1 * b[0](i,j,k);
41  sten[1] = Real(0.0);
42  sten[2] += bf2 * b[0](i,j,k);
43  }
44  }
45 
46  // xhi
47  if (i == boxhi.x) {
48  int cdir = Orientation(Direction::x, Orientation::high);
49  if (msk[cdir](i+1,j,k) > 0) {
50  Real bf1, bf2;
51  detail::comp_bf(bf1, bf2, sb, dx[0], bctype[cdir], bcl[cdir], bho);
52  sten[0] += bf1 * b[0](i+1,j,k);
53  sten[1] += bf2 * b[0](i+1,j,k);
54  sten[2] = Real(0.0);
55  }
56  }
57 
58  // ylo
59  if (j == boxlo.y) {
60  int cdir = Orientation(Direction::y, Orientation::low);
61  if (msk[cdir](i,j-1,k) > 0) {
62  Real bf1, bf2;
63  detail::comp_bf(bf1, bf2, sb, dx[1], bctype[cdir], bcl[cdir], bho);
64  sten[0] += bf1 * b[1](i,j,k);
65  sten[3] = Real(0.0);
66  sten[4] += bf2 * b[1](i,j,k);
67  }
68  }
69 
70  // yhi
71  if (j == boxhi.y) {
72  int cdir = Orientation(Direction::y, Orientation::high);
73  if (msk[cdir](i,j+1,k) > 0) {
74  Real bf1, bf2;
75  detail::comp_bf(bf1, bf2, sb, dx[1], bctype[cdir], bcl[cdir], bho);
76  sten[0] += bf1 * b[1](i,j+1,k);
77  sten[3] += bf2 * b[1](i,j+1,k);
78  sten[4] = Real(0.0);
79  }
80  }
81 
82  diaginv(i,j,k) = Real(1.0) / sten[0];
83  sten[0] = Real(1.0);
84  for (int m = 1; m < 2*AMREX_SPACEDIM+1; ++m) {
85  sten[m] *= diaginv(i,j,k);
86  }
87 }
88 
89 template <typename Int>
90 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
91 void habec_ijmat (GpuArray<Real,2*AMREX_SPACEDIM+1>& sten, Array4<Int> const& ncols,
92  Array4<Real> const& diaginv, int i, int j, int k,
93  Array4<Int const> const& cell_id,
94  Real sa, Array4<Real const> const& a,
95  Real sb, GpuArray<Real,AMREX_SPACEDIM> const& dx,
96  GpuArray<Array4<Real const>, AMREX_SPACEDIM> const& b,
97  GpuArray<int,AMREX_SPACEDIM*2> const& bctype,
98  GpuArray<Real,AMREX_SPACEDIM*2> const& bcl, int bho,
99  Array4<int const> const& osm)
100 {
101  if (!osm || osm(i,j,k) != 0) {
102  sten[1] = -(sb / (dx[0]*dx[0])) * b[0](i,j,k);
103  sten[2] = -(sb / (dx[0]*dx[0])) * b[0](i+1,j,k);
104  sten[3] = -(sb / (dx[1]*dx[1])) * b[1](i,j,k);
105  sten[4] = -(sb / (dx[1]*dx[1])) * b[1](i,j+1,k);
106  sten[0] = -(sten[1] + sten[2] + sten[3] + sten[4]);
107  if (sa != Real(0.0)) {
108  sten[0] += sa * a(i,j,k);
109  }
110 
111  // xlo
112  if (cell_id(i-1,j,k) < 0) {
113  int cdir = Orientation(Direction::x, Orientation::low);
114  Real bf1, bf2;
115  detail::comp_bf(bf1, bf2, sb, dx[0], bctype[cdir], bcl[cdir], bho);
116  sten[0] += bf1 * b[0](i,j,k);
117  sten[1] = Real(0.0);
118  sten[2] += bf2 * b[0](i,j,k);
119  }
120 
121  // xhi
122  if (cell_id(i+1,j,k) < 0) {
123  int cdir = Orientation(Direction::x, Orientation::high);
124  Real bf1, bf2;
125  detail::comp_bf(bf1, bf2, sb, dx[0], bctype[cdir], bcl[cdir], bho);
126  sten[0] += bf1 * b[0](i+1,j,k);
127  sten[1] += bf2 * b[0](i+1,j,k);
128  sten[2] = Real(0.0);
129  }
130 
131  // ylo
132  if (cell_id(i,j-1,k) < 0) {
133  int cdir = Orientation(Direction::y, Orientation::low);
134  Real bf1, bf2;
135  detail::comp_bf(bf1, bf2, sb, dx[1], bctype[cdir], bcl[cdir], bho);
136  sten[0] += bf1 * b[1](i,j,k);
137  sten[3] = Real(0.0);
138  sten[4] += bf2 * b[1](i,j,k);
139  }
140 
141  // yhi
142  if (cell_id(i,j+1,k) < 0) {
143  int cdir = Orientation(Direction::y, Orientation::high);
144  Real bf1, bf2;
145  detail::comp_bf(bf1, bf2, sb, dx[1], bctype[cdir], bcl[cdir], bho);
146  sten[0] += bf1 * b[1](i,j+1,k);
147  sten[3] += bf2 * b[1](i,j+1,k);
148  sten[4] = Real(0.0);
149  }
150  } else {
151  sten[0] = Real(1.0);
152  for (int m = 1; m < 2*AMREX_SPACEDIM+1; ++m) {
153  sten[m] = Real(0.0);
154  }
155  }
156 
157  diaginv(i,j,k) = Real(1.0) / sten[0];
158  sten[0] = Real(1.0);
159  for (int m = 1; m < 2*AMREX_SPACEDIM+1; ++m) {
160  sten[m] *= diaginv(i,j,k);
161  }
162 
163  ncols(i,j,k) = 0;
164  for (int m = 0; m < 2*AMREX_SPACEDIM+1; ++m) {
165  ncols(i,j,k) += (sten[m] != Real(0.0));
166  }
167 }
168 
169 template <typename Int>
170 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
171 void habec_cols (GpuArray<Int,2*AMREX_SPACEDIM+1>& sten, int i, int j, int /*k*/,
172  Array4<Int const> const& cell_id)
173 {
174  sten[0] = cell_id(i ,j ,0);
175  sten[1] = cell_id(i-1,j ,0);
176  sten[2] = cell_id(i+1,j ,0);
177  sten[3] = cell_id(i ,j-1,0);
178  sten[4] = cell_id(i ,j+1,0);
179 }
180 
181 #ifdef AMREX_USE_EB
182 
183 template <typename Int>
184 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
185 void habec_cols_eb (GpuArray<Int,9>& sten, int i, int j, int /*k*/,
186  Array4<Int const> const& cell_id)
187 {
188  sten[0] = cell_id(i-1,j-1,0);
189  sten[1] = cell_id(i ,j-1,0);
190  sten[2] = cell_id(i+1,j-1,0);
191  sten[3] = cell_id(i-1,j ,0);
192  sten[4] = cell_id(i ,j ,0);
193  sten[5] = cell_id(i+1,j ,0);
194  sten[6] = cell_id(i-1,j+1,0);
195  sten[7] = cell_id(i ,j+1,0);
196  sten[8] = cell_id(i+1,j+1,0);
197 }
198 
199 template <typename Int>
200 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
201 void habec_ijmat_eb (GpuArray<Real,9>& sten, Array4<Int> const& ncols,
202  Array4<Real> const& diaginv, int i, int j, int k,
203  Array4<Int const> const& cell_id,
204  Real sa, Array4<Real const> const& a,
205  Real sb, GpuArray<Real,AMREX_SPACEDIM> const& dx,
206  GpuArray<Array4<Real const>, AMREX_SPACEDIM> const& b,
207  GpuArray<int,AMREX_SPACEDIM*2> const& bctype,
208  GpuArray<Real,AMREX_SPACEDIM*2> const& bcl, int bho,
209  Array4<const EBCellFlag> const& flag,
210  Array4<Real const> const& vfrc,
211  Array4<Real const> const& apx,
212  Array4<Real const> const& apy,
213  Array4<Real const> const& fcx,
214  Array4<Real const> const& fcy,
215  Array4<Real const> const& ba,
216  Array4<Real const> const& bcen,
217  Array4<Real const> const& beb)
218 {
219  for (int m = 0; m < 9; ++m) {
220  sten[m] = Real(0.0);
221  }
222 
223  auto& mat_tmp = reinterpret_cast<Array2D<Real,-1,1,-1,1>&>(sten);
224 
225  GpuArray<Real,AMREX_SPACEDIM> fac{sb/(dx[0]*dx[0]), sb/(dx[1]*dx[1])};
226 
227  if (flag(i,j,k).isRegular())
228  {
229  mat_tmp(0,0) = fac[0]*(b[0](i,j,k)+b[0](i+1,j,k))
230  + fac[1]*(b[1](i,j,k)+b[1](i,j+1,k));
231  mat_tmp(-1, 0) = -fac[0]*b[0](i,j,k);
232  mat_tmp( 1, 0) = -fac[0]*b[0](i+1,j,k);
233  mat_tmp( 0,-1) = -fac[1]*b[1](i,j,k);
234  mat_tmp( 0, 1) = -fac[1]*b[1](i,j+1,k);
235 
236  if (cell_id(i-1,j,k) < 0) {
237  int cdir = Orientation(Direction::x, Orientation::low);
238  Real bf1, bf2;
239  detail::comp_bf(bf1, bf2, sb, dx[0], bctype[cdir], bcl[cdir], bho);
240  mat_tmp(0,0) += bf1*b[0](i,j,k);
241  mat_tmp(-1,0) = Real(0.0);
242  mat_tmp(1,0) += bf2*b[0](i,j,k);
243  }
244 
245  if (cell_id(i+1,j,k) < 0) {
246  int cdir = Orientation(Direction::x, Orientation::high);
247  Real bf1, bf2;
248  detail::comp_bf(bf1, bf2, sb, dx[0], bctype[cdir], bcl[cdir], bho);
249  mat_tmp(0,0) += bf1*b[0](i+1,j,k);
250  mat_tmp(1,0) = Real(0.0);
251  mat_tmp(-1,0) += bf2*b[0](i+1,j,k);
252  }
253 
254  if (cell_id(i,j-1,k) < 0) {
255  int cdir = Orientation(Direction::y, Orientation::low);
256  Real bf1, bf2;
257  detail::comp_bf(bf1, bf2, sb, dx[1], bctype[cdir], bcl[cdir], bho);
258  mat_tmp(0,0) += bf1*b[1](i,j,k);
259  mat_tmp(0,-1) = Real(0.0);
260  mat_tmp(0,1) += bf2*b[1](i,j,k);
261  }
262 
263  if (cell_id(i,j+1,k) < 0) {
264  int cdir = Orientation(Direction::y, Orientation::high);
265  Real bf1, bf2;
266  detail::comp_bf(bf1, bf2, sb, dx[1], bctype[cdir], bcl[cdir], bho);
267  mat_tmp(0,0) += bf1*b[1](i,j+1,k);
268  mat_tmp(0,1) = Real(0.0);
269  mat_tmp(0,-1) += bf2*b[1](i,j+1,k);
270  }
271 
272  if (sa != Real(0.0)) {
273  mat_tmp(0,0) += sa*a(i,j,k);
274  }
275  }
276  else if (flag(i,j,k).isSingleValued())
277  {
278  int cdir = Orientation(Direction::x, Orientation::low);
279  Real area = apx(i,j,k);
280 
281  if (area > Real(0.0))
282  {
283  int joff, jj;
284  Real fracy;
285  Real f = fac[0];
286 
287  if (area != Real(1.0)) {
288  joff = static_cast<int>(std::copysign(Real(1.0), fcx(i,j,k)));
289  jj = j+joff;
290  if (cell_id(i-1,jj,k) < 0 && cell_id(i,jj,k) < 0) {
291  fracy = Real(0.0);
292  } else {
293  fracy = std::abs(fcx(i,j,k));
294  }
295  } else {
296  joff = 0;
297  jj = j;
298  fracy = Real(0.0);
299  }
300 
301  Real bf1, bf2, bflo;
302  detail::comp_bflo(bf1, bf2, bflo, sb, dx[0], bctype[cdir], bcl[cdir], bho);
303 
304  Real tmp = (Real(1.0)-fracy)*area*b[0](i,j,k);
305  if (cell_id(i-1,j,k) >= 0) {
306  mat_tmp( 0,0) += tmp*f;
307  mat_tmp(-1,0) -= tmp*f;
308  } else if (cell_id(i+1,j,k) < 0 || apx(i+1,j,k) == Real(0.0)) {
309  mat_tmp(0,0) += tmp*(f+bflo);
310  } else {
311  mat_tmp(0,0) += tmp*(f+bf1);
312  mat_tmp(1,0) += tmp* bf2;
313  }
314 
315  if (fracy > Real(0.0)) {
316  if (cell_id(i-1,jj,k) >= 0 && cell_id(i,jj,k) >= 0) {
317  mat_tmp(-1,joff) -= fracy*area*b[0](i,jj,k)*f;
318  mat_tmp( 0,joff) += fracy*area*b[0](i,jj,k)*f;
319  } else if (cell_id(i+1,jj,k) < 0 || apx(i+1,jj,k) == Real(0.0)) {
320  mat_tmp( 0,joff) += tmp*(f+bflo);
321  } else {
322  mat_tmp( 0,joff) += tmp*(f+bf1);
323  mat_tmp( 1,joff) += tmp* bf2;
324  }
325  }
326  }
327 
328  cdir = Orientation(Direction::x, Orientation::high);
329  area = apx(i+1,j,k);
330  if (area > Real(0.0))
331  {
332  int joff, jj;
333  Real fracy;
334  Real f = fac[0];
335 
336  if (area != Real(1.0)) {
337  joff = static_cast<int>(std::copysign(Real(1.0), fcx(i+1,j,k)));
338  jj = j+joff;
339  if (cell_id(i,jj,k) < 0 && cell_id(i+1,jj,k) < 0) {
340  fracy = Real(0.0);
341  } else {
342  fracy = std::abs(fcx(i+1,j,k));
343  }
344  } else {
345  joff = 0;
346  jj = j;
347  fracy = Real(0.0);
348  }
349 
350  Real bf1, bf2, bflo;
351  detail::comp_bflo(bf1, bf2, bflo, sb, dx[0], bctype[cdir], bcl[cdir], bho);
352 
353  Real tmp = (Real(1.0)-fracy)*area*b[0](i+1,j,k);
354  if (cell_id(i+1,j,k) >= 0) {
355  mat_tmp(0,0) += tmp*f;
356  mat_tmp(1,0) -= tmp*f;
357  } else if (cell_id(i-1,j,k) < 0 || apx(i,j,k) == Real(0.0)) {
358  mat_tmp(0,0) += tmp*(f+bflo);
359  } else {
360  mat_tmp( 0,0) += tmp*(f+bf1);
361  mat_tmp(-1,0) += tmp* bf2;
362  }
363 
364  if (fracy > Real(0.0)) {
365  if (cell_id(i,jj,k) >= 0 && cell_id(i+1,jj,k) >= 0) {
366  mat_tmp(0,joff) += fracy*area*b[0](i+1,jj,k)*f;
367  mat_tmp(1,joff) -= fracy*area*b[0](i+1,jj,k)*f;
368  } else if (cell_id(i-1,jj,k) < 0 || apx(i,jj,k) == Real(0.0)) {
369  mat_tmp(0,joff) += tmp*(f+bflo);
370  } else{
371  mat_tmp( 0,joff) += tmp*(f+bf1);
372  mat_tmp(-1,joff) += tmp* bf2;
373  }
374  }
375  }
376 
377  cdir = Orientation(Direction::y, Orientation::low);
378  area = apy(i,j,k);
379  if (area > Real(0.0))
380  {
381  int ioff, ii;
382  Real fracx ;
383  Real f = fac[1];
384 
385  if (area != Real(1.0)) {
386  ioff = static_cast<int>(std::copysign(Real(1.0), fcy(i,j,k)));
387  ii = i+ioff;
388  if (cell_id(ii,j-1,k) < 0 && cell_id(ii,j,k) < 0) {
389  fracx = Real(0.0);
390  } else {
391  fracx = std::abs(fcy(i,j,k));
392  }
393  } else {
394  ioff = 0;
395  ii = i;
396  fracx = Real(0.0);
397  }
398 
399  Real bf1, bf2, bflo;
400  detail::comp_bflo(bf1, bf2, bflo, sb, dx[1], bctype[cdir], bcl[cdir], bho);
401 
402  Real tmp = (Real(1.0)-fracx)*area*b[1](i,j,k);
403  if (cell_id(i,j-1,k) >= 0) {
404  mat_tmp(0, 0) += tmp*f;
405  mat_tmp(0,-1) -= tmp*f;
406  } else if (cell_id(i,j+1,k) < 0 || apy(i,j+1,k) == Real(0.0)) {
407  mat_tmp(0,0) += tmp*(f+bflo);
408  } else {
409  mat_tmp(0,0) += tmp*(f+bf1);
410  mat_tmp(0,1) += tmp* bf2;
411  }
412 
413  if (fracx > Real(0.0)) {
414  if (cell_id(ii,j-1,k) >= 0 && cell_id(ii,j,k) >= 0) {
415  mat_tmp(ioff,-1) -= fracx*area*b[1](ii,j,k)*f;
416  mat_tmp(ioff, 0) += fracx*area*b[1](ii,j,k)*f;
417  } else if (cell_id(ii,j+1,k) < 0 || apy(ii,j+1,k) == Real(0.0)) {
418  mat_tmp(ioff,0) += tmp*(f+bflo);
419  } else {
420  mat_tmp(ioff,0) += tmp*(f+bf1);
421  mat_tmp(ioff,1) += tmp* bf2;
422  }
423  }
424  }
425 
426  cdir = Orientation(Direction::y, Orientation::high);
427  area = apy(i,j+1,k);
428  if (area > Real(0.0))
429  {
430  int ioff, ii;
431  Real fracx;
432  Real f = fac[1];
433 
434  if (area != Real(1.0)) {
435  ioff = static_cast<int>(std::copysign(Real(1.0), fcy(i,j+1,k)));
436  ii = i+ioff;
437  if (cell_id(ii,j,k) < 0 && cell_id(ii,j+1,k) < 0) {
438  fracx = Real(0.0);
439  } else {
440  fracx = std::abs(fcy(i,j+1,k));
441  }
442  } else {
443  ioff = 0;
444  ii = i;
445  fracx = Real(0.0);
446  }
447 
448  Real bf1, bf2, bflo;
449  detail::comp_bflo(bf1, bf2, bflo, sb, dx[1], bctype[cdir], bcl[cdir], bho);
450 
451  Real tmp = (Real(1.0)-fracx)*area*b[1](i,j+1,k);
452  if (cell_id(i,j+1,k) >= 0) {
453  mat_tmp(0,0) += tmp*f;
454  mat_tmp(0,1) -= tmp*f;
455  } else if (cell_id(i,j-1,k) < 0 || apy(i,j,k) == Real(0.0)) {
456  mat_tmp(0,0) += tmp*(f+bflo);
457  } else {
458  mat_tmp(0, 0) += tmp*(f+bf1);
459  mat_tmp(0,-1) += tmp* bf2;
460  }
461 
462  if (fracx > Real(0.0)) {
463  if (cell_id(ii,j,k) >= 0 && cell_id(ii,j+1,k) >= 0) {
464  mat_tmp(ioff,0) += fracx*area*b[1](ii,j+1,k)*f;
465  mat_tmp(ioff,1) -= fracx*area*b[1](ii,j+1,k)*f;
466  } else if (cell_id(ii,j-1,k) < 0 || apy(ii,j,k) == Real(0.0)) {
467  mat_tmp(ioff,0) += tmp*(f+bflo);
468  } else {
469  mat_tmp(ioff, 0) += tmp*(f+bf1);
470  mat_tmp(ioff,-1) += tmp* bf2;
471  }
472  }
473  }
474 
475  if (beb) {
476  Real anorm = std::sqrt((apx(i,j,k) - apx(i+1,j,k))*(apx(i,j,k) - apx(i+1,j,k))
477  + (apy(i,j,k) - apy(i,j+1,k))*(apy(i,j,k) - apy(i,j+1,k)));
478 
479  Real anorminv = Real(1.0)/anorm;
480  Real anrmx = (apx(i,j,k) - apx(i+1,j,k))*anorminv;
481  Real anrmy = (apy(i,j,k) - apy(i,j+1,k))*anorminv;
482  Real sx = std::copysign(Real(1.0),anrmx);
483  Real sy = std::copysign(Real(1.0),anrmy);
484  Real bctx = bcen(i,j,k,0);
485  Real bcty = bcen(i,j,k,1);
486  Real dg = get_dx_eb(vfrc(i,j,k)) / amrex::max(std::abs(anrmx), std::abs(anrmy));
487  Real gx = bctx - dg*anrmx;
488  Real gy = bcty - dg*anrmy;
489  int ioff = -static_cast<int>(sx);
490  int joff = -static_cast<int>(sy);
491  Array1D<Real,0,3> phig1{Real(1.0) + gx*sx + gy*sy + gx*gy*sx*sy,
492  - gx*sx - gx*gy*sx*sy,
493  - gy*sy - gx*gy*sx*sy,
494  + gx*gy*sx*sy};
495 
496  Array1D<Real,0,3> feb;
497  for (int ii=0; ii<4; ii++){
498  feb(ii) = -phig1(ii) * (ba(i,j,k) * beb(i,j,k) / dg);
499  }
500  mat_tmp(0 , 0 ) -= feb(0)*fac[0];
501  mat_tmp(ioff, 0 ) -= feb(1)*fac[0];
502  mat_tmp(0 ,joff) -= feb(2)*fac[0];
503  mat_tmp(ioff,joff) -= feb(3)*fac[0];
504  }
505 
506  for (int jj=-1; jj<=1; jj++) {
507  for (int ii=-1; ii<=1; ii++) {
508  mat_tmp(ii,jj) *= Real(1.0)/vfrc(i,j,k);
509  }}
510 
511  if (sa != Real(0.0)) {
512  mat_tmp(0,0) += sa*a(i,j,k);
513  }
514  }
515  else
516  {
517  for (int jj=-1; jj<=1; jj++) {
518  for (int ii=-1; ii<=1; ii++) {
519  mat_tmp(ii,jj) = Real(0.0);
520  }}
521  mat_tmp(0,0) = Real(1.0);
522  }
523 
524  diaginv(i,j,k) = Real(1.0) / mat_tmp(0,0);
525  for (int jj=-1; jj<=1; jj++) {
526  for (int ii=-1; ii<=1; ii++) {
527  mat_tmp(ii,jj) *= diaginv(i,j,k);
528  }}
529  mat_tmp(0,0) = Real(1.0);
530 
531  ncols(i,j,k) = 0;
532  for (int m = 0; m < 9; ++m) {
533  ncols(i,j,k) += (sten[m] != Real(0.0));
534  }
535 }
536 
537 #endif
538 
539 }
540 #endif
AMReX_EBFabFactory.H
AMReX_EBMultiFabUtil.H
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_MLLinOp_K.H
AMReX_MultiCutFab.H
amrex::Orientation
Encapsulation of the Orientation of the Faces of a Box.
Definition: AMReX_Orientation.H:29
amrex::Orientation::low
@ low
Definition: AMReX_Orientation.H:34
amrex::Orientation::high
@ high
Definition: AMReX_Orientation.H:34
amrex::SundialsUserFun::f
static int f(amrex::Real t, N_Vector y_data, N_Vector y_rhs, void *user_data)
Definition: AMReX_SundialsIntegrator.H:44
amrex::detail::comp_bf
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void comp_bf(Real &bf1, Real &bf2, Real sb, Real h, int bct, Real bcl, int bho)
Definition: AMReX_Habec_K.H:11
amrex::detail::comp_bflo
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void comp_bflo(Real &bf1, Real &bf2, Real &bflo, Real sb, Real h, int bct, Real bcl, int bho)
Definition: AMReX_Habec_K.H:34
amrex
Definition: AMReX_Amr.cpp:49
amrex::habec_ijmat
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void habec_ijmat(GpuArray< Real, 2 *AMREX_SPACEDIM+1 > &sten, Array4< Int > const &ncols, Array4< Real > const &diaginv, int i, int j, int k, Array4< Int const > const &cell_id, Real sa, Array4< Real const > const &a, Real sb, GpuArray< Real, AMREX_SPACEDIM > const &dx, GpuArray< Array4< Real const >, AMREX_SPACEDIM > const &b, GpuArray< int, AMREX_SPACEDIM *2 > const &bctype, GpuArray< Real, AMREX_SPACEDIM *2 > const &bcl, int bho, Array4< int const > const &osm)
Definition: AMReX_Habec_2D_K.H:91
amrex::max
AMREX_GPU_HOST_DEVICE constexpr AMREX_FORCE_INLINE const T & max(const T &a, const T &b) noexcept
Definition: AMReX_Algorithm.H:35
amrex::abs
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::habec_cols
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void habec_cols(GpuArray< Int, 2 *AMREX_SPACEDIM+1 > &sten, int i, int j, int, Array4< Int const > const &cell_id)
Definition: AMReX_Habec_2D_K.H:171
amrex::habec_mat
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void habec_mat(GpuArray< Real, 2 *AMREX_SPACEDIM+1 > &sten, int i, int j, int k, Dim3 const &boxlo, Dim3 const &boxhi, Real sa, Array4< Real const > const &a, Real sb, GpuArray< Real, AMREX_SPACEDIM > const &dx, GpuArray< Array4< Real const >, AMREX_SPACEDIM > const &b, GpuArray< int, AMREX_SPACEDIM *2 > const &bctype, GpuArray< Real, AMREX_SPACEDIM *2 > const &bcl, int bho, GpuArray< Array4< int const >, AMREX_SPACEDIM *2 > const &msk, Array4< Real > const &diaginv)
Definition: AMReX_Habec_2D_K.H:15
amrex::sqrt
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE GpuComplex< T > sqrt(const GpuComplex< T > &a_z) noexcept
Return the square root of a complex number.
Definition: AMReX_GpuComplex.H:373
amrex::Array4
Definition: AMReX_Array4.H:61
amrex::Dim3
Definition: AMReX_Dim3.H:12
amrex::Dim3::x
int x
Definition: AMReX_Dim3.H:12
amrex::Dim3::y
int y
Definition: AMReX_Dim3.H:12
amrex::GpuArray
Definition: AMReX_Array.H:34