/************************************************************************************* Grid physics library, www.github.com/paboyle/Grid Source file: MultiRHSDeflation.h Copyright (C) 2023 Author: Peter Boyle 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 */ #pragma once NAMESPACE_BEGIN(Grid); /* Need helper object for BLAS accelerated mrhs projection i) MultiRHS Deflation Import Evecs -> nev x vol x internal Import vector of Lattice objects -> nrhs x vol x internal => Cij (nrhs x Nev) via GEMM. => Guess (nrhs x vol x internal) = C x evecs (via GEMM) Export ii) MultiRHS block projection Import basis -> nblock x nbasis x (block x internal) Import vector of fine lattice objects -> nblock x nrhs x (block x internal) => coarse_(nrhs x nbasis )^block = via batched GEMM iii) Alternate interface: Import higher dim Lattice object-> vol x nrhs layout */ template class MultiRHSDeflation { public: typedef typename Field::scalar_type scalar; typedef typename Field::scalar_object scalar_object; int nev; std::vector eval; GridBase *grid; uint64_t vol; uint64_t words; deviceVector BLAS_E; // nev x vol -- the eigenbasis (up to a 1/sqrt(lambda)) deviceVector BLAS_R; // nrhs x vol -- the sources deviceVector BLAS_G; // nrhs x vol -- the guess deviceVector BLAS_C; // nrhs x nev -- the coefficients MultiRHSDeflation(){}; ~MultiRHSDeflation(){ Deallocate(); }; void Deallocate(void) { nev=0; grid=nullptr; vol=0; words=0; BLAS_E.resize(0); BLAS_R.resize(0); BLAS_C.resize(0); BLAS_G.resize(0); } void Allocate(int _nev,GridBase *_grid) { nev=_nev; grid=_grid; vol = grid->lSites(); words = sizeof(scalar_object)/sizeof(scalar); eval.resize(nev); BLAS_E.resize (vol * words * nev ); std::cout << GridLogMessage << " Allocate for "< &evec,std::vector &_eval) { ImportEigenBasis(evec,_eval,0,evec.size()); } // Could use to import a batch of eigenvectors void ImportEigenBasis(std::vector &evec,std::vector &_eval, int _ev0, int _nev) { assert(_ev0+_nev<=evec.size()); Allocate(_nev,evec[0].Grid()); // Imports a sub-batch of eigenvectors, _ev0, ..., _ev0+_nev-1 for(int e=0;e &source,std::vector & guess) { int nrhs = source.size(); assert(source.size()==guess.size()); assert(grid == guess[0].Grid()); conformable(guess[0],source[0]); int64_t vw = vol * words; RealD t0 = usecond(); BLAS_R.resize(nrhs * vw); // cost free if size doesn't change BLAS_G.resize(nrhs * vw); // cost free if size doesn't change BLAS_C.resize(nev * nrhs);// cost free if size doesn't change ///////////////////////////////////////////// // Copy in the multi-rhs sources ///////////////////////////////////////////// // for(int r=0;r