From cd2fb68905bef09c00ac09d73ac8082685d7b74e Mon Sep 17 00:00:00 2001 From: Peter Boyle Date: Tue, 30 Jun 2015 15:01:26 +0100 Subject: [PATCH] big commit fixing nocompiles in defective C++11 compilers (gcc, icpc). stared getting to near the bleeding edge I guess --- tests/Test_GaugeAction.cc | 2 +- tests/Test_cf_coarsen_support.cc | 2 +- tests/Test_cf_cr_unprec.cc | 2 +- tests/Test_dwf_cg_prec.cc | 2 +- tests/Test_dwf_cg_schur.cc | 2 +- tests/Test_dwf_cg_unprec.cc | 2 +- tests/Test_dwf_cr_unprec.cc | 2 +- tests/Test_dwf_hdcr.cc | 218 ++++++++++++++++++++++++++----- tests/Test_gamma.cc | 8 +- tests/Test_main.cc | 30 ++--- tests/Test_nersc_io.cc | 4 +- tests/Test_quenched_update.cc | 2 +- tests/Test_wilson_cg_prec.cc | 2 +- tests/Test_wilson_cg_unprec.cc | 6 - tests/Test_wilson_cr_unprec.cc | 4 - 15 files changed, 218 insertions(+), 70 deletions(-) diff --git a/tests/Test_GaugeAction.cc b/tests/Test_GaugeAction.cc index 56634341..1ce8306e 100644 --- a/tests/Test_GaugeAction.cc +++ b/tests/Test_GaugeAction.cc @@ -53,7 +53,7 @@ int main (int argc, char ** argv) readNerscConfiguration(Umu,header,file); for(int mu=0;mu(Umu,mu); + U[mu] = PeekIndex(Umu,mu); } // Painful ; fix syntactical niceness : to check reader diff --git a/tests/Test_cf_coarsen_support.cc b/tests/Test_cf_coarsen_support.cc index 9d81a798..4a4cf3e9 100644 --- a/tests/Test_cf_coarsen_support.cc +++ b/tests/Test_cf_coarsen_support.cc @@ -41,7 +41,7 @@ int main (int argc, char ** argv) std::vector U(4,UGrid); for(int mu=0;mu(Umu,mu); + U[mu] = PeekIndex(Umu,mu); } RealD mass=0.1; diff --git a/tests/Test_cf_cr_unprec.cc b/tests/Test_cf_cr_unprec.cc index 9c1aea60..e74fc85e 100644 --- a/tests/Test_cf_cr_unprec.cc +++ b/tests/Test_cf_cr_unprec.cc @@ -38,7 +38,7 @@ int main (int argc, char ** argv) std::vector U(4,UGrid); for(int mu=0;mu(Umu,mu); + U[mu] = PeekIndex(Umu,mu); } RealD mass=0.1; diff --git a/tests/Test_dwf_cg_prec.cc b/tests/Test_dwf_cg_prec.cc index 32a0c3ae..28848ff9 100644 --- a/tests/Test_dwf_cg_prec.cc +++ b/tests/Test_dwf_cg_prec.cc @@ -38,7 +38,7 @@ int main (int argc, char ** argv) std::vector U(4,UGrid); for(int mu=0;mu(Umu,mu); + U[mu] = PeekIndex(Umu,mu); } RealD mass=0.1; diff --git a/tests/Test_dwf_cg_schur.cc b/tests/Test_dwf_cg_schur.cc index d080045b..8415c489 100644 --- a/tests/Test_dwf_cg_schur.cc +++ b/tests/Test_dwf_cg_schur.cc @@ -38,7 +38,7 @@ int main (int argc, char ** argv) std::vector U(4,UGrid); for(int mu=0;mu(Umu,mu); + U[mu] = PeekIndex(Umu,mu); } RealD mass=0.1; diff --git a/tests/Test_dwf_cg_unprec.cc b/tests/Test_dwf_cg_unprec.cc index 7a13be43..f4f63c50 100644 --- a/tests/Test_dwf_cg_unprec.cc +++ b/tests/Test_dwf_cg_unprec.cc @@ -38,7 +38,7 @@ int main (int argc, char ** argv) std::vector U(4,UGrid); for(int mu=0;mu(Umu,mu); + U[mu] = PeekIndex(Umu,mu); } RealD mass=0.1; diff --git a/tests/Test_dwf_cr_unprec.cc b/tests/Test_dwf_cr_unprec.cc index 6cd2b32a..2cf7065c 100644 --- a/tests/Test_dwf_cr_unprec.cc +++ b/tests/Test_dwf_cr_unprec.cc @@ -43,7 +43,7 @@ int main (int argc, char ** argv) std::vector U(4,UGrid); for(int mu=0;mu(Umu,mu); + U[mu] = PeekIndex(Umu,mu); } ConjugateResidual MCR(1.0e-8,10000); diff --git a/tests/Test_dwf_hdcr.cc b/tests/Test_dwf_hdcr.cc index 44b05714..8d1a5fce 100644 --- a/tests/Test_dwf_hdcr.cc +++ b/tests/Test_dwf_hdcr.cc @@ -1,11 +1,12 @@ #include #include +#include using namespace std; using namespace Grid; using namespace Grid::QCD; -template +template class MultiGridPreconditioner : public LinearFunction< Lattice > { public: @@ -21,19 +22,22 @@ public: Aggregates & _Aggregates; CoarseOperator & _CoarseOperator; + Matrix & _Matrix; FineOperator & _FineOperator; // Constructor - MultiGridPreconditioner(Aggregates &Agg, CoarseOperator &Coarse, FineOperator &Fine) + MultiGridPreconditioner(Aggregates &Agg, CoarseOperator &Coarse, FineOperator &Fine,Matrix &FineMatrix) : _Aggregates(Agg), _CoarseOperator(Coarse), - _FineOperator(Fine) + _FineOperator(Fine), + _Matrix(FineMatrix) { } - +#if 0 void operator()(const FineField &in, FineField & out) { FineField Min(in._grid); + FineField tmp(in._grid); CoarseVector Csrc(_CoarseOperator.Grid()); CoarseVector Ctmp(_CoarseOperator.Grid()); @@ -44,17 +48,40 @@ public: _Aggregates.PromoteFromSubspace(Csrc,out); std::cout<<"Completeness: "< MCR(1.0e-2,1000); - ConjugateGradient CG(1.0e-2,10000); + // Build some solvers + ConjugateGradient fCG(1.0e-1,1000); + ConjugateGradient CG(1.0e-8,100000); + + //////////////////////////////////////////////////////////////////////// + // ADEF2: [PTM+Q] in = [1 - Q A] M in + Q in = Min + Q [ in -A Min] + //////////////////////////////////////////////////////////////////////// - // [PTM+Q] in = [1 - Q A] M in + Q in = Min + Q [ in -A Min] // Smoothing step, followed by coarse grid correction + MdagMLinearOperator MdagMOp(_Matrix); - MCR(_FineOperator,in,Min); - _FineOperator.Op(Min,out); - out = in -out; // out = in - A Min + Min=in; + std::cout<< " Preconditioner in " << norm2(in)< MdagMOp(_CoarseOperator); + fCG(MdagMOp,tmp,out); + + _FineOperator.Op(out,tmp); + + std::cout<< " Preconditioner in " << norm2(in)< MdagMOp(_CoarseOperator); HermitianLinearOperator HermOp(_CoarseOperator); Csol=zero; _Aggregates.ProjectToSubspace (Csrc,out); @@ -63,6 +90,121 @@ public: _Aggregates.PromoteFromSubspace(Csol,out); out = Min + out;; + */ + + } +#endif + + //////////////////////////////////////////////////////////////////////// + // ADEF2: [PTM+Q] in = [1 - Q A] M in + Q in = Min + Q [ in -A Min] + // ADEF1: [MP+Q ] in =M [1 - A Q] in + Q in + //////////////////////////////////////////////////////////////////////// +#if 0 + void operator()(const FineField &in, FineField & out) { + + CoarseVector Csrc(_CoarseOperator.Grid()); + CoarseVector Ctmp(_CoarseOperator.Grid()); + CoarseVector Csol(_CoarseOperator.Grid()); + + ConjugateGradient CG(1.0e-10,100000); + ConjugateGradient fCG(3.0e-2,1000); + + HermitianLinearOperator HermOp(_CoarseOperator); + MdagMLinearOperator MdagMOp(_CoarseOperator); + MdagMLinearOperator fMdagMOp(_Matrix); + + FineField tmp(in._grid); + FineField res(in._grid); + FineField Min(in._grid); + + // Monitor completeness of low mode space + _Aggregates.ProjectToSubspace (Csrc,in); + _Aggregates.PromoteFromSubspace(Csrc,out); + std::cout<<"Coarse Grid Preconditioner\nCompleteness in: "< CG(1.0e-10,100000); + ConjugateGradient fCG(1.0e-3,1000); + + HermitianLinearOperator HermOp(_CoarseOperator); + MdagMLinearOperator MdagMOp(_CoarseOperator); + MdagMLinearOperator fMdagMOp(_Matrix); + + FineField tmp(in._grid); + FineField res(in._grid); + FineField Qin(in._grid); + + // Monitor completeness of low mode space + // _Aggregates.ProjectToSubspace (Csrc,in); + // _Aggregates.PromoteFromSubspace(Csrc,out); + // std::cout<<"Coarse Grid Preconditioner\nCompleteness in: "< clatt = GridDefaultLatt(); for(int d=0;d seeds4({1,2,3,4}); @@ -96,7 +239,9 @@ int main (int argc, char ** argv) GridParallelRNG RNG4(UGrid); RNG4.SeedFixedIntegers(seeds4); GridParallelRNG CRNG(Coarse5d);CRNG.SeedFixedIntegers(cseeds); - LatticeFermion src(FGrid); gaussian(RNG5,src); + Gamma g5(Gamma::Gamma5); + + LatticeFermion src(FGrid); gaussian(RNG5,src);// src=src+g5*src; LatticeFermion result(FGrid); result=zero; LatticeFermion ref(FGrid); ref=zero; LatticeFermion tmp(FGrid); @@ -112,7 +257,7 @@ int main (int argc, char ** argv) // SU3::HotConfiguration(RNG4,Umu); // Umu=zero; - RealD mass=0.04; + RealD mass=0.01; RealD M5=1.8; std::cout << "**************************************************"<< std::endl; @@ -120,7 +265,7 @@ int main (int argc, char ** argv) std::cout << "**************************************************"<< std::endl; DomainWallFermion Ddwf(Umu,*FGrid,*FrbGrid,*UGrid,*UrbGrid,mass,M5); - const int nbasis = 4; + const int nbasis = 6; typedef Aggregation Subspace; typedef CoarsenedMatrix CoarseOperator; @@ -132,7 +277,12 @@ int main (int argc, char ** argv) MdagMLinearOperator HermDefOp(Ddwf); Subspace Aggregates(Coarse5d,FGrid); Aggregates.CreateSubspace(RNG5,HermDefOp); - +// for(int i=0;i PosdefLdop(LDOp); - ConjugateGradient CG(1.0e-6,10000); + ConjugateGradient CG(1.0e-6,100000); CG(PosdefLdop,c_src,c_res); - std::cout << "**************************************************"<< std::endl; - std::cout << "Solving indef-MCR on coarse space "<< std::endl; - std::cout << "**************************************************"<< std::endl; - HermitianLinearOperator HermIndefLdop(LDOp); - ConjugateResidual MCR(1.0e-6,10000); + // std::cout << "**************************************************"<< std::endl; + // std::cout << "Solving indef-MCR on coarse space "<< std::endl; + // std::cout << "**************************************************"<< std::endl; + // HermitianLinearOperator HermIndefLdop(LDOp); + // ConjugateResidual MCR(1.0e-6,100000); //MCR(HermIndefLdop,c_src,c_res); std::cout << "**************************************************"<< std::endl; std::cout << "Building deflation preconditioner "<< std::endl; std::cout << "**************************************************"<< std::endl; - MultiGridPreconditioner Precon(Aggregates, LDOp,HermIndefOp); + MultiGridPreconditioner Precon(Aggregates, LDOp,HermIndefOp,Ddwf); + TrivialPrecon simple; std::cout << "**************************************************"<< std::endl; - std::cout << "Building a one level PGCR "<< std::endl; + std::cout << "Unprec CG "<< std::endl; std::cout << "**************************************************"<< std::endl; - TrivialPrecon simple; - PrecGeneralisedConjugateResidual GCR(1.0e-6,10000,simple,8,64); - GCR(HermIndefOp,src,result); + // TrivialPrecon simple; + ConjugateGradient fCG(1.0e-8,100000); + fCG(HermDefOp,src,result); + + std::cout << "**************************************************"<< std::endl; + std::cout << "Testing GCR on indef matrix "<< std::endl; + std::cout << "**************************************************"<< std::endl; + // PrecGeneralisedConjugateResidual UPGCR(1.0e-8,100000,simple,8,128); + // UPGCR(HermIndefOp,src,result); std::cout << "**************************************************"<< std::endl; std::cout << "Building a two level PGCR "<< std::endl; std::cout << "**************************************************"<< std::endl; - PrecGeneralisedConjugateResidual PGCR(1.0e-6,10000,Precon,8,64); + PrecGeneralisedConjugateResidual PGCR(1.0e-8,100000,Precon,8,128); + std::cout<<"checking norm src "<(cVec,scalar,1); + pokeIndex<2>(cVec,scalar,1); scalar=transpose(scalar); - scalar=transposeIndex<1>(scalar); - scalar=traceIndex<1>(scalar); - scalar=peekIndex<1>(cVec,0); + scalar=TransposeIndex(scalar); + scalar=TraceIndex(scalar); + scalar=PeekIndex(cVec,0); scalar=trace(scalar); scalar=localInnerProduct(cVec,cVec); @@ -230,8 +230,8 @@ int main (int argc, char ** argv) cm = ProjectOnGroup(cm); std::cout << cm << " " << std::endl; - det = Determinant(cm); - std::cout << "determinant: " << det << std::endl; + // det = Determinant(cm); + // std::cout << "determinant: " << det << std::endl; // Foo = Foo+scalar; // LatticeColourMatrix+Scalar @@ -257,16 +257,16 @@ int main (int argc, char ** argv) SpinMatrix s_m; SpinColourMatrix sc_m; - s_m = traceIndex<1>(sc_m); // Map to traceColour - c_m = traceIndex<2>(sc_m); // map to traceSpin + s_m = TensorIndexRecursion::traceIndex(sc_m); // Map to traceColour + c_m = TensorIndexRecursion::traceIndex(sc_m); // map to traceSpin - c = traceIndex<2>(s_m); - c = traceIndex<1>(c_m); + c = TensorIndexRecursion::traceIndex(s_m); + c = TensorIndexRecursion::traceIndex(c_m); - s_m = peekIndex<1>(scm,0,0); - c_m = peekIndex<2>(scm,1,2); + s_m = TensorIndexRecursion::peekIndex(scm,0,0); + c_m = TensorIndexRecursion::peekIndex(scm,1,2); // c_m = peekSpin(scm,1,2); - c_m = peekIdiot(scm,1,2); + // c_m = peekIdiot(scm,1,2); printf("c. Level %d\n",c_m.TensorLevel); printf("c. Level %d\n",c_m().TensorLevel); @@ -277,7 +277,7 @@ int main (int argc, char ** argv) c = scm()(1,1)(1,2); scm()(1,1)(2,1) = c; - pokeIndex<1> (c_m,c,0,0); + pokeIndex (c_m,c,0,0); } FooBar = Bar; @@ -346,7 +346,7 @@ int main (int argc, char ** argv) LatticeGaugeField U(&Fine); // LatticeColourMatrix Uy = peekLorentz(U,1); - LatticeColourMatrix Uy = peekDumKopf(U,1); + // LatticeColourMatrix Uy = peekDumKopf(U,1); flops = ncall*1.0*volume*(8*Nc*Nc*Nc); bytes = ncall*1.0*volume*Nc*Nc *2*3*sizeof(Grid::Real); diff --git a/tests/Test_nersc_io.cc b/tests/Test_nersc_io.cc index b1b9d7ae..f4afb6d0 100644 --- a/tests/Test_nersc_io.cc +++ b/tests/Test_nersc_io.cc @@ -10,7 +10,7 @@ int main (int argc, char ** argv) Grid_init(&argc,&argv); - std::vector simd_layout = GridDefaultSimd(4,vComplexF::Nsimd()); + std::vector simd_layout = GridDefaultSimd(4,vComplex::Nsimd()); std::vector mpi_layout = GridDefaultMpi(); std::vector latt_size ({16,16,16,32}); std::vector clatt_size ({4,4,4,8}); @@ -31,7 +31,7 @@ int main (int argc, char ** argv) readNerscConfiguration(Umu,header,file); for(int mu=0;mu(Umu,mu); + U[mu] = PeekIndex(Umu,mu); } // Painful ; fix syntactical niceness diff --git a/tests/Test_quenched_update.cc b/tests/Test_quenched_update.cc index 2d5ea2ee..887a27a8 100644 --- a/tests/Test_quenched_update.cc +++ b/tests/Test_quenched_update.cc @@ -63,7 +63,7 @@ int main (int argc, char ** argv) ColourWilsonLoops::Staple(staple,Umu,mu); staple = adj(staple); - link = peekIndex(Umu,mu); + link = PeekIndex(Umu,mu); for( int subgroup=0;subgroup U(4,&Grid); for(int mu=0;mu(Umu,mu); + U[mu] = PeekIndex(Umu,mu); } RealD mass=0.5; diff --git a/tests/Test_wilson_cg_unprec.cc b/tests/Test_wilson_cg_unprec.cc index 13996017..4e42af42 100644 --- a/tests/Test_wilson_cg_unprec.cc +++ b/tests/Test_wilson_cg_unprec.cc @@ -34,16 +34,10 @@ int main (int argc, char ** argv) LatticeFermion result(&Grid); result=zero; LatticeGaugeField Umu(&Grid); random(pRNG,Umu); - std::vector U(4,&Grid); - double volume=1; for(int mu=0;mu(Umu,mu); - } RealD mass=0.5; WilsonFermion Dw(Umu,Grid,RBGrid,mass); diff --git a/tests/Test_wilson_cr_unprec.cc b/tests/Test_wilson_cr_unprec.cc index 7fbf8c8f..1aa04632 100644 --- a/tests/Test_wilson_cr_unprec.cc +++ b/tests/Test_wilson_cr_unprec.cc @@ -42,10 +42,6 @@ int main (int argc, char ** argv) volume=volume*latt_size[mu]; } - for(int mu=0;mu(Umu,mu); - } - RealD mass=0.5; WilsonFermion Dw(Umu,Grid,RBGrid,mass);