1
1
#include " Naive.hpp"
2
2
#include " equillibrium.hpp"
3
+ #include " reducer.hpp"
3
4
4
5
#include < cnpy.h>
5
6
#include < cstdlib>
@@ -255,7 +256,7 @@ void Naive::run(Dim N, Dim saveInterval) {
255
256
256
257
// / f_pred from Viriato
257
258
auto GM_Nonlinear_K_Loop = g.cBufMXY ();
258
- Real sumAParRelError = 0 ;
259
+ SumReducer< Real> sumAParRelError = 0 ;
259
260
g.for_each_kxky ([&](Dim kx, Dim ky) {
260
261
GM_Nonlinear_K_Loop (kx, ky, A_PAR) = nonlinear::A (
261
262
bracketAParPhiG2Ne_K_Loop (kx, ky), bracketUEParPhi_K_Loop (kx, ky), g.kPerp2 (kx, ky));
@@ -493,8 +494,7 @@ Naive::Buf::R_XY Naive::getMoment(Dim m) const {
493
494
Real Naive::getTimestep (DxDy<View::R_XY> dPhi, DxDy<View::R_XY> dNE, DxDy<View::R_XY> dAPar) {
494
495
// compute flows
495
496
DxDy<View::R_XY> ve, b;
496
- Real vyMax{0 }, vxMax{0 }, bxMax{0 }, byMax{0 };
497
- bPerpMax = 0 ;
497
+ MaxReducer<Real> vyMax{0 }, vxMax{0 }, bxMax{0 }, byMax{0 }, bPerpMaxRed{0 };
498
498
499
499
// Note that this is minus in Viriato, but we don't care because we're taking the absolute value
500
500
// anyway.
@@ -504,17 +504,20 @@ Real Naive::getTimestep(DxDy<View::R_XY> dPhi, DxDy<View::R_XY> dNE, DxDy<View::
504
504
b.DY = dAPar.DX ;
505
505
506
506
g.for_each_xy ([&](Dim x, Dim y) {
507
- bxMax = std::max (bxMax, std::abs (b.DX (x, y)));
508
- byMax = std::max (byMax, std::abs (b.DY (x, y)));
509
- bPerpMax = std::max (bPerpMax, std::sqrt (b.DX (x, y) * b.DX (x, y) + b.DY (x, y) * b.DY (x, y)));
510
- vxMax = std::max (vxMax, std::abs (ve.DX (x, y)));
511
- vyMax = std::max (vyMax, std::abs (ve.DY (x, y)));
507
+ bxMax = std::max<Real>(bxMax, std::abs (b.DX (x, y)));
508
+ byMax = std::max<Real>(byMax, std::abs (b.DY (x, y)));
509
+ bPerpMaxRed =
510
+ std::max<Real>(bPerpMaxRed, std::sqrt (b.DX (x, y) * b.DX (x, y) + b.DY (x, y) * b.DY (x, y)));
511
+ vxMax = std::max<Real>(vxMax, std::abs (ve.DX (x, y)));
512
+ vyMax = std::max<Real>(vyMax, std::abs (ve.DY (x, y)));
512
513
if (rhoI >= smallRhoI) {
513
- vxMax = std::max (vxMax, rhoS * rhoS * std::abs (dNE.DX (x, y)));
514
- vyMax = std::max (vyMax, rhoS * rhoS * std::abs (dNE.DY (x, y)));
514
+ vxMax = std::max<Real> (vxMax, rhoS * rhoS * std::abs (dNE.DX (x, y)));
515
+ vyMax = std::max<Real> (vyMax, rhoS * rhoS * std::abs (dNE.DY (x, y)));
515
516
}
516
517
});
517
518
519
+ bPerpMax = bPerpMaxRed;
520
+
518
521
Real kperpDum2 = std::pow (g.ky_ (g.KY / 2 ), 2 ) + std::pow (Real (g.KX ), 2 );
519
522
Real omegaKaw;
520
523
if (rhoI < smallRhoI) {
0 commit comments