QUDA  1.0.0
staggered_invertmsrc_test.cpp
Go to the documentation of this file.
1 #include <stdlib.h>
2 #include <stdio.h>
3 #include <time.h>
4 #include <math.h>
5 
6 #include <test_util.h>
7 #include <dslash_util.h>
8 #include <blas_reference.h>
10 #include <quda.h>
11 #include <string.h>
12 #include "misc.h"
13 #include <gauge_field.h>
14 #include <blas_quda.h>
15 
16 #if defined(QMP_COMMS)
17 #include <qmp.h>
18 #elif defined(MPI_COMMS)
19 #include <mpi.h>
20 #endif
21 
22 #define MAX(a,b) ((a)>(b)?(a):(b))
23 #define mySpinorSiteSize 6
24 
25 extern void usage(char** argv);
26 void *qdp_fatlink[4];
27 void *qdp_longlink[4];
28 
29 void *fatlink;
30 void *longlink;
31 
32 #ifdef MULTI_GPU
34 #endif
35 
36 extern int device;
37 extern int niter;
38 extern int Nsrc; // number of spinors to apply to simultaneously
39 
40 
41 
43 extern QudaPrecision prec;
45 
52 
55 
56 extern double tol; // tolerance for inverter
57 extern double tol_hq; // heavy-quark tolerance for inverter
58 extern int test_type;
59 extern int xdim;
60 extern int ydim;
61 extern int zdim;
62 extern int tdim;
63 extern int gridsize_from_cmdline[];
64 
65 // Dirac operator type
67 
69 extern double mass; // the mass of the Dirac operator
70 
71 extern double mass;
72 
73 static void end();
74 
75 template<typename Float>
76 void constructSpinorField(Float *res) {
77  for(int i = 0; i < Vh; i++) {
78  for (int s = 0; s < 1; s++) {
79  for (int m = 0; m < 3; m++) {
80  res[i*(1*3*2) + s*(3*2) + m*(2) + 0] = rand() / (Float)RAND_MAX;
81  res[i*(1*3*2) + s*(3*2) + m*(2) + 1] = rand() / (Float)RAND_MAX;
82  }
83  }
84  }
85 }
86 
87 
88 static void
90  int X1, int X2, int X3, int X4,
93  double mass, double tol, int maxiter, double reliable_delta,
94  double tadpole_coeff
95  )
96 {
97  gaugeParam->X[0] = X1;
98  gaugeParam->X[1] = X2;
99  gaugeParam->X[2] = X3;
100  gaugeParam->X[3] = X4;
101 
102  gaugeParam->cpu_prec = cpu_prec;
103  gaugeParam->cuda_prec = prec;
104  gaugeParam->reconstruct = link_recon;
105  gaugeParam->cuda_prec_sloppy = prec_sloppy;
107  gaugeParam->gauge_fix = QUDA_GAUGE_FIXED_NO;
108  gaugeParam->anisotropy = 1.0;
109  gaugeParam->tadpole_coeff = tadpole_coeff;
110 
113 
114  gaugeParam->scale = dslash_type == QUDA_STAGGERED_DSLASH ? 1.0 : -1.0/(24.0*tadpole_coeff*tadpole_coeff);
115 
116  gaugeParam->t_boundary = QUDA_ANTI_PERIODIC_T;
117  gaugeParam->gauge_order = QUDA_MILC_GAUGE_ORDER;
118  gaugeParam->ga_pad = X1*X2*X3/2;
119 
120  inv_param->verbosity = QUDA_VERBOSE;
121  inv_param->mass = mass;
122 
123  // outer solver parameters
124  inv_param->inv_type = inv_type;
125  inv_param->tol = tol;
126  inv_param->tol_restart = 1e-3; //now theoretical background for this parameter...
127  inv_param->maxiter = niter;
128  inv_param->reliable_delta = 0;//1e-1;
129  inv_param->use_sloppy_partial_accumulator = false;
130  inv_param->pipeline = false;
131 
132  inv_param->Ls = 1;
133 
134 
135  if(tol_hq == 0 && tol == 0){
136  errorQuda("qudaInvert: requesting zero residual\n");
137  exit(1);
138  }
139  // require both L2 relative and heavy quark residual to determine convergence
140  inv_param->residual_type = static_cast<QudaResidualType_s>(0);
141  inv_param->residual_type = (tol != 0) ? static_cast<QudaResidualType_s> ( inv_param->residual_type | QUDA_L2_RELATIVE_RESIDUAL) : inv_param->residual_type;
142  inv_param->residual_type = (tol_hq != 0) ? static_cast<QudaResidualType_s> (inv_param->residual_type | QUDA_HEAVY_QUARK_RESIDUAL) : inv_param->residual_type;
143 
144  inv_param->tol_hq = tol_hq; // specify a tolerance for the residual for heavy quark residual
145 
146  inv_param->Nsteps = 2;
147 
148 
149  //inv_param->inv_type = QUDA_GCR_INVERTER;
150  //inv_param->gcrNkrylov = 10;
151 
152  // domain decomposition preconditioner parameters
154  inv_param->tol_precondition = 1e-1;
155  inv_param->maxiter_precondition = 10;
156  inv_param->verbosity_precondition = QUDA_SILENT;
158 
160  inv_param->solve_type = QUDA_NORMOP_PC_SOLVE;
161  inv_param->matpc_type = QUDA_MATPC_EVEN_EVEN;
162  inv_param->dagger = QUDA_DAG_NO;
164 
165  inv_param->cpu_prec = cpu_prec;
166  inv_param->cuda_prec = prec;
167  inv_param->cuda_prec_sloppy = prec_sloppy;
169  inv_param->gamma_basis = QUDA_DEGRAND_ROSSI_GAMMA_BASIS; // this is meaningless, but must be thus set
170  inv_param->dirac_order = QUDA_DIRAC_ORDER;
171 
172  inv_param->dslash_type = dslash_type;
173 
174  inv_param->sp_pad = X1*X2*X3/2;
176 
179 }
180 
181 
182  int
184 {
187 
188  set_params(&gaugeParam, &inv_param,
189  xdim, ydim, zdim, tdim,
191  link_recon, link_recon_sloppy, mass, tol, 500, 1e-3,
192  0.8);
193 
194  // this must be before the FaceBuffer is created (this is because it allocates pinned memory - FIXME)
195  initQuda(device);
196 
197  setDims(gaugeParam.X);
198  dw_setDims(gaugeParam.X,1); // so we can use 5-d indexing from dwf
200 
201  size_t gSize = (gaugeParam.cpu_prec == QUDA_DOUBLE_PRECISION) ? sizeof(double) : sizeof(float);
202  for (int dir = 0; dir < 4; dir++) {
203  qdp_fatlink[dir] = malloc(V*gaugeSiteSize*gSize);
204  qdp_longlink[dir] = malloc(V*gaugeSiteSize*gSize);
205  }
206  fatlink = malloc(4*V*gaugeSiteSize*gSize);
207  longlink = malloc(4*V*gaugeSiteSize*gSize);
208 
210  &gaugeParam, dslash_type);
211 
212  for(int dir=0; dir<4; ++dir){
213  for(int i=0; i<V; ++i){
214  for(int j=0; j<gaugeSiteSize; ++j){
215  if(gaugeParam.cpu_prec == QUDA_DOUBLE_PRECISION){
216  ((double*)fatlink)[(i*4 + dir)*gaugeSiteSize + j] = ((double*)qdp_fatlink[dir])[i*gaugeSiteSize + j];
217  ((double*)longlink)[(i*4 + dir)*gaugeSiteSize + j] = ((double*)qdp_longlink[dir])[i*gaugeSiteSize + j];
218  }else{
219  ((float*)fatlink)[(i*4 + dir)*gaugeSiteSize + j] = ((float*)qdp_fatlink[dir])[i*gaugeSiteSize + j];
220  ((float*)longlink)[(i*4 + dir)*gaugeSiteSize + j] = ((float*)qdp_longlink[dir])[i*gaugeSiteSize + j];
221  }
222  }
223  }
224  }
225 
226 
228  csParam.nColor=3;
229  csParam.nSpin=1;
230  csParam.nDim=5;
231  for (int d = 0; d < 4; d++) csParam.x[d] = gaugeParam.X[d];
232  csParam.x[0] /= 2;
233  csParam.x[4] = 1;
234 
235  csParam.setPrecision(inv_param.cpu_prec);
236  csParam.pad = 0;
240  csParam.gammaBasis = inv_param.gamma_basis;
241  csParam.create = QUDA_ZERO_FIELD_CREATE;
242  in = new cpuColorSpinorField(csParam);
243  out = new cpuColorSpinorField(csParam);
244  ref = new cpuColorSpinorField(csParam);
245  tmp = new cpuColorSpinorField(csParam);
246 
247  if (inv_param.cpu_prec == QUDA_SINGLE_PRECISION){
248  constructSpinorField((float*)in->V());
249  }else{
250  constructSpinorField((double*)in->V());
251  }
252 
253 #ifdef MULTI_GPU
254  int tmp_value = MAX(ydim*zdim*tdim/2, xdim*zdim*tdim/2);
255  tmp_value = MAX(tmp_value, xdim*ydim*tdim/2);
256  tmp_value = MAX(tmp_value, xdim*ydim*zdim/2);
257 
258  int fat_pad = tmp_value;
259  int link_pad = 3*tmp_value;
260 
261  // FIXME: currently assume staggered is SU(3)
262  gaugeParam.type = dslash_type == QUDA_STAGGERED_DSLASH ?
264  gaugeParam.reconstruct = QUDA_RECONSTRUCT_NO;
265  GaugeFieldParam cpuFatParam(fatlink, gaugeParam);
267  cpuFat = new cpuGaugeField(cpuFatParam);
268  ghost_fatlink = (void**)cpuFat->Ghost();
269 
270  gaugeParam.type = QUDA_ASQTAD_LONG_LINKS;
271  GaugeFieldParam cpuLongParam(longlink, gaugeParam);
272  cpuLongParam.ghostExchange = QUDA_GHOST_EXCHANGE_PAD;
273  cpuLong = new cpuGaugeField(cpuLongParam);
274  ghost_longlink = (void**)cpuLong->Ghost();
275 
276 
277 #else
278  int fat_pad = 0;
279  int link_pad = 0;
280 #endif
281 
282  gaugeParam.type = dslash_type == QUDA_STAGGERED_DSLASH ?
284  gaugeParam.ga_pad = fat_pad;
286  gaugeParam.reconstruct = link_recon;
288  } else {
289  gaugeParam.reconstruct= gaugeParam.reconstruct_sloppy = QUDA_RECONSTRUCT_NO;
290  }
292  loadGaugeQuda(fatlink, &gaugeParam);
293 
295  gaugeParam.type = QUDA_ASQTAD_LONG_LINKS;
296  gaugeParam.ga_pad = link_pad;
297  gaugeParam.reconstruct= link_recon;
299  loadGaugeQuda(longlink, &gaugeParam);
300  }
301 
302  double time0 = -((double)clock()); // Start the timer
303 
304  double nrm2=0;
305  double src2=0;
306  int ret = 0;
307 
308 
309 
310  switch(test_type){
311  case 0: //even
313  inv_param.inv_type = QUDA_GCR_INVERTER;
314  inv_param.gcrNkrylov = 50;
315  }else if(inv_type == QUDA_PCG_INVERTER){
316  inv_param.inv_type = QUDA_PCG_INVERTER;
317  }
318  inv_param.matpc_type = QUDA_MATPC_EVEN_EVEN;
319  #define NUM_SRC 20
320  inv_param.num_src=Nsrc; // number of spinors to apply to simultaneously
321  void* outArray[NUM_SRC];
322  void* inArray[NUM_SRC];
323  int len;
324 
325  cpuColorSpinorField* spinorOutArray[NUM_SRC];
326  cpuColorSpinorField* spinorInArray[NUM_SRC];
327  spinorOutArray[0] = out;
328  spinorInArray[0] = in;
329  // in = new cpuColorSpinorField(csParam);
330  // out = new cpuColorSpinorField(csParam);
331  // ref = new cpuColorSpinorField(csParam);
332  // tmp = new cpuColorSpinorField(csParam);
333 
334  if (inv_param.cpu_prec == QUDA_SINGLE_PRECISION){
335  constructSpinorField((float*)in->V());
336  }else{
337  constructSpinorField((double*)in->V());
338  }
339 
340  for(int i=1;i < inv_param.num_src; i++){
341  spinorOutArray[i] = new cpuColorSpinorField(csParam);
342  spinorInArray[i] = new cpuColorSpinorField(csParam);
343  if (inv_param.cpu_prec == QUDA_SINGLE_PRECISION){
344  constructSpinorField((float*)spinorInArray[i]->V());
345  }else{
346  constructSpinorField((double*)spinorInArray[i]->V());
347  }
348  }
349 
350  for(int i=0;i < inv_param.num_src; i++){
351  outArray[i] = spinorOutArray[i]->V();
352  inArray[i] = spinorInArray[i]->V();
353  // inv_param.offset[i] = 4*masses[i]*masses[i];
354  }
355  invertMultiSrcQuda(outArray, inArray, &inv_param);
356 
357  time0 += clock();
358  time0 /= CLOCKS_PER_SEC;
359 
360 
361 
362 #ifdef MULTI_GPU
364  out, mass, 0, inv_param.cpu_prec, gaugeParam.cpu_prec, tmp, QUDA_EVEN_PARITY);
365 #else
366  matdagmat(ref->V(), qdp_fatlink, qdp_longlink, out->V(), mass, 0, inv_param.cpu_prec, gaugeParam.cpu_prec, tmp->V(), QUDA_EVEN_PARITY);
367 #endif
368 
369  mxpy(in->V(), ref->V(), Vh*mySpinorSiteSize, inv_param.cpu_prec);
370  nrm2 = norm_2(ref->V(), Vh*mySpinorSiteSize, inv_param.cpu_prec);
371  src2 = norm_2(in->V(), Vh*mySpinorSiteSize, inv_param.cpu_prec);
372 
373  for(int i=1; i < inv_param.num_src;i++) delete spinorOutArray[i];
374  for(int i=1; i < inv_param.num_src;i++) delete spinorInArray[i];
375 
376 
377  break;
378 
379  case 1: //odd
381  inv_param.inv_type = QUDA_GCR_INVERTER;
382  inv_param.gcrNkrylov = 50;
383  }else if(inv_type == QUDA_PCG_INVERTER){
384  inv_param.inv_type = QUDA_PCG_INVERTER;
385  }
386 
387  inv_param.matpc_type = QUDA_MATPC_ODD_ODD;
388  invertQuda(out->V(), in->V(), &inv_param);
389  time0 += clock(); // stop the timer
390  time0 /= CLOCKS_PER_SEC;
391 
392 #ifdef MULTI_GPU
394  out, mass, 0, inv_param.cpu_prec, gaugeParam.cpu_prec, tmp, QUDA_ODD_PARITY);
395 #else
396  matdagmat(ref->V(), qdp_fatlink, qdp_longlink, out->V(), mass, 0, inv_param.cpu_prec, gaugeParam.cpu_prec, tmp->V(), QUDA_ODD_PARITY);
397 #endif
398  mxpy(in->V(), ref->V(), Vh*mySpinorSiteSize, inv_param.cpu_prec);
399  nrm2 = norm_2(ref->V(), Vh*mySpinorSiteSize, inv_param.cpu_prec);
400  src2 = norm_2(in->V(), Vh*mySpinorSiteSize, inv_param.cpu_prec);
401 
402  break;
403 
404  case 2: //full spinor
405 
406  errorQuda("full spinor not supported\n");
407  break;
408 
409  case 3: //multi mass CG, even
410  case 4:
411 
412 #define NUM_OFFSETS 12
413 
414  {
415  double masses[NUM_OFFSETS] ={0.06, 0.061, 0.064, 0.070, 0.077, 0.081, 0.1, 0.11, 0.12, 0.13, 0.14, 0.205};
416  inv_param.num_offset = NUM_OFFSETS;
417  // these can be set independently
418  for (int i=0; i<inv_param.num_offset; i++) {
419  inv_param.tol_offset[i] = inv_param.tol;
420  inv_param.tol_hq_offset[i] = inv_param.tol_hq;
421  }
422  void* outArray[NUM_OFFSETS];
423  int len;
424 
425  cpuColorSpinorField* spinorOutArray[NUM_OFFSETS];
426  spinorOutArray[0] = out;
427  for(int i=1;i < inv_param.num_offset; i++){
428  spinorOutArray[i] = new cpuColorSpinorField(csParam);
429  }
430 
431  for(int i=0;i < inv_param.num_offset; i++){
432  outArray[i] = spinorOutArray[i]->V();
433  inv_param.offset[i] = 4*masses[i]*masses[i];
434  }
435 
436  len=Vh;
437 
438  if (test_type == 3) {
439  inv_param.matpc_type = QUDA_MATPC_EVEN_EVEN;
440  } else {
441  inv_param.matpc_type = QUDA_MATPC_ODD_ODD;
442  }
443 
444  invertMultiShiftQuda(outArray, in->V(), &inv_param);
445 
446  cudaDeviceSynchronize();
447  time0 += clock(); // stop the timer
448  time0 /= CLOCKS_PER_SEC;
449 
450  printfQuda("done: total time = %g secs, compute time = %g, %i iter / %g secs = %g gflops\n",
451  time0, inv_param.secs, inv_param.iter, inv_param.secs,
452  inv_param.gflops/inv_param.secs);
453 
454 
455  printfQuda("checking the solution\n");
457  if (inv_param.solve_type == QUDA_NORMOP_SOLVE){
458  //parity = QUDA_EVENODD_PARITY;
459  errorQuda("full parity not supported\n");
460  }else if (inv_param.matpc_type == QUDA_MATPC_EVEN_EVEN){
461  parity = QUDA_EVEN_PARITY;
462  }else if (inv_param.matpc_type == QUDA_MATPC_ODD_ODD){
463  parity = QUDA_ODD_PARITY;
464  }else{
465  errorQuda("ERROR: invalid spinor parity \n");
466  exit(1);
467  }
468  for(int i=0;i < inv_param.num_offset;i++){
469  printfQuda("%dth solution: mass=%f, ", i, masses[i]);
470 #ifdef MULTI_GPU
472  spinorOutArray[i], masses[i], 0, inv_param.cpu_prec,
473  gaugeParam.cpu_prec, tmp, parity);
474 #else
475  matdagmat(ref->V(), qdp_fatlink, qdp_longlink, outArray[i], masses[i], 0, inv_param.cpu_prec, gaugeParam.cpu_prec, tmp->V(), parity);
476 #endif
477 
478  mxpy(in->V(), ref->V(), len*mySpinorSiteSize, inv_param.cpu_prec);
479  double nrm2 = norm_2(ref->V(), len*mySpinorSiteSize, inv_param.cpu_prec);
480  double src2 = norm_2(in->V(), len*mySpinorSiteSize, inv_param.cpu_prec);
481  double hqr = sqrt(blas::HeavyQuarkResidualNorm(*spinorOutArray[i], *ref).z);
482  double l2r = sqrt(nrm2/src2);
483 
484  printfQuda("Shift %d residuals: (L2 relative) tol %g, QUDA = %g, host = %g; (heavy-quark) tol %g, QUDA = %g, host = %g\n",
485  i, inv_param.tol_offset[i], inv_param.true_res_offset[i], l2r,
486  inv_param.tol_hq_offset[i], inv_param.true_res_hq_offset[i], hqr);
487 
488  //emperical, if the cpu residue is more than 1 order the target accuracy, the it fails to converge
489  if (sqrt(nrm2/src2) > 10*inv_param.tol_offset[i]){
490  ret |=1;
491  }
492  }
493 
494  for(int i=1; i < inv_param.num_offset;i++) delete spinorOutArray[i];
495  }
496  break;
497 
498  default:
499  errorQuda("Unsupported test type");
500 
501  }//switch
502 
503  if (test_type <=2){
504 
505  double hqr = sqrt(blas::HeavyQuarkResidualNorm(*out, *ref).z);
506  double l2r = sqrt(nrm2/src2);
507 
508  printfQuda("Residuals: (L2 relative) tol %g, QUDA = %g, host = %g; (heavy-quark) tol %g, QUDA = %g, host = %g\n",
509  inv_param.tol, inv_param.true_res, l2r, inv_param.tol_hq, inv_param.true_res_hq, hqr);
510 
511  printfQuda("done: total time = %g secs, compute time = %g secs, %i iter / %g secs = %g gflops, \n",
512  time0, inv_param.secs, inv_param.iter, inv_param.secs,
513  inv_param.gflops/inv_param.secs);
514  }
515 
516  end();
517  return ret;
518 }
519 
520 
521 
522  static void
523 end(void)
524 {
525  for(int i=0;i < 4;i++){
526  free(qdp_fatlink[i]);
527  free(qdp_longlink[i]);
528  }
529 
530  free(fatlink);
531  free(longlink);
532 
533  delete in;
534  delete out;
535  delete ref;
536  delete tmp;
537 
538  if (cpuFat) delete cpuFat;
539  if (cpuLong) delete cpuLong;
540 
541  endQuda();
542 }
543 
544 
545  void
547 {
548  printfQuda("running the following test:\n");
549 
550  printfQuda("prec sloppy_prec link_recon sloppy_link_recon test_type S_dimension T_dimension\n");
551  printfQuda("%s %s %s %s %s %d/%d/%d %d \n",
555 
556  printfQuda("Grid partition info: X Y Z T\n");
557  printfQuda(" %d %d %d %d\n",
558  dimPartitioned(0),
559  dimPartitioned(1),
560  dimPartitioned(2),
561  dimPartitioned(3));
562 
563  return ;
564 
565 }
566 
567  void
568 usage_extra(char** argv )
569 {
570  printfQuda("Extra options:\n");
571  printfQuda(" --test <0/1> # Test method\n");
572  printfQuda(" 0: Even even spinor CG inverter\n");
573  printfQuda(" 1: Odd odd spinor CG inverter\n");
574  printfQuda(" 3: Even even spinor multishift CG inverter\n");
575  printfQuda(" 4: Odd odd spinor multishift CG inverter\n");
576  printfQuda(" --cpu_prec <double/single/half> # Set CPU precision\n");
577 
578  return ;
579 }
580 int main(int argc, char** argv)
581 {
582  for (int i = 1; i < argc; i++) {
583 
584  if(process_command_line_option(argc, argv, &i) == 0){
585  continue;
586  }
587 
588 
589 
590  if( strcmp(argv[i], "--cpu_prec") == 0){
591  if (i+1 >= argc){
592  usage(argv);
593  }
594  cpu_prec= get_prec(argv[i+1]);
595  i++;
596  continue;
597  }
598 
599  printf("ERROR: Invalid option:%s\n", argv[i]);
600  usage(argv);
601  }
602 
604  prec_sloppy = prec;
605  }
608  }
609 
610  if(inv_type != QUDA_CG_INVERTER){
611  if(test_type != 0 && test_type != 1) errorQuda("Preconditioning is currently not supported in multi-shift solver solvers");
612  }
613 
614 
615  // initialize QMP/MPI, QUDA comms grid and RNG (test_util.cpp)
616  initComms(argc, argv, gridsize_from_cmdline);
617 
619 
620  printfQuda("dslash_type = %d\n", dslash_type);
621 
622  int ret = invert_test();
623 
624  // finalize the communications layer
625  finalizeComms();
626 
627  return ret;
628 }
int maxiter_precondition
Definition: quda.h:292
static size_t gSize
double secs
Definition: quda.h:251
void * qdp_longlink[4]
int dimPartitioned(int dim)
Definition: test_util.cpp:1776
QudaDiracFieldOrder dirac_order
Definition: quda.h:219
QudaMassNormalization mass_normalization
Definition: quda.h:208
double tol_hq_offset[QUDA_MAX_MULTI_SHIFT]
Definition: quda.h:182
QudaReconstructType reconstruct_sloppy
Definition: quda.h:53
double anisotropy
Definition: quda.h:38
QudaGhostExchange ghostExchange
Definition: lattice_field.h:76
QudaGaugeParam gaugeParam
Definition: covdev_test.cpp:36
void setPrecision(QudaPrecision precision, QudaPrecision ghost_precision=QUDA_INVALID_PRECISION, bool force_native=false)
int invert_test(void)
void invertMultiShiftQuda(void **_hp_x, void *_hp_b, QudaInvertParam *param)
void endQuda(void)
QudaSolveType solve_type
Definition: quda.h:205
QudaVerbosity verbosity_precondition
Definition: quda.h:286
cpuColorSpinorField * out
enum QudaPrecision_s QudaPrecision
QudaReconstructType link_recon
Definition: test_util.cpp:1605
int ga_pad
Definition: quda.h:63
int xdim
Definition: test_util.cpp:1615
void dw_setDims(int *X, const int L5)
Definition: test_util.cpp:187
QudaGaugeFixed gauge_fix
Definition: quda.h:61
QudaInverterType inv_type_precondition
Definition: quda.h:270
QudaLinkType type
Definition: quda.h:42
void invertQuda(void *h_x, void *h_b, QudaInvertParam *param)
static int X2
Definition: face_gauge.cpp:42
#define errorQuda(...)
Definition: util_quda.h:121
double tol
Definition: quda.h:121
QudaDslashType dslash_type
Definition: quda.h:102
QudaInverterType inv_type
Definition: quda.h:103
void ** ghost_fatlink
QudaPrecision cuda_prec
Definition: quda.h:214
int niter
Definition: test_util.cpp:1629
__host__ __device__ ValueType sqrt(ValueType x)
Definition: complex_quda.h:120
void loadGaugeQuda(void *h_gauge, QudaGaugeParam *param)
int device
Definition: test_util.cpp:1602
QudaPrecision cpu_prec
Definition: quda.h:213
int process_command_line_option(int argc, char **argv, int *idx)
Definition: test_util.cpp:2019
#define mySpinorSiteSize
QudaPrecision prec_sloppy
Definition: test_util.cpp:1609
#define NUM_SRC
QudaDagType dagger
Definition: quda.h:207
void matdagmat_mg4dir(cpuColorSpinorField *out, void **link, void **ghostLink, cpuColorSpinorField *in, int dagger_bit, int mu, QudaPrecision sPrecision, QudaPrecision gPrecision, cpuColorSpinorField *tmp, QudaParity parity)
void finalizeComms()
Definition: test_util.cpp:128
double reliable_delta
Definition: test_util.cpp:1658
cpuColorSpinorField * in
QudaGaugeFieldOrder gauge_order
Definition: quda.h:43
double true_res
Definition: quda.h:126
const char * get_prec_str(QudaPrecision prec)
Definition: misc.cpp:701
QudaSiteSubset siteSubset
Definition: lattice_field.h:71
const char * get_test_type(int t)
Definition: misc.cpp:796
void * qdp_fatlink[4]
void setDims(int *)
Definition: test_util.cpp:151
QudaFieldLocation input_location
Definition: quda.h:99
double true_res_hq_offset[QUDA_MAX_MULTI_SHIFT]
Definition: quda.h:191
double reliable_delta
Definition: quda.h:129
cpuColorSpinorField * tmp
int gridsize_from_cmdline[]
Definition: test_util.cpp:49
int Nsrc
Definition: test_util.cpp:1627
QudaUseInitGuess use_init_guess
Definition: quda.h:231
QudaSolutionType solution_type
Definition: quda.h:204
QudaPrecision cpu_prec
int x[QUDA_MAX_DIM]
Definition: lattice_field.h:67
void display_test_info()
double scale
Definition: quda.h:40
void initQuda(int device)
int tdim
Definition: test_util.cpp:1618
QudaFieldLocation output_location
Definition: quda.h:100
QudaInvertParam inv_param
Definition: covdev_test.cpp:37
void ** ghost_longlink
QudaPrecision cuda_prec_sloppy
Definition: quda.h:215
QudaVerbosity verbosity
Definition: quda.h:244
void setSpinorSiteSize(int n)
Definition: test_util.cpp:211
ColorSpinorParam csParam
Definition: pack_test.cpp:24
double tol_offset[QUDA_MAX_MULTI_SHIFT]
Definition: quda.h:179
double true_res_offset[QUDA_MAX_MULTI_SHIFT]
Definition: quda.h:185
QudaInvertParam newQudaInvertParam(void)
double gflops
Definition: quda.h:250
const char * get_recon_str(QudaReconstructType recon)
Definition: misc.cpp:768
int main(int argc, char **argv)
QudaPrecision cuda_prec_precondition
Definition: quda.h:58
double tol_hq
Definition: quda.h:123
int test_type
Definition: test_util.cpp:1636
cpuColorSpinorField * ref
double true_res_hq
Definition: quda.h:127
void matdagmat(void *out, void **link, void *in, int dagger_bit, int mu, QudaPrecision sPrecision, QudaPrecision gPrecision, void *tmp, QudaParity parity)
QudaGammaBasis gamma_basis
Definition: quda.h:221
QudaPrecision prec
Definition: test_util.cpp:1608
QudaPrecision cuda_prec_sloppy
Definition: quda.h:52
QudaReconstructType link_recon_sloppy
Definition: test_util.cpp:1606
const void ** Ghost() const
Definition: gauge_field.h:323
double tol_precondition
Definition: quda.h:289
double3 HeavyQuarkResidualNorm(ColorSpinorField &x, ColorSpinorField &r)
Definition: reduce_quda.cu:809
double offset[QUDA_MAX_MULTI_SHIFT]
Definition: quda.h:176
int use_sloppy_partial_accumulator
Definition: quda.h:132
enum QudaParity_s QudaParity
QudaInverterType inv_type
Definition: test_util.cpp:1640
QudaReconstructType reconstruct
Definition: quda.h:50
QudaPrecision cuda_prec
Definition: quda.h:49
int X[4]
Definition: quda.h:36
double mass
Definition: quda.h:105
int gcrNkrylov
Definition: quda.h:259
void usage(char **argv)
Definition: test_util.cpp:1783
static int X3
Definition: face_gauge.cpp:42
int V
Definition: test_util.cpp:27
double norm_2(void *v, int len, QudaPrecision precision)
static int X1
Definition: face_gauge.cpp:42
int ydim
Definition: test_util.cpp:1616
cpuGaugeField * cpuFat
void construct_fat_long_gauge_field(void **fatlink, void **longlink, int type, QudaPrecision precision, QudaGaugeParam *param, QudaDslashType dslash_type)
Definition: test_util.cpp:1062
static void end()
QudaResidualType_s
Definition: enum_quda.h:186
static void set_params(QudaGaugeParam *gaugeParam, QudaInvertParam *inv_param, int X1, int X2, int X3, int X4, QudaPrecision cpu_prec, QudaPrecision prec, QudaPrecision prec_sloppy, QudaReconstructType link_recon, QudaReconstructType link_recon_sloppy, double mass, double tol, int maxiter, double reliable_delta, double tadpole_coeff)
double tadpole_coeff
Definition: quda.h:39
QudaPrecision cuda_prec_precondition
Definition: quda.h:217
double tol_restart
Definition: quda.h:122
enum QudaReconstructType_s QudaReconstructType
Main header file for the QUDA library.
void invertMultiSrcQuda(void **_hp_x, void **_hp_b, QudaInvertParam *param)
QudaDslashType dslash_type
Definition: test_util.cpp:1621
__shared__ float s[]
double tol_hq
Definition: test_util.cpp:1657
double tol
Definition: test_util.cpp:1656
cpuGaugeField * cpuLong
#define printfQuda(...)
Definition: util_quda.h:115
QudaTboundary t_boundary
Definition: quda.h:45
void usage_extra(char **argv)
enum QudaDslashType_s QudaDslashType
void * longlink
void mxpy(ColorSpinorField &x, ColorSpinorField &y)
Definition: blas_quda.h:34
QudaResidualType residual_type
Definition: quda.h:320
int num_offset
Definition: quda.h:169
void * fatlink
#define NUM_OFFSETS
#define MAX(a, b)
void initComms(int argc, char **argv, int *const commDims)
Definition: test_util.cpp:88
int zdim
Definition: test_util.cpp:1617
QudaParity parity
Definition: covdev_test.cpp:54
QudaMatPCType matpc_type
Definition: quda.h:206
enum QudaInverterType_s QudaInverterType
QudaPrecision get_prec(QIO_Reader *infile)
Definition: qio_field.cpp:69
QudaPrecision cpu_prec
Definition: quda.h:47
void constructSpinorField(Float *res)
static int X4
Definition: face_gauge.cpp:42
#define gaugeSiteSize
Definition: face_gauge.cpp:34
QudaGaugeParam newQudaGaugeParam(void)
QudaPreserveSource preserve_source
Definition: quda.h:211
int Vh
Definition: test_util.cpp:28