portelli/Grid
1
0
Fork 0
mirror of https://github.com/paboyle/Grid.git synced 2025-06-15 14:27:06 +01:00
Code Releases Activity

Added a mixed precision multishift algorithm for which the matrix multiplies are performed in single precision but the search directions are accumulated in double precision.

Browse Source 
A reliable update step is performed at a tunable frequency to correct the residual. A final mixed-prec single-shift solve is performed on each pole to perform cleanup if necessary.
A test is provided to demonstrate the algorithm.
This commit is contained in:
Christopher Kelly
2021-01-06 12:21:30 -05:00
parent 1fb41a4300
commit 1b84f59273
3 changed files with 504 additions and 0 deletions
Download Patch File Download Diff File Expand all files Collapse all files

93
tests/solver/Test_dwf_multishift_mixedprec.cc Normal file
View File

@ -0,0 +1,93 @@
/*************************************************************************************
Grid physics library, www.github.com/paboyle/Grid
Source file: ./tests/Test_dwf_multishift_mixedprec.cc
Copyright (C) 2015
Author: Christopher Kelly <ckelly@bnl.gov>
This program is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
This program is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License along
with this program; if not, write to the Free Software Foundation, Inc.,
51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
See the full license in the file "LICENSE" in the top level distribution directory
*************************************************************************************/
/* END LEGAL */
#include <Grid/Grid.h>
using namespace Grid;
int main (int argc, char ** argv)
{
Grid_init(&argc, &argv);
const int Ls = 16;
GridCartesian* UGrid_d = SpaceTimeGrid::makeFourDimGrid(GridDefaultLatt(), GridDefaultSimd(Nd, vComplexD::Nsimd()), GridDefaultMpi());
GridRedBlackCartesian* UrbGrid_d = SpaceTimeGrid::makeFourDimRedBlackGrid(UGrid_d);
GridCartesian* FGrid_d = SpaceTimeGrid::makeFiveDimGrid(Ls, UGrid_d);
GridRedBlackCartesian* FrbGrid_d = SpaceTimeGrid::makeFiveDimRedBlackGrid(Ls, UGrid_d);
GridCartesian* UGrid_f = SpaceTimeGrid::makeFourDimGrid(GridDefaultLatt(), GridDefaultSimd(Nd, vComplexF::Nsimd()), GridDefaultMpi());
GridRedBlackCartesian* UrbGrid_f = SpaceTimeGrid::makeFourDimRedBlackGrid(UGrid_f);
GridCartesian* FGrid_f = SpaceTimeGrid::makeFiveDimGrid(Ls, UGrid_f);
GridRedBlackCartesian* FrbGrid_f = SpaceTimeGrid::makeFiveDimRedBlackGrid(Ls, UGrid_f);
std::vector<int> seeds4({1, 2, 3, 4});
std::vector<int> seeds5({5, 6, 7, 8});
GridParallelRNG RNG5(FGrid_d);
RNG5.SeedFixedIntegers(seeds5);
GridParallelRNG RNG4(UGrid_d);
RNG4.SeedFixedIntegers(seeds4);
LatticeFermionD src_d(FGrid_d);
random(RNG5, src_d);
LatticeGaugeFieldD Umu_d(UGrid_d);
SU<Nc>::HotConfiguration(RNG4, Umu_d);
LatticeGaugeFieldF Umu_f(UGrid_f);
precisionChange(Umu_f, Umu_d);
std::cout << GridLogMessage << "Lattice dimensions: " << GridDefaultLatt() << " Ls: " << Ls << std::endl;
RealD mass = 0.01;
RealD M5 = 1.8;
DomainWallFermionD Ddwf_d(Umu_d, *FGrid_d, *FrbGrid_d, *UGrid_d, *UrbGrid_d, mass, M5);
DomainWallFermionF Ddwf_f(Umu_f, *FGrid_f, *FrbGrid_f, *UGrid_f, *UrbGrid_f, mass, M5);
LatticeFermionD src_o_d(FrbGrid_d);
pickCheckerboard(Odd, src_o_d, src_d);
SchurDiagMooeeOperator<DomainWallFermionD, LatticeFermionD> HermOpEO_d(Ddwf_d);
SchurDiagMooeeOperator<DomainWallFermionF, LatticeFermionF> HermOpEO_f(Ddwf_f);
AlgRemez remez(1e-4, 64, 50);
int order = 15;
remez.generateApprox(order, 1, 2); //sqrt
MultiShiftFunction shifts(remez, 1e-10, false);
int relup_freq = 50;
double t1=usecond();
ConjugateGradientMultiShiftMixedPrec<LatticeFermionD,LatticeFermionF> mcg(10000, shifts, FrbGrid_f, HermOpEO_f, relup_freq);
std::vector<LatticeFermionD> results_o_d(order, FrbGrid_d);
mcg(HermOpEO_d, src_o_d, results_o_d);
double t2=usecond();
std::cout<<GridLogMessage << "Test: Total usec = "<< (t2-t1)<<std::endl;
Grid_finalize();
}