@@ -73,37 +73,29 @@ inline void correct(const Grid_t& grid, E1& efield, const Int3& efield_ib,
7373 Int3 lx, rx, ly, ry, lz, rz;
7474 detail::find_limits (grid, p, lx, rx, ly, ry, lz, rz);
7575
76- if (!grid.isInvar (0 )) {
77- Int3 l = lx, r = rx;
78- auto ex = efield.view (_all, _all, _all, 0 , p);
79- auto res = mf.view (_all, _all, _all, 0 , p);
80- ex.view (_s (l[0 ], r[0 ]), _s (l[1 ], r[1 ]), _s (l[2 ], r[2 ])) =
81- ex.view (_s (l[0 ], r[0 ]), _s (l[1 ], r[1 ]), _s (l[2 ], r[2 ])) +
82- (res.view (_s (l[0 ] + 1 , r[0 ] + 1 ), _s (l[1 ], r[1 ]), _s (l[2 ], r[2 ])) -
83- res.view (_s (l[0 ], r[0 ]), _s (l[1 ], r[1 ]), _s (l[2 ], r[2 ]))) *
84- fac[0 ];
85- }
76+ Int3 ls[3 ] = {lx, ly, lz};
77+ Int3 rs[3 ] = {rx, ry, rz};
8678
87- {
88- Int3 l = ly, r = ry;
89- auto ey = efield.view (_all, _all, _all, 1 , p);
90- auto res = mf.view (_all, _all, _all, 0 , p);
91- ey.view (_s (l[0 ], r[0 ]), _s (l[1 ], r[1 ]), _s (l[2 ], r[2 ])) =
92- ey.view (_s (l[0 ], r[0 ]), _s (l[1 ], r[1 ]), _s (l[2 ], r[2 ])) +
93- (res.view (_s (l[0 ], r[0 ]), _s (l[1 ] + 1 , r[1 ] + 1 ), _s (l[2 ], r[2 ])) -
94- res.view (_s (l[0 ], r[0 ]), _s (l[1 ], r[1 ]), _s (l[2 ], r[2 ]))) *
95- fac[1 ];
96- }
79+ auto res = mf.view (_all, _all, _all, 0 , p);
80+ for (int d = 0 ; d < 3 ; d++) {
81+ if (grid.isInvar (d)) {
82+ continue ;
83+ }
84+
85+ Int3 l = ls[d];
86+ Int3 r = rs[d];
87+ auto e_comp = efield.view (_all, _all, _all, d, p);
88+
89+ gt::gslice s1x = _s (l[0 ], r[0 ]);
90+ gt::gslice s1y = _s (l[1 ], r[1 ]);
91+ gt::gslice s1z = _s (l[2 ], r[2 ]);
92+
93+ gt::gslice s2[3 ] = {_s (l[0 ], r[0 ]), _s (l[1 ], r[1 ]), _s (l[2 ], r[2 ])};
94+ s2[d] = _s (l[d] + 1 , r[d] + 1 );
9795
98- {
99- Int3 l = lz, r = rz;
100- auto ez = efield.view (_all, _all, _all, 2 , p);
101- auto res = mf.view (_all, _all, _all, 0 , p);
102- ez.view (_s (l[0 ], r[0 ]), _s (l[1 ], r[1 ]), _s (l[2 ], r[2 ])) =
103- ez.view (_s (l[0 ], r[0 ]), _s (l[1 ], r[1 ]), _s (l[2 ], r[2 ])) +
104- (res.view (_s (l[0 ], r[0 ]), _s (l[1 ], r[1 ]), _s (l[2 ] + 1 , r[2 ] + 1 )) -
105- res.view (_s (l[0 ], r[0 ]), _s (l[1 ], r[1 ]), _s (l[2 ], r[2 ]))) *
106- fac[2 ];
96+ e_comp.view (s1x, s1y, s1z) =
97+ e_comp.view (s1x, s1y, s1z) +
98+ (res.view (s2[0 ], s2[1 ], s2[2 ]) - res.view (s1x, s1y, s1z)) * fac[d];
10799 }
108100 }
109101}
0 commit comments