/home/runner/work/kynema-sgf/kynema-sgf/src/fvm/nonLinearSum.H Source File

Kynema-SGF API: /home/runner/work/kynema-sgf/kynema-sgf/src/fvm/nonLinearSum.H Source File
Kynema-SGF API v0.1.0
CFD solver for wind plant simulations
Loading...
Searching...
No Matches
nonLinearSum.H
Go to the documentation of this file.
1#ifndef NONLINEARSUM_H
2#define NONLINEARSUM_H
3
4#include "src/fvm/fvm_utils.H"
5
6namespace kynema_sgf::fvm {
7
12template <typename FTypeIn, typename FTypeOut>
14{
15 NonLinearSum(FTypeOut& strphi, const FTypeIn& phi)
16 : m_strphi(strphi), m_phi(phi)
17 {
18 AMREX_ALWAYS_ASSERT(AMREX_SPACEDIM == m_phi.num_comp());
19 AMREX_ALWAYS_ASSERT(m_phi.num_grow() > amrex::IntVect(0));
20 }
21
22 template <typename Stencil>
23 void apply(const int lev) const
24 {
25 const auto& geom = m_phi.repo().mesh().Geom(lev);
26 const auto& idx = geom.InvCellSizeArray();
27 const auto dlo = geom.Domain().smallEnd();
28 const auto dhi = geom.Domain().bigEnd();
29
30 auto& out_mf = m_strphi(lev);
31 auto const& in_arrs = m_phi(lev).const_arrays();
32 auto const& out_arrs = out_mf.arrays();
33
34 amrex::ParallelFor(
35 out_mf, amrex::IntVect(0),
36 [=] AMREX_GPU_DEVICE(int nbx, int i, int j, int k) {
37 if (!Stencil::applies_to(i, j, k, dlo, dhi)) {
38 return;
39 }
40 auto const& phi = in_arrs[nbx];
41 auto const& strphi = out_arrs[nbx];
42 const amrex::Real ux = (Stencil::c00 * phi(i + 1, j, k, 0) +
43 Stencil::c01 * phi(i, j, k, 0) +
44 Stencil::c02 * phi(i - 1, j, k, 0)) *
45 idx[0];
46 const amrex::Real vx = (Stencil::c00 * phi(i + 1, j, k, 1) +
47 Stencil::c01 * phi(i, j, k, 1) +
48 Stencil::c02 * phi(i - 1, j, k, 1)) *
49 idx[0];
50 const amrex::Real wx = (Stencil::c00 * phi(i + 1, j, k, 2) +
51 Stencil::c01 * phi(i, j, k, 2) +
52 Stencil::c02 * phi(i - 1, j, k, 2)) *
53 idx[0];
54 const amrex::Real uy = (Stencil::c10 * phi(i, j + 1, k, 0) +
55 Stencil::c11 * phi(i, j, k, 0) +
56 Stencil::c12 * phi(i, j - 1, k, 0)) *
57 idx[1];
58 const amrex::Real vy = (Stencil::c10 * phi(i, j + 1, k, 1) +
59 Stencil::c11 * phi(i, j, k, 1) +
60 Stencil::c12 * phi(i, j - 1, k, 1)) *
61 idx[1];
62 const amrex::Real wy = (Stencil::c10 * phi(i, j + 1, k, 2) +
63 Stencil::c11 * phi(i, j, k, 2) +
64 Stencil::c12 * phi(i, j - 1, k, 2)) *
65 idx[1];
66 const amrex::Real uz = (Stencil::c20 * phi(i, j, k + 1, 0) +
67 Stencil::c21 * phi(i, j, k, 0) +
68 Stencil::c22 * phi(i, j, k - 1, 0)) *
69 idx[2];
70 const amrex::Real vz = (Stencil::c20 * phi(i, j, k + 1, 1) +
71 Stencil::c21 * phi(i, j, k, 1) +
72 Stencil::c22 * phi(i, j, k - 1, 1)) *
73 idx[2];
74 const amrex::Real wz = (Stencil::c20 * phi(i, j, k + 1, 2) +
75 Stencil::c21 * phi(i, j, k, 2) +
76 Stencil::c22 * phi(i, j, k - 1, 2)) *
77 idx[2];
78 // N11
79 strphi(i, j, k, 0) = (ux * ux + uy * vx + uz * wx) -
80 (ux * ux + uy * uy + uz * uz) +
81 (3 * (ux * ux + vx * vx + wx * wx)) +
82 (ux * ux + vx * uy + wx * uz);
83 // N12
84 strphi(i, j, k, 1) = (ux * uy + uy * vy + uz * wy) -
85 (ux * vx + uy * vy + uz * vz) +
86 (3 * (ux * uy + vx * vy + wx * wy)) +
87 (ux * vx + vx * vy + wx * vz);
88 // N13
89 strphi(i, j, k, 2) = (ux * uz + uy * vz + uz * wz) -
90 (ux * wx + uy * wy + uz * wz) +
91 (3 * (ux * uz + vx * vz + wx * wz)) +
92 (ux * wx + vx * wy + wx * wz);
93 // N21
94 strphi(i, j, k, 3) = (vx * ux + vy * vx + vz * wx) -
95 (vx * ux + vy * uy + vz * uz) +
96 (3 * (uy * ux + vy * vx + wy * wx)) +
97 (uy * ux + vy * uy + wy * uz);
98 // N22
99 strphi(i, j, k, 4) = (vx * uy + vy * vy + vz * wy) -
100 (vx * vx + vy * vy + vz * vz) +
101 (3 * (uy * vx + vy * vy + wy * wy)) +
102 (uy * vx + vy * vy + wy * vz);
103 // N23
104 strphi(i, j, k, 5) = (vx * uz + vy * vz + vz * wz) -
105 (vx * wx + vy * wy + vz * wz) +
106 (3 * (uy * wx + vy * vz + wy * wz)) +
107 (uy * wx + vy * wy + wy * wz);
108 // N31
109 strphi(i, j, k, 6) = (wx * ux + wy * vx + wz * wx) -
110 (wx * ux + wy * uy + wz * uz) +
111 (3 * (uz * ux + vz * vx + wz * wx)) +
112 (uz * ux + vz * uy + wz * uz);
113 // N32
114 strphi(i, j, k, 7) = (wx * uy + wy * vy + wz * wy) -
115 (wx * vx + wy * vy + wz * vz) +
116 (3 * (uz * uy + vz * vy + wz * wy)) +
117 (uz * vx + vz * vy + wz * vz);
118 // N33
119 strphi(i, j, k, 8) = (wx * uz + wy * vz + wz * wz) -
120 (wx * wx + wy * wy + wz * wz) +
121 (3 * (uz * uz + vz * vz + wz * wz)) +
122 (uz * wx + vz * wy + wz * wz);
123 });
124 }
125
126 FTypeOut& m_strphi;
127 const FTypeIn& m_phi;
128};
129
136template <typename FTypeIn, typename FTypeOut>
137void nonlinearsum(FTypeOut& strphi, const FTypeIn& phi)
138{
139 BL_PROFILE("kynema-sgf::fvm::nonlinearsum");
140 NonLinearSum<FTypeIn, FTypeOut> str(strphi, phi);
141 impl::apply(str, phi);
142}
143
149template <typename FType>
150std::unique_ptr<ScratchField> nonlinearsum(const FType& phi)
151{
152 const std::string gname = phi.name() + "_nonlinearsum";
153 auto strphi = phi.repo().create_scratch_field(gname, 9);
154 nonlinearsum(*strphi, phi);
155 return strphi;
156}
157
158} // namespace kynema_sgf::fvm
159
160#endif /* NONLINEARSUM_H */
void nonlinearsum(FTypeOut &strphi, const FTypeIn &phi)
Definition nonLinearSum.H:137
AMREX_INLINE void apply(const FvmOp &fvmop, const FType &fld)
Definition fvm_utils.H:20
Definition SchemeTraits.H:6
Definition nonLinearSum.H:14
const FTypeIn & m_phi
Definition nonLinearSum.H:127
NonLinearSum(FTypeOut &strphi, const FTypeIn &phi)
Definition nonLinearSum.H:15
FTypeOut & m_strphi
Definition nonLinearSum.H:126
void apply(const int lev) const
Definition nonLinearSum.H:23