From f69008edf1766d6cb29401caccfb8cd6274ae4d0 Mon Sep 17 00:00:00 2001 From: Daniel Richtmann Date: Tue, 3 Apr 2018 17:26:49 +0200 Subject: [PATCH] WilsonMG: Add functionality to report timings to MG preconditioner --- tests/solver/Test_wilsonclover_mg.cc | 133 ++++++++++++++++++++++++--- 1 file changed, 121 insertions(+), 12 deletions(-) diff --git a/tests/solver/Test_wilsonclover_mg.cc b/tests/solver/Test_wilsonclover_mg.cc index 5d03e4d6..064b4d6e 100644 --- a/tests/solver/Test_wilsonclover_mg.cc +++ b/tests/solver/Test_wilsonclover_mg.cc @@ -157,6 +157,7 @@ public: virtual void setup() = 0; virtual void operator()(Field const &in, Field &out) = 0; virtual void runChecks() = 0; + virtual void reportTimings() = 0; }; template @@ -180,16 +181,27 @@ public: // Member Data ///////////////////////////////////////////// - int _CurrentLevel; - int _NextCoarserLevel; - MultiGridParams & _MultiGridParams; - LevelInfo & _LevelInfo; - FineMatrix & _FineMatrix; - FineMatrix & _SmootherMatrix; - Aggregates _Aggregates; - CoarseMatrix _CoarseMatrix; + int _CurrentLevel; + int _NextCoarserLevel; + + MultiGridParams &_MultiGridParams; + LevelInfo & _LevelInfo; + + FineMatrix & _FineMatrix; + FineMatrix & _SmootherMatrix; + Aggregates _Aggregates; + CoarseMatrix _CoarseMatrix; + std::unique_ptr _NextPreconditionerLevel; + GridStopWatch _SetupTotalTimer; + GridStopWatch _SetupNextLevelTimer; + GridStopWatch _SolveTotalTimer; + GridStopWatch _SolveRestrictionTimer; + GridStopWatch _SolveProlongationTimer; + GridStopWatch _SolveSmootherTimer; + GridStopWatch _SolveNextLevelTimer; + ///////////////////////////////////////////// // Member Functions ///////////////////////////////////////////// @@ -203,12 +215,17 @@ public: , _SmootherMatrix(SmootherMat) , _Aggregates(_LevelInfo.Grids[_NextCoarserLevel], _LevelInfo.Grids[_CurrentLevel], 0) , _CoarseMatrix(*_LevelInfo.Grids[_NextCoarserLevel]) { + _NextPreconditionerLevel = std::unique_ptr(new NextPreconditionerLevel(_MultiGridParams, _LevelInfo, _CoarseMatrix, _CoarseMatrix)); + + resetTimers(); } void setup() { + _SetupTotalTimer.Start(); + Gamma g5(Gamma::Algebra::Gamma5); MdagMLinearOperator fineMdagMOp(_FineMatrix); @@ -228,7 +245,11 @@ public: _CoarseMatrix.CoarsenOperator(_LevelInfo.Grids[_CurrentLevel], fineMdagMOp, _Aggregates); + _SetupNextLevelTimer.Start(); _NextPreconditionerLevel->setup(); + _SetupNextLevelTimer.Stop(); + + _SetupTotalTimer.Stop(); } virtual void operator()(Lattice const &in, Lattice &out) { @@ -245,6 +266,8 @@ public: void vCycle(Lattice const &in, Lattice &out) { + _SolveTotalTimer.Start(); + RealD inputNorm = norm2(in); CoarseVector coarseSrc(_LevelInfo.Grids[_NextCoarserLevel]); @@ -265,16 +288,26 @@ public: MdagMLinearOperator fineMdagMOp(_FineMatrix); MdagMLinearOperator fineSmootherMdagMOp(_SmootherMatrix); + _SolveRestrictionTimer.Start(); _Aggregates.ProjectToSubspace(coarseSrc, in); + _SolveRestrictionTimer.Stop(); + + _SolveNextLevelTimer.Start(); (*_NextPreconditionerLevel)(coarseSrc, coarseSol); + _SolveNextLevelTimer.Stop(); + + _SolveProlongationTimer.Start(); _Aggregates.PromoteFromSubspace(coarseSol, out); + _SolveProlongationTimer.Stop(); fineMdagMOp.Op(out, fineTmp); fineTmp = in - fineTmp; auto r = norm2(fineTmp); auto residualAfterCoarseGridCorrection = std::sqrt(r / inputNorm); + _SolveSmootherTimer.Start(); fineFGMRES(fineSmootherMdagMOp, in, out); + _SolveSmootherTimer.Stop(); fineMdagMOp.Op(out, fineTmp); fineTmp = in - fineTmp; @@ -284,10 +317,14 @@ public: std::cout << GridLogMG << " Level " << _CurrentLevel << ": V-cycle: Input norm = " << std::sqrt(inputNorm) << " Coarse residual = " << residualAfterCoarseGridCorrection << " Post-Smoother residual = " << residualAfterPostSmoother << std::endl; + + _SolveTotalTimer.Stop(); } void kCycle(Lattice const &in, Lattice &out) { + _SolveTotalTimer.Start(); + RealD inputNorm = norm2(in); CoarseVector coarseSrc(_LevelInfo.Grids[_NextCoarserLevel]); @@ -315,16 +352,26 @@ public: MdagMLinearOperator fineSmootherMdagMOp(_SmootherMatrix); MdagMLinearOperator coarseMdagMOp(_CoarseMatrix); + _SolveRestrictionTimer.Start(); _Aggregates.ProjectToSubspace(coarseSrc, in); + _SolveRestrictionTimer.Stop(); + + _SolveNextLevelTimer.Start(); coarseFGMRES(coarseMdagMOp, coarseSrc, coarseSol); + _SolveNextLevelTimer.Stop(); + + _SolveProlongationTimer.Start(); _Aggregates.PromoteFromSubspace(coarseSol, out); + _SolveProlongationTimer.Stop(); fineMdagMOp.Op(out, fineTmp); fineTmp = in - fineTmp; auto r = norm2(fineTmp); auto residualAfterCoarseGridCorrection = std::sqrt(r / inputNorm); + _SolveSmootherTimer.Start(); fineFGMRES(fineSmootherMdagMOp, in, out); + _SolveSmootherTimer.Stop(); fineMdagMOp.Op(out, fineTmp); fineTmp = in - fineTmp; @@ -334,6 +381,8 @@ public: std::cout << GridLogMG << " Level " << _CurrentLevel << ": K-cycle: Input norm = " << std::sqrt(inputNorm) << " Coarse residual = " << residualAfterCoarseGridCorrection << " Post-Smoother residual = " << residualAfterPostSmoother << std::endl; + + _SolveTotalTimer.Stop(); } void runChecks() { @@ -477,6 +526,34 @@ public: _NextPreconditionerLevel->runChecks(); } + + void reportTimings() { + + // clang-format off + std::cout << GridLogMG << " Level " << _CurrentLevel << ": Time elapsed: Setup total " << _SetupTotalTimer.Elapsed() << std::endl; + std::cout << GridLogMG << " Level " << _CurrentLevel << ": Time elapsed: Setup next level " << _SetupNextLevelTimer.Elapsed() << std::endl; + std::cout << GridLogMG << " Level " << _CurrentLevel << ": Time elapsed: Solve total " << _SolveTotalTimer.Elapsed() << std::endl; + std::cout << GridLogMG << " Level " << _CurrentLevel << ": Time elapsed: Solve restriction " << _SolveRestrictionTimer.Elapsed() << std::endl; + std::cout << GridLogMG << " Level " << _CurrentLevel << ": Time elapsed: Solve prolongation " << _SolveProlongationTimer.Elapsed() << std::endl; + std::cout << GridLogMG << " Level " << _CurrentLevel << ": Time elapsed: Solve smoother " << _SolveSmootherTimer.Elapsed() << std::endl; + std::cout << GridLogMG << " Level " << _CurrentLevel << ": Time elapsed: Solve next level " << _SolveNextLevelTimer.Elapsed() << std::endl; + // clang-format on + + _NextPreconditionerLevel->reportTimings(); + } + + void resetTimers() { + + _SetupTotalTimer.Reset(); + _SetupNextLevelTimer.Reset(); + _SolveTotalTimer.Reset(); + _SolveRestrictionTimer.Reset(); + _SolveProlongationTimer.Reset(); + _SolveSmootherTimer.Reset(); + _SolveNextLevelTimer.Reset(); + + _NextPreconditionerLevel->resetTimers(); + } }; // Specialization for the coarsest level @@ -494,11 +571,16 @@ public: // Member Data ///////////////////////////////////////////// - int _CurrentLevel; + int _CurrentLevel; + MultiGridParams &_MultiGridParams; LevelInfo & _LevelInfo; - FineMatrix & _FineMatrix; - FineMatrix & _SmootherMatrix; + + FineMatrix &_FineMatrix; + FineMatrix &_SmootherMatrix; + + GridStopWatch _SolveTotalTimer; + GridStopWatch _SolveSmootherTimer; ///////////////////////////////////////////// // Member Functions @@ -509,12 +591,17 @@ public: , _MultiGridParams(mgParams) , _LevelInfo(LvlInfo) , _FineMatrix(FineMat) - , _SmootherMatrix(SmootherMat) {} + , _SmootherMatrix(SmootherMat) { + + resetTimers(); + } void setup() {} virtual void operator()(Lattice const &in, Lattice &out) { + _SolveTotalTimer.Start(); + conformable(_LevelInfo.Grids[_CurrentLevel], in._grid); conformable(in, out); @@ -527,10 +614,28 @@ public: MdagMLinearOperator fineMdagMOp(_FineMatrix); + _SolveSmootherTimer.Start(); fineFGMRES(fineMdagMOp, in, out); + _SolveSmootherTimer.Stop(); + + _SolveTotalTimer.Stop(); } void runChecks() {} + + void reportTimings() { + + // clang-format off + std::cout << GridLogMG << " Level " << _CurrentLevel << ": Time elapsed: Solve total " << _SolveTotalTimer.Elapsed() << std::endl; + std::cout << GridLogMG << " Level " << _CurrentLevel << ": Time elapsed: Solve smoother " << _SolveSmootherTimer.Elapsed() << std::endl; + // clang-format on + } + + void resetTimers() { + + _SolveTotalTimer.Reset(); + _SolveSmootherTimer.Reset(); + } }; template @@ -680,6 +785,8 @@ int main(int argc, char **argv) { std::cout << std::endl; } + MGPreconDw->reportTimings(); + std::cout << GridLogMessage << "**************************************************" << std::endl; std::cout << GridLogMessage << "Testing Multigrid for Wilson Clover" << std::endl; std::cout << GridLogMessage << "**************************************************" << std::endl; @@ -701,5 +808,7 @@ int main(int argc, char **argv) { std::cout << std::endl; } + MGPreconDwc->reportTimings(); + Grid_finalize(); }