CHROMA
multi_syssolver_mdagm_cg_wilson_quda_w.cc
Go to the documentation of this file.
1 /*! \file
2  * \brief Solve a MdagM*psi=chi linear system by CG2
3  */
4 
9 #include "quda.h"
10 
11 #include <cstdlib>
12 
13 namespace Chroma
14 {
15 
16  //! CG2 system solver namespace
17  namespace MdagMMultiSysSolverCGQudaWilsonEnv
18  {
19  //! Callback function
21  const std::string& path,
22  Handle< FermState< LatticeFermion, multi1d<LatticeColorMatrix>, multi1d<LatticeColorMatrix> > > state,
24  {
25  return new MdagMMultiSysSolverCGQudaWilson(A, state,SysSolverQUDAWilsonParams(xml_in, path));
26  }
27 
28  //! Name to be used
29  const std::string name("MULTI_CG_QUDA_WILSON_INVERTER");
30 
31  //! Local registration flag
32  static bool registered = false;
33 
34  //! Register all the factories
35  bool registerAll()
36  {
37  bool success = true;
38  if (! registered)
39  {
41  registered = true;
42  }
43  return success;
44  }
45  }
46 
47  SystemSolverResults_t
48  MdagMMultiSysSolverCGQudaWilson::qudaInvertMulti(const T& chi_s,
49  multi1d<T>& psi_s,
50  const multi1d<Real> shifts) const{
51 
52  SystemSolverResults_t ret;
53 
54 
55  void *spinorIn=(void *)&(chi_s.elem(rb[1].start()).elem(0).elem(0).real());
56 
57 
58  void** spinorOut = new void*[shifts.size()];
59 
60  // This is a hangover from mallocing. I expect now new will throw an exception
61  if (spinorOut == nullptr ) {
62  QDPIO::cerr << "Couldn't allocate spinorOut" << std::endl;
63  QDP_abort(1);
64  }
65 
66  double resid_sq;
67 
68  psi_s.resize( shifts.size());
69 
70  if ( shifts.size() > QUDA_MAX_MULTI_SHIFT ) {
71  QDPIO::cerr << "You want more shifts than QUDA_MAX_MULTI_SHIFT" << std::endl;
72  QDPIO::cerr << "Requested : " << shifts.size() << " QUDA_MAX_MULTI_SHIFT=" << QUDA_MAX_MULTI_SHIFT << std::endl;
73  QDP_abort(1);
74  }
75 
76  for(int s=0; s < shifts.size(); s++) {
77  psi_s[s][ rb[1] ] = zero;
78  spinorOut[s] = (void *)&(psi_s[s].elem(rb[1].start()).elem(0).elem(0).real());
79  quda_inv_param.offset[s] = toDouble(shifts[s]);
80  }
81  quda_inv_param.num_offset = shifts.size();
82 
83  for (int i=0; i< quda_inv_param.num_offset; i++) quda_inv_param.tol_offset[i] = quda_inv_param.tol;
84  // Do the solve here
85  StopWatch swatch1;
86  swatch1.reset();
87  swatch1.start();
88  QDPIO::cout << "CALLING QUDA SOLVER" << std::endl << std::flush ;
89  invertMultiShiftQuda(spinorOut, spinorIn, (QudaInvertParam*)&quda_inv_param);
90  swatch1.stop();
91 
92  // Tidy Up
93  delete [] spinorOut;
94 
95 
96  QDPIO::cout << "QUDA_"<<solver_string<<"_WILSON_SOLVER: time="<< quda_inv_param.secs <<" s" ;
97  QDPIO::cout << "\tPerformance="<< quda_inv_param.gflops/quda_inv_param.secs<<" GFLOPS" ;
98  QDPIO::cout << "\tTotal Time (incl. load gauge)=" << swatch1.getTimeInSeconds() <<" s"<<std::endl;
99 
100  ret.n_count =quda_inv_param.iter;
101 
102  return ret;
103 
104  }
105 
106 
107 }
108 
Support class for fermion actions and linear operators.
Definition: state.h:94
Class for counted reference semantics.
Definition: handle.h:33
static T & Instance()
Definition: singleton.h:432
Register MdagM system solvers.
Solve a MdagM*psi=chi linear system by CG2 using CG.
Factory for producing system solvers for MdagM*psi = chi.
MdagMMultiSystemSolver< LatticeFermion > * createFerm(XMLReader &xml_in, const std::string &path, Handle< FermState< LatticeFermion, multi1d< LatticeColorMatrix >, multi1d< LatticeColorMatrix > > > state, Handle< LinearOperator< LatticeFermion > > A)
Callback function.
const std::string name("MULTI_CG_QUDA_WILSON_INVERTER")
Name to be used.
Asqtad Staggered-Dirac operator.
Definition: klein_gord.cc:10
int i
Definition: pbg5p_w.cc:55
A(A, psi, r, Ncb, PLUS)
const WilsonTypeFermAct< multi1d< LatticeFermion > > Handle< const ConnectState > state
Definition: pbg5p_w.cc:28
Double zero
Definition: invbicg.cc:106
multi1d< LatticeFermion > s(Ncb)
::std::string string
Definition: gtest.h:1979
LatticeFermion T
Definition: t_clover.cc:11