@@ -169,34 +169,20 @@ MLCGSolverT<MF>::solve_bicgstab (MF& sol, const MF& rhs, RT eps_rel, RT eps_abs)
169169 return ret;
170170 }
171171
172+ RT rho = dotxy (rh,r);
173+ PrintDots (Post_Calc," Dots(rh,r)" ,rh,r);
174+ p.LocalCopy (r,0 ,0 ,ncomp,nghost); // This is the true initialization. Move this to here to avoid iter==1 check every iteration
175+ PrintDots (Post_Calc, " Dots(p,p)" ,p,p);
176+
172177 for (; iter <= maxiter; ++iter)
173178 {
174179 Half_Iter = TitleHalf_Iter + std::to_string (iter);
175180 Iteration = TitleIteration + std::to_string (iter);
176- Post_Calc = TitlePost_Calc + std::to_string (iter-1 );
177-
178- const RT rho = dotxy (rh,r);
179- PrintDots (Post_Calc," Dots(rh,r)" ,rh,r);
180-
181+ Post_Calc = TitlePost_Calc + std::to_string (iter);
181182 if ( rho == 0 )
182183 {
183184 ret = 1 ; break ;
184185 }
185- if ( iter == 1 )
186- {
187- p.LocalCopy (r,0 ,0 ,ncomp,nghost);
188- PrintDots (Post_Calc, " Dots(p,p)" ,p,p);
189- }
190- else
191- {
192- const RT beta = (rho/rho_1)*(alpha/omega);
193- // amrex::Print() << Iteration << " beta: " << beta << '\n';
194- PrintDots (Post_Calc, " Dots(p,p) before " ,p,p);
195- MF::Saxpy (p, -omega, v, 0 , 0 , ncomp, nghost); // p += -omega*v
196- PrintDots (Post_Calc, " Dots(p,p) after Saxpy" ,p,p);
197- MF::Xpay (p, beta, r, 0 , 0 , ncomp, nghost); // p = r + beta*p
198- PrintDots (Post_Calc, " Dots(p,p) after Xpay" ,p,p);
199- }
200186 ph.LocalCopy (p,0 ,0 ,ncomp,nghost);
201187 PrintDots (Half_Iter," Dots(ph,ph) before Lp.apply" ,ph,ph);
202188 Lp.apply (amrlev, mglev, v, ph, MLLinOpT<MF>::BCMode::Homogeneous, MLLinOpT<MF>::StateMode::Correction);
@@ -286,6 +272,15 @@ MLCGSolverT<MF>::solve_bicgstab (MF& sol, const MF& rhs, RT eps_rel, RT eps_abs)
286272 ret = 4 ; break ;
287273 }
288274 rho_1 = rho;
275+ rho = dotxy (rh,r);
276+ PrintDots (Post_Calc," Dots(rh,r)" ,rh,r);
277+ const RT beta = (rho/rho_1)*(alpha/omega);
278+ // amrex::Print() << Iteration << " beta: " << beta << '\n';
279+ PrintDots (Post_Calc, " Dots(p,p) before " ,p,p);
280+ MF::Saxpy (p, -omega, v, 0 , 0 , ncomp, nghost); // p += -omega*v
281+ PrintDots (Post_Calc, " Dots(p,p) after Saxpy" ,p,p);
282+ MF::Xpay (p, beta, r, 0 , 0 , ncomp, nghost); // p = r + beta*p
283+ PrintDots (Post_Calc, " Dots(p,p) after Xpay" ,p,p);
289284 }
290285
291286 if ( verbose > 0 )
0 commit comments