CHROMA
syssolver_mdagm_clover_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 "io/aniso_io.h"
10 
11 
12 #include "handle.h"
15 #include "meas/glue/mesplq.h"
16 // QUDA Headers
17 #include <quda.h>
18 // #include <util_quda.h>
19 
20 
21 //#undef BUILD_QUDA_DEVIFACE_GAUGE
22 //#undef BUILD_QUDA_DEVIFACE_SPINOR
23 //#undef BUILD_QUDA_DEVIFACE_CLOVER
24 
25 
26 namespace Chroma
27 {
28  namespace MdagMSysSolverQUDACloverEnv
29  {
30 
31  //! Anonymous namespace
32  namespace
33  {
34  //! Name to be used
35  const std::string name("QUDA_CLOVER_INVERTER");
36 
37  //! Local registration flag
38  bool registered = false;
39  }
40 
41 
42 
44  const std::string& path,
45  Handle< FermState< LatticeFermion, multi1d<LatticeColorMatrix>, multi1d<LatticeColorMatrix> > > state,
46 
48  {
49  return new MdagMSysSolverQUDAClover(A, state,SysSolverQUDACloverParams(xml_in, path));
50  }
51 
52  //! Register all the factories
53  bool registerAll()
54  {
55  bool success = true;
56  if (! registered)
57  {
59  registered = true;
60  }
61  return success;
62  }
63  }
64 
66  MdagMSysSolverQUDAClover::qudaInvert(const CloverTermT<T, U>& clover,
67  const CloverTermT<T, U>& invclov,
68  const T& chi_s,
69  T& psi_s) const{
70 
72 
73  void *spinorIn;
74  void *spinorOut;
75 
76 #ifdef BUILD_QUDA_DEVIFACE_SPINOR
77  std::vector<QDPCache::ArgKey> ids;
78 #endif
79 
80  // Solve A_oo - D A^{-1}_ee D -- chroma conventions.
81  // No need to transform source
82 #ifndef BUILD_QUDA_DEVIFACE_SPINOR
83  spinorIn =(void *)&(chi_s.elem(rb[1].start()).elem(0).elem(0).real());
84 #else
85  //spinorIn = GetMemoryPtr( chi_s.getId() );
86  ids.push_back(chi_s.getId());
87 #endif
88 
89 #ifndef BUILD_QUDA_DEVIFACE_SPINOR
90  spinorOut =(void *)&(psi_s.elem(rb[1].start()).elem(0).elem(0).real());
91 #else
92  ids.push_back(psi_s.getId());
93  auto dev_ptr = GetMemoryPtr(ids);
94  spinorIn = dev_ptr[0];
95  spinorOut = dev_ptr[1];
96  //void* spinorOut = GetMemoryPtr( psi_s.getId() );
97 #endif
98 
99  // Do the solve here
100  StopWatch swatch1;
101  swatch1.reset();
102  swatch1.start();
103  invertQuda(spinorOut, spinorIn, (QudaInvertParam*)&quda_inv_param);
104  swatch1.stop();
105 
106 
107 
108  QDPIO::cout << "QUDA_"<<solver_string<<"_CLOVER_SOLVER: time="<< quda_inv_param.secs <<" s" ;
109  QDPIO::cout << "\tPerformance="<< quda_inv_param.gflops/quda_inv_param.secs<<" GFLOPS" ;
110  QDPIO::cout << "\tTotal Time (incl. load gauge)=" << swatch1.getTimeInSeconds() <<" s"<<std::endl;
111 
112  ret.n_count =quda_inv_param.iter;
113 
114  return ret;
115 
116  }
117 
118 
119 }
120 
121 
122 
123 // DEAD Test Code
124 #if 0
125  void* spinorIn =(void *)&(mod_chi.elem(rb[1].start()).elem(0).elem(0).real());
126  void* spinorOut =(void *)&(psi_s.elem(rb[0].start()).elem(0).elem(0).real());
127 
128  // OK Here I have a chance to test directly
129  // Even Target Checkerboard, No Dagger
130  psi_s = zero;
131  quda_inv_param.dslash_type = QUDA_CLOVER_WILSON_DSLASH; // Sets Clover Matrix
132  dslashQuda(spinorOut, spinorIn, &quda_inv_param, 0, 0);
133 
134 
135 
136  // Need to create a simple ferm state from the links_single...
138  const AnisoParam_t& aniso = invParam.CloverParams.anisoParam;
139  QDPWilsonDslashT<T,Q,Q> qdp_dslash(pstate, aniso);
140 
141  T tmp,psi2;
142  tmp=zero;
143  psi2=zero;
144  // qdp_dslash.apply(psi2, mod_chi, PLUS, 0);
145  qdp_dslash.apply(tmp, mod_chi, PLUS, 0);
146  invclov.apply(psi2,tmp, PLUS, 0);
147 
148 
149  T r=zero;
150  r = psi2 - psi_s;
151 
152  QDPIO::cout << "CB=0" << std::endl;
153  QDPIO::cout << "Dslash Test: || r || = " << sqrt(norm2(r,rb[0])) << std::endl;
154  // QDPIO::cout << "Dslash Test: || r ||/|| psi || = " << sqrt(norm2(r,rb[0])/norm2(psi_s, rb[0])) << std::endl;
155 
156  QDPIO::cout << "CB=1: Should be zero" << std::endl;
157  QDPIO::cout << "Dslash Test: || r || = " << sqrt(norm2(r,rb[1])) << std::endl;
158  //QDPIO::cout << "Dslash Test: || r ||/|| psi || = " << sqrt(norm2(r,rb[1])/norm2(psi_s, rb[1])) << std::endl;
159 
160  const int* tab = rb[0].siteTable().slice();
161  for(int i=0; i < rb[0].numSiteTable(); i++) {
162  int j = tab[i];
163  bool printSite=false;
164 
165  for(int spin=0; spin < 4; spin++) {
166  for(int col=0; col < 3; col++) {
167  if( (fabs(r.elem(j).elem(spin).elem(col).real()) > 1.0e-5 )
168  || (fabs(r.elem(j).elem(spin).elem(col).imag()) > 1.0e-5 )) {
169  printSite=true;
170  }
171  }
172  }
173  if( printSite ) {
174 
175  for(int spin=0; spin < 4; spin++) {
176  for(int col=0; col < 3; col++) {
177  QDPIO::cout << "Site= " << j << " Spin= "<< spin << " Col= " << col << " spinor = ( "
178  << psi2.elem(j).elem(spin).elem(col).real() << " , "
179  << psi2.elem(j).elem(spin).elem(col).imag() << " )" << std::endl;
180  }
181  }
182  QDPIO::cout << std::endl;
183  }
184  }
185  QDP_abort(1);
186 #endif
Anisotropy parameters.
Support class for fermion actions and linear operators.
Definition: state.h:94
Class for counted reference semantics.
Definition: handle.h:33
Periodic version of FermState.
void apply(T &chi, const T &psi, enum PlusMinus isign, int cb) const
General Wilson-Dirac dslash.
Definition: lwldslash_w.h:48
static T & Instance()
Definition: singleton.h:432
Class for counted reference semantics.
unsigned j
Definition: ldumul_w.cc:35
unsigned i
Definition: ldumul_w.cc:34
Wilson Dslash linear operator.
Double tmp
Definition: meslate.cc:60
static bool registered
Local registration flag.
const std::string name
Name to be used.
MdagMSystemSolver< LatticeFermion > * createFerm(XMLReader &xml_in, const std::string &path, Handle< FermState< LatticeFermion, multi1d< LatticeColorMatrix >, multi1d< LatticeColorMatrix > > > state, Handle< LinearOperator< LatticeFermion > > A)
bool registerAll()
Register all the factories.
Asqtad Staggered-Dirac operator.
Definition: klein_gord.cc:10
LinOpSysSolverMGProtoClover::T T
@ PLUS
Definition: chromabase.h:45
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
::std::string string
Definition: gtest.h:1979
multi1d< LatticeFermion > r(Ncb)
Periodic ferm state and a creator.
Parameters for anisotropy.
Definition: aniso_io.h:24
Holds return info from SystemSolver call.
Definition: syssolver.h:17
Register MdagM system solvers.
Solve a MdagM*psi=chi linear system by BiCGStab.
Factory for producing system solvers for MdagM*psi = chi.