Skip to content

Commit 166e971

Browse files
ProExpertProgProExpertProg
committed
Update exp_* signatures, save hyper, use dt member
1 parent 0fe404f commit 166e971

File tree

2 files changed

+44
-50
lines changed

2 files changed

+44
-50
lines changed

‎cpp/lib/Naive.cpp‎

Lines changed: 27 additions & 40 deletions
Original file line numberDiff line numberDiff line change
@@ -62,7 +62,6 @@ void Naive::run(Dim N, Dim saveInterval) {
6262

6363
bool divergent = false, repeat = false, noInc = false;
6464
int divergentCount = 0, repeatCount = 0;
65-
HyperCoefficients hyper{};
6665

6766
bool saved = false;
6867
// Manually increment t only if not diverging
@@ -116,16 +115,14 @@ void Naive::run(Dim N, Dim saveInterval) {
116115
g.for_each_kxky([&](Dim kx, Dim ky) {
117116
GM_Nonlinear_K(kx, ky, N_E) =
118117
nonlinear::N(bracketPhiNE_K(kx, ky), bracketAParUEKPar_K(kx, ky));
119-
GM_K_Star(kx, ky, N_E) =
120-
exp_nu(kx, ky, hyper.nu_2, dt) * moments_K(kx, ky, N_E) +
121-
dt / 2.0 * (1 + exp_nu(kx, ky, hyper.nu_2, dt)) * GM_Nonlinear_K(kx, ky, N_E);
118+
GM_K_Star(kx, ky, N_E) = exp_nu(kx, ky) * moments_K(kx, ky, N_E) +
119+
dt / 2.0 * (1 + exp_nu(kx, ky)) * GM_Nonlinear_K(kx, ky, N_E);
122120

123121
GM_Nonlinear_K(kx, ky, A_PAR) =
124122
nonlinear::A(bracketAParPhiG2Ne_K(kx, ky), bracketUEParPhi_K(kx, ky), g.kPerp2(kx, ky));
125-
GM_K_Star(kx, ky, A_PAR) =
126-
exp_eta(kx, ky, hyper.eta2, dt) * moments_K(kx, ky, A_PAR) +
127-
dt / 2.0 * (1 + exp_eta(kx, ky, hyper.eta2, dt)) * GM_Nonlinear_K(kx, ky, A_PAR) +
128-
(1.0 - exp_eta(kx, ky, hyper.eta2, dt)) * aParEq_K(kx, ky);
123+
GM_K_Star(kx, ky, A_PAR) = exp_eta(kx, ky) * moments_K(kx, ky, A_PAR) +
124+
dt / 2.0 * (1 + exp_eta(kx, ky)) * GM_Nonlinear_K(kx, ky, A_PAR) +
125+
(1.0 - exp_eta(kx, ky)) * aParEq_K(kx, ky);
129126
});
130127

131128
if (g.M > 2) {
@@ -150,17 +147,14 @@ void Naive::run(Dim N, Dim saveInterval) {
150147
g.for_each_kxky([&](Dim kx, Dim ky) {
151148
GM_Nonlinear_K(kx, ky, G_MIN) = nonlinear::G2(
152149
bracketPhiG2_K(kx, ky), bracketAParG3_K(kx, ky), bracketAParUEKPar_K(kx, ky));
153-
GM_K_Star(kx, ky, G_MIN) =
154-
exp_nu(kx, ky, hyper.nu_2, dt) * moments_K(kx, ky, G_MIN) +
155-
dt / 2.0 * (1 + exp_nu(kx, ky, hyper.nu_2, dt)) * GM_Nonlinear_K(kx, ky, G_MIN);
150+
GM_K_Star(kx, ky, G_MIN) = exp_nu(kx, ky) * moments_K(kx, ky, G_MIN) +
151+
dt / 2.0 * (1 + exp_nu(kx, ky)) * GM_Nonlinear_K(kx, ky, G_MIN);
156152

157153
GM_Nonlinear_K(kx, ky, LAST) =
158154
nonlinear::GLast(bracketPhiGLast_K(kx, ky), bracketTotalGLast_K(kx, ky));
159155
GM_K_Star(kx, ky, LAST) =
160-
exp_gm(LAST, hyper.nu_ei, dt) * exp_nu(kx, ky, hyper.nu_g, dt) *
161-
moments_K(kx, ky, LAST) +
162-
dt / 2.0 * (1 + exp_gm(LAST, hyper.nu_ei, dt) * exp_nu(kx, ky, hyper.nu_g, dt)) *
163-
GM_Nonlinear_K(kx, ky, LAST);
156+
exp_gm(LAST) * exp_nu_g(kx, ky) * moments_K(kx, ky, LAST) +
157+
dt / 2.0 * (1 + exp_gm(LAST) * exp_nu_g(kx, ky)) * GM_Nonlinear_K(kx, ky, LAST);
164158
});
165159

166160
auto dGMinusPlus = g.dBufXY();
@@ -179,9 +173,8 @@ void Naive::run(Dim N, Dim saveInterval) {
179173
GM_Nonlinear_K(kx, ky, m) =
180174
nonlinear::GM(m, bracketPhiGM_K(kx, ky), bracketAParGMMinusPlus_K(kx, ky));
181175
GM_K_Star(kx, ky, m) =
182-
exp_gm(m, hyper.nu_ei, dt) * exp_nu(kx, ky, hyper.nu_g, dt) * moments_K(kx, ky, m) +
183-
dt / 2.0 * (1 + exp_gm(m, hyper.nu_ei, dt) * exp_nu(kx, ky, hyper.nu_g, dt)) *
184-
GM_Nonlinear_K(kx, ky, m);
176+
exp_gm(m) * exp_nu_g(kx, ky) * moments_K(kx, ky, m) +
177+
dt / 2.0 * (1 + exp_gm(m) * exp_nu_g(kx, ky)) * GM_Nonlinear_K(kx, ky, m);
185178
});
186179
}
187180
}
@@ -251,13 +244,12 @@ void Naive::run(Dim N, Dim saveInterval) {
251244
GM_Nonlinear_K_Loop(kx, ky, A_PAR) = nonlinear::A(
252245
bracketAParPhiG2Ne_K_Loop(kx, ky), bracketUEParPhi_K_Loop(kx, ky), g.kPerp2(kx, ky));
253246
// TODO(OPT) reuse star
254-
momentsNew_K(kx, ky, A_PAR) =
255-
1.0 / (1.0 + semiImplicitOperator(kx, ky) / 4.0) *
256-
(exp_eta(kx, ky, hyper.eta2, dt) * moments_K(kx, ky, A_PAR) +
257-
dt / 2.0 * exp_eta(kx, ky, hyper.eta2, dt) * GM_Nonlinear_K(kx, ky, A_PAR) +
258-
dt / 2.0 * GM_Nonlinear_K_Loop(kx, ky, A_PAR) +
259-
(1.0 - exp_eta(kx, ky, hyper.eta2, dt)) * aParEq_K(kx, ky) +
260-
semiImplicitOperator(kx, ky) / 4.0 * guessAPar_K(kx, ky));
247+
momentsNew_K(kx, ky, A_PAR) = 1.0 / (1.0 + semiImplicitOperator(kx, ky) / 4.0) *
248+
(exp_eta(kx, ky) * moments_K(kx, ky, A_PAR) +
249+
dt / 2.0 * exp_eta(kx, ky) * GM_Nonlinear_K(kx, ky, A_PAR) +
250+
dt / 2.0 * GM_Nonlinear_K_Loop(kx, ky, A_PAR) +
251+
(1.0 - exp_eta(kx, ky)) * aParEq_K(kx, ky) +
252+
semiImplicitOperator(kx, ky) / 4.0 * guessAPar_K(kx, ky));
261253
ueKPar_K_New(kx, ky) = -g.kPerp2(kx, ky) * momentsNew_K(kx, ky, A_PAR);
262254

263255
sumAParRelError += std::norm(momentsNew_K(kx, ky, A_PAR) - moments_K(kx, ky, A_PAR));
@@ -285,10 +277,9 @@ void Naive::run(Dim N, Dim saveInterval) {
285277
GM_Nonlinear_K_Loop(kx, ky, N_E) =
286278
nonlinear::N(bracketPhiNE_K_Loop(kx, ky), bracketAParUEKPar_K_Loop(kx, ky));
287279
// TODO(OPT) reuse star
288-
momentsNew_K(kx, ky, N_E) =
289-
exp_nu(kx, ky, hyper.nu_2, dt) * moments_K(kx, ky, N_E) +
290-
dt / 2.0 * exp_nu(kx, ky, hyper.nu_2, dt) * GM_Nonlinear_K(kx, ky, N_E) +
291-
dt / 2.0 * GM_Nonlinear_K_Loop(kx, ky, N_E);
280+
momentsNew_K(kx, ky, N_E) = exp_nu(kx, ky) * moments_K(kx, ky, N_E) +
281+
dt / 2.0 * exp_nu(kx, ky) * GM_Nonlinear_K(kx, ky, N_E) +
282+
dt / 2.0 * GM_Nonlinear_K_Loop(kx, ky, N_E);
292283

293284
phi_K_New(kx, ky) =
294285
(kx | ky) == 0 ? 0 : nonlinear::phi(momentsNew_K(kx, ky, N_E), g.kPerp2(kx, ky));
@@ -307,10 +298,9 @@ void Naive::run(Dim N, Dim saveInterval) {
307298
nonlinear::G2(bracketPhiG2_K_Loop(kx, ky), bracketAParG3_K_Loop(kx, ky),
308299
bracketAParUEKPar_K_Loop(kx, ky));
309300
// TODO(OPT) reuse star
310-
momentsNew_K(kx, ky, G_MIN) =
311-
exp_nu(kx, ky, hyper.nu_2, dt) * moments_K(kx, ky, G_MIN) +
312-
dt / 2.0 * exp_nu(kx, ky, hyper.nu_2, dt) * GM_Nonlinear_K(kx, ky, G_MIN) +
313-
dt / 2.0 * GM_Nonlinear_K_Loop(kx, ky, G_MIN);
301+
momentsNew_K(kx, ky, G_MIN) = exp_nu(kx, ky) * moments_K(kx, ky, G_MIN) +
302+
dt / 2.0 * exp_nu(kx, ky) * GM_Nonlinear_K(kx, ky, G_MIN) +
303+
dt / 2.0 * GM_Nonlinear_K_Loop(kx, ky, G_MIN);
314304
});
315305
DerivateNewMoment(G_MIN);
316306

@@ -332,9 +322,8 @@ void Naive::run(Dim N, Dim saveInterval) {
332322
bracketAParGMMinusPlus_K_Loop(kx, ky));
333323
// TODO(OPT) reuse star
334324
momentsNew_K(kx, ky, m) =
335-
exp_gm(m, hyper.nu_ei, dt) * exp_nu(kx, ky, hyper.nu_g, dt) * moments_K(kx, ky, m) +
336-
dt / 2.0 * exp_gm(m, hyper.nu_ei, dt) * exp_nu(kx, ky, hyper.nu_g, dt) *
337-
GM_Nonlinear_K(kx, ky, m) +
325+
exp_gm(m) * exp_nu_g(kx, ky) * moments_K(kx, ky, m) +
326+
dt / 2.0 * exp_gm(m) * exp_nu_g(kx, ky) * GM_Nonlinear_K(kx, ky, m) +
338327
dt / 2.0 * GM_Nonlinear_K_Loop(kx, ky, m);
339328
});
340329

