diff --git a/Grid/algorithms/iterative/AdefGeneric.h b/Grid/algorithms/iterative/AdefGeneric.h index 9ae1e6113..82eb68668 100644 --- a/Grid/algorithms/iterative/AdefGeneric.h +++ b/Grid/algorithms/iterative/AdefGeneric.h @@ -53,7 +53,22 @@ class TwoLevelCG : public LinearFunction // Fine operator, Smoother, CoarseSolver LinearOperatorBase &_FineLinop; LinearFunction &_Smoother; - + + GridStopWatch ProjectTimer; + GridStopWatch PromoteTimer; + GridStopWatch CoarseTimer; + GridStopWatch SmoothTimer; + GridStopWatch MatrixTimer; + GridStopWatch M3Timer; + GridStopWatch LinalgTimer; + + int64_t M3Calls; + int64_t SmoothCalls; + int64_t MatrixCalls; + int64_t ProjectCalls; + int64_t CoarseCalls; + int64_t PromoteCalls; + // more most opertor functions TwoLevelCG(RealD tol, Integer maxit, @@ -103,12 +118,20 @@ class TwoLevelCG : public LinearFunction RealD tn; GridStopWatch HDCGTimer; + ProjectTimer.Reset(); + PromoteTimer.Reset(); + CoarseTimer.Reset(); + SmoothTimer.Reset(); + MatrixTimer.Reset(); + M3Timer.Reset(); + LinalgTimer.Reset(); + M3Calls = SmoothCalls = MatrixCalls = ProjectCalls = CoarseCalls = PromoteCalls = 0; HDCGTimer.Start(); ////////////////////////// // x0 = Vstart -- possibly modify guess ////////////////////////// Vstart(x,src); - + // r0 = b -A x0 _FineLinop.HermOp(x,mmp[0]); axpy (r, -1.0,mmp[0], src); // Recomputes r=src-Ax0 @@ -145,33 +168,40 @@ class TwoLevelCG : public LinearFunction int peri_kp = (k+1) % mmax; rtz=rtzp; + M3Timer.Start(); d= PcgM3(p[peri_k],mmp[peri_k]); + M3Timer.Stop(); + M3Calls++; a = rtz/d; - + // Memorise this pAp[peri_k] = d; - + + LinalgTimer.Start(); axpy(x,a,p[peri_k],x); RealD rn = axpy_norm(r,-a,mmp[peri_k],r); + LinalgTimer.Stop(); // Compute z = M x PcgM1(r,z); - + { RealD n1,n2; n1=norm2(r); n2=norm2(z); std::cout << GridLogMessage<<"HDCG::fPcg iteration "< z + b p + // Standard search direction p -> z + b p b = (rtzp)/rtz; int northog; @@ -202,8 +232,25 @@ class TwoLevelCG : public LinearFunction if ( rn <= rsq ) { HDCGTimer.Stop(); - std::cout< double { + return (n > 0) ? sw.useconds() * 1e-3 / n : 0.0; + }; + std::cout< CoarseField PleftProj(this->coarsegrid); CoarseField PleftMss_proj(this->coarsegrid); - GridStopWatch SmootherTimer; - GridStopWatch MatrixTimer; - SmootherTimer.Start(); + this->SmoothTimer.Start(); this->_Smoother(in,Min); - SmootherTimer.Stop(); + this->SmoothTimer.Stop(); + this->SmoothCalls++; - MatrixTimer.Start(); + this->MatrixTimer.Start(); this->_FineLinop.HermOp(Min,out); - MatrixTimer.Stop(); + this->MatrixTimer.Stop(); + this->MatrixCalls++; axpy(tmp,-1.0,out,in); // tmp = in - A Min - GridStopWatch ProjTimer; - GridStopWatch CoarseTimer; - GridStopWatch PromTimer; - ProjTimer.Start(); - this->_Aggregates.ProjectToSubspace(PleftProj,tmp); - ProjTimer.Stop(); - CoarseTimer.Start(); + this->ProjectTimer.Start(); + this->_Aggregates.ProjectToSubspace(PleftProj,tmp); + this->ProjectTimer.Stop(); + this->ProjectCalls++; + this->CoarseTimer.Start(); this->_CoarseSolver(PleftProj,PleftMss_proj); // Ass^{-1} [in - A Min]_s - CoarseTimer.Stop(); - PromTimer.Start(); - this->_Aggregates.PromoteFromSubspace(PleftMss_proj,tmp);// tmp = Q[in - A Min] - PromTimer.Stop(); - std::cout << GridLogPerformance << "PcgM1 breakdown "<CoarseTimer.Stop(); + this->CoarseCalls++; + this->PromoteTimer.Start(); + this->_Aggregates.PromoteFromSubspace(PleftMss_proj,tmp);// tmp = Q[in - A Min] + this->PromoteTimer.Stop(); + this->PromoteCalls++; axpy(out,1.0,Min,tmp); // Min+tmp }