@@ -363,10 +352,8 @@ void Naive::run(Dim N, Dim saveInterval) {
363352
nonlinear::GLast(bracketPhiGLast_K_Loop(kx, ky), bracketTotalGLast_K_Loop(kx, ky));
364353
// TODO(OPT) reuse star
365354
momentsNew_K(kx, ky, LAST) =
366-
exp_gm(LAST, hyper.nu_ei, dt) * exp_nu(kx, ky, hyper.nu_g, dt) *
367-
moments_K(kx, ky, LAST) +
368-
dt / 2.0 * exp_gm(LAST, hyper.nu_ei, dt) * exp_nu(kx, ky, hyper.nu_g, dt) *
369-
GM_Nonlinear_K(kx, ky, LAST) +
355+
exp_gm(LAST) * exp_nu_g(kx, ky) * moments_K(kx, ky, LAST) +
356+
dt / 2.0 * exp_gm(LAST) * exp_nu_g(kx, ky) * GM_Nonlinear_K(kx, ky, LAST) +
370357
dt / 2.0 * GM_Nonlinear_K_Loop(kx, ky, LAST);
371358
});
372359
DerivateNewMoment(LAST);

‎cpp/lib/Naive.hpp‎

Lines changed: 17 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -42,14 +42,14 @@ class Naive : public ahr::HermiteRunner {
4242
Brackets br{g, tf, hlFilter};
4343

4444
private:
45-
4645
using View = Grid::View;
4746
using Buf = Grid::Buf;
4847

4948
const Real XYNorm{1.0 / Real(g.X) / Real(g.Y)}; ///< Normalization factor for FFT
5049

51-
Real dt{-1}; ///< timestep
52-
Real elapsedT{0.0}; ///< total time elapsed
50+
Real dt{-1}; ///< timestep
51+
Real elapsedT{0.0}; ///< total time elapsed
52+
HyperCoefficients hyper{}; ///< Hyper-coefficients, updated every dt
5353

5454
void fftHL(View::R_XY in, View::C_XY out); ///< FFT with Hou-Li Filter
5555

@@ -97,17 +97,24 @@ class Naive : public ahr::HermiteRunner {
9797
// TODO other file/class
9898
// =================
9999

100-
[[nodiscard]] Real exp_nu(Dim kx, Dim ky, Real nu2, Real dt) const {
101-
return std::exp(-(nu * g.kPerp2(kx, ky) + nu2 * std::pow(g.kPerp2(kx, ky), hyper_order)) * dt);
100+
[[nodiscard]] Real exp_nu(Dim kx, Dim ky) const {
101+
return std::exp(
102+
-(nu * g.kPerp2(kx, ky) + hyper.nu_2 * std::pow(g.kPerp2(kx, ky), hyper_order)) * dt);
103+
}
104+
105+
[[nodiscard]] Real exp_nu_g(Dim kx, Dim ky) const {
106+
return std::exp(
107+
-(nu * g.kPerp2(kx, ky) + hyper.nu_g * std::pow(g.kPerp2(kx, ky), hyper_order)) * dt);
102108
}
103109

104-
[[nodiscard]] Real exp_gm(Dim m, Real hyper_nuei, Real dt) const {
105-
return std::exp(-(Real(m) * nu_ei + std::pow(m, 2 * hyper_morder) * hyper_nuei) * dt);
110+
[[nodiscard]] Real exp_gm(Dim m) const {
111+
return std::exp(-(Real(m) * nu_ei + std::pow(m, 2 * hyper_morder) * hyper.nu_ei) * dt);
106112
}
107113

108-
[[nodiscard]] Real exp_eta(Dim kx, Dim ky, Real res2, Real dt) const {
109-
return std::exp(-(res * g.kPerp2(kx, ky) + res2 * std::pow(g.kPerp2(kx, ky), hyper_order)) *
110-
dt / (1.0 + g.kPerp2(kx, ky) * de * de));
114+
[[nodiscard]] Real exp_eta(Dim kx, Dim ky) const {
115+
return std::exp(
116+
-(res * g.kPerp2(kx, ky) + hyper.eta2 * std::pow(g.kPerp2(kx, ky), hyper_order)) * dt /
117+
(1.0 + g.kPerp2(kx, ky) * de * de));
111118
}
112119

113120
/// getTimestep calculates flows and magnetic fields to determine a dt.

0 commit comments

Comments
 (0)