QUDA  1.0.0
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros
quda_milc_interface.h
Go to the documentation of this file.
1 #ifndef _QUDA_MILC_INTERFACE_H
2 #define _QUDA_MILC_INTERFACE_H
3 
4 #include <enum_quda.h>
5 #include <quda.h>
6 
16 #ifdef __cplusplus
17 extern "C" {
18 #endif
19 
23  typedef struct {
24  void *site;
25  void *link;
26  size_t link_offset;
27  void *mom;
28  size_t mom_offset;
29  size_t size;
31 
35  typedef struct {
36  int max_iter;
39  double boundary_phase[4];
43  double tadpole;
46  double naik_epsilon;
48 
53  typedef struct {
55  int nev;
58  double tol_restart;
59 
62  double inc_tol;
63  double eigenval_tol;
64 
67 
70 
71  char *vec_infile;
72  char *vec_outfile;
73 
74  } QudaEigArgs_t;
75 
76 
80  typedef struct {
81  const int* latsize;
82  const int* machsize;
83  int device;
84  } QudaLayout_t;
85 
86 
90  typedef struct {
93  } QudaInitArgs_t; // passed to the initialization struct
94 
95 
99  typedef struct {
104  double force_filter;
106 
107 
111  typedef struct {
115 
121  void qudaSetMPICommHandle(void *mycomm);
122 
128  void qudaInit(QudaInitArgs_t input);
129 
136 
140  void qudaFinalize();
141 
147  void* qudaAllocatePinned(size_t bytes);
148 
153  void qudaFreePinned(void *ptr);
154 
161  void qudaHisqParamsInit(QudaHisqParams_t hisq_params);
162 
175  void qudaLoadKSLink(int precision,
176  QudaFatLinkArgs_t fatlink_args,
177  const double act_path_coeff[6],
178  void* inlink,
179  void* fatlink,
180  void* longlink);
181 
194  void qudaLoadUnitarizedLink(int precision,
195  QudaFatLinkArgs_t fatlink_args,
196  const double path_coeff[6],
197  void* inlink,
198  void* fatlink,
199  void* ulink);
200 
201 
214  void qudaDslash(int external_precision,
215  int quda_precision,
216  QudaInvertArgs_t inv_args,
217  const void* const milc_fatlink,
218  const void* const milc_longlink,
219  void* source,
220  void* solution,
221  int* num_iters);
222 
245  void qudaDDInvert(int external_precision,
246  int quda_precision,
247  double mass,
248  QudaInvertArgs_t inv_args,
249  double target_residual,
250  double target_fermilab_residual,
251  const int * const domain_overlap,
252  const void* const fatlink,
253  const void* const longlink,
254  void* source,
255  void* solution,
256  double* const final_residual,
257  double* const final_fermilab_residual,
258  int* num_iters);
259 
260 
261 
282  void qudaInvert(int external_precision,
283  int quda_precision,
284  double mass,
285  QudaInvertArgs_t inv_args,
286  double target_residual,
287  double target_fermilab_residual,
288  const void* const milc_fatlink,
289  const void* const milc_longlink,
290  void* source,
291  void* solution,
292  double* const final_resid,
293  double* const final_rel_resid,
294  int* num_iters);
295 
317  void qudaInvertMsrc(int external_precision,
318  int quda_precision,
319  double mass,
320  QudaInvertArgs_t inv_args,
321  double target_residual,
322  double target_fermilab_residual,
323  const void* const fatlink,
324  const void* const longlink,
325  void** sourceArray,
326  void** solutionArray,
327  double* const final_residual,
328  double* const final_fermilab_residual,
329  int* num_iters,
330  int num_src);
331 
357  int external_precision,
358  int precision,
359  int num_offsets,
360  double* const offset,
361  QudaInvertArgs_t inv_args,
362  const double* target_residual,
363  const double* target_fermilab_residual,
364  const void* const milc_fatlink,
365  const void* const milc_longlink,
366  void* source,
367  void** solutionArray,
368  double* const final_residual,
369  double* const final_fermilab_residual,
370  int* num_iters);
371 
400  void qudaEigCGInvert(
401  int external_precision,
402  int quda_precision,
403  double mass,
404  QudaInvertArgs_t inv_args,
405  double target_residual,
406  double target_fermilab_residual,
407  const void* const fatlink,
408  const void* const longlink,
409  void* source,
410  void* solution,
411  QudaEigArgs_t eig_args,
412  const int rhs_idx,//current rhs
413  const int last_rhs_flag,//is this the last rhs to solve?
414  double* const final_residual,
415  double* const final_fermilab_residual,
416  int *num_iters);
417 
440  void qudaCloverInvert(int external_precision,
441  int quda_precision,
442  double kappa,
443  double clover_coeff,
444  QudaInvertArgs_t inv_args,
445  double target_residual,
446  double target_fermilab_residual,
447  const void* milc_link,
448  void* milc_clover,
449  void* milc_clover_inv,
450  void* source,
451  void* solution,
452  double* const final_residual,
453  double* const final_fermilab_residual,
454  int* num_iters);
455 
485  int external_precision,
486  int quda_precision,
487  double kappa,
488  double clover_coeff,
489  QudaInvertArgs_t inv_args,
490  double target_residual,
491  double target_fermilab_residual,
492  const void* milc_link,
493  void* milc_clover,
494  void* milc_clover_inv,
495  void* source,
496  void* solution,
497  QudaEigArgs_t eig_args,
498  const int rhs_idx,//current rhs
499  const int last_rhs_flag,//is this the last rhs to solve?
500  double* const final_residual,
501  double* const final_fermilab_residual,
502  int *num_iters);
503 
504 
513  void qudaLoadGaugeField(int external_precision,
514  int quda_precision,
515  QudaInvertArgs_t inv_args,
516  const void* milc_link) ;
517 
521  void qudaFreeGaugeField();
522 
541  void qudaLoadCloverField(int external_precision,
542  int quda_precision,
543  QudaInvertArgs_t inv_args,
544  void* milc_clover,
545  void* milc_clover_inv,
548  double clover_coeff,
549  int compute_trlog,
550  double *trlog) ;
551 
555  void qudaFreeCloverField();
556 
583  void qudaCloverMultishiftInvert(int external_precision,
584  int quda_precision,
585  int num_offsets,
586  double* const offset,
587  double kappa,
588  double clover_coeff,
589  QudaInvertArgs_t inv_args,
590  const double* target_residual,
591  const void* milc_link,
592  void* milc_clover,
593  void* milc_clover_inv,
594  void* source,
595  void** solutionArray,
596  double* const final_residual,
597  int* num_iters
598  );
599 
618  void qudaHisqForce(int precision,
619  int num_terms,
620  int num_naik_terms,
621  double dt,
622  double** coeff,
623  void** quark_field,
624  const double level2_coeff[6],
625  const double fat7_coeff[6],
626  const void* const w_link,
627  const void* const v_link,
628  const void* const u_link,
629  void* const milc_momentum);
630 
642  void qudaGaugeForce(int precision,
643  int num_loop_types,
644  double milc_loop_coeff[3],
645  double eb3,
647 
656  void qudaUpdateU(int precision,
657  double eps,
659 
670  double qudaMomAction(int precision,
671  void *momentum);
672 
684  void qudaRephase(int prec, void *gauge, int flag, double i_mu);
685 
694  void qudaUnitarizeSU3(int prec, double tol, QudaMILCSiteArg_t *arg);
695 
714  void qudaCloverForce(void *mom, double dt, void **x, void **p, double *coeff, double kappa,
715  double ck, int nvec, double multiplicity, void *gauge, int precision,
716  QudaInvertArgs_t inv_args);
717 
729  void qudaCloverTrace(void* out,
730  void* dummy,
731  int mu,
732  int nu);
733 
734 
750  void qudaCloverDerivative(void* out,
751  void* gauge,
752  void* oprod,
753  int mu,
754  int nu,
755  double coeff,
756  int precision,
757  int parity,
758  int conjugate);
759 
760 
770  void* qudaCreateExtendedGaugeField(void* gauge,
771  int geometry,
772  int precision);
773 
783  void* qudaResidentExtendedGaugeField(void* gauge,
784  int geometry,
785  int precision);
786 
795  void* qudaCreateGaugeField(void* gauge,
796  int geometry,
797  int precision);
798 
805  void qudaSaveGaugeField(void* gauge,
806  void* inGauge);
807 
813  void qudaDestroyGaugeField(void* gauge);
814 
815 
828  void qudaGaugeFixingOVR( const int precision,
829  const unsigned int gauge_dir,
830  const int Nsteps,
831  const int verbose_interval,
832  const double relax_boost,
833  const double tolerance,
834  const unsigned int reunit_interval,
835  const unsigned int stopWtheta,
836  void* milc_sitelink
837  );
838 
839 
852  void qudaGaugeFixingFFT( int precision,
853  unsigned int gauge_dir,
854  int Nsteps,
855  int verbose_interval,
856  double alpha,
857  unsigned int autotune,
858  double tolerance,
859  unsigned int stopWtheta,
860  void* milc_sitelink
861  );
862 
863  /* The below declarations are for removed functions from prior versions of QUDA. */
864 
869  void qudaAsqtadForce(int precision,
870  const double act_path_coeff[6],
871  const void* const one_link_src[4],
872  const void* const naik_src[4],
873  const void* const link,
874  void* const milc_momentum);
875 
880  void qudaComputeOprod(int precision,
881  int num_terms,
882  int num_naik_terms,
883  double** coeff,
884  double scale,
885  void** quark_field,
886  void* oprod[3]);
887 
888 #ifdef __cplusplus
889 }
890 #endif
891 
892 
893 #endif // _QUDA_MILC_INTERFACE_H
void qudaHisqParamsInit(QudaHisqParams_t hisq_params)
double mu
Definition: test_util.cpp:1648
void qudaDDInvert(int external_precision, int quda_precision, double mass, QudaInvertArgs_t inv_args, double target_residual, double target_fermilab_residual, const int *const domain_overlap, const void *const fatlink, const void *const longlink, void *source, void *solution, double *const final_residual, double *const final_fermilab_residual, int *num_iters)
enum QudaPrecision_s QudaPrecision
void qudaUnitarizeSU3(int prec, double tol, QudaMILCSiteArg_t *arg)
void * qudaCreateExtendedGaugeField(void *gauge, int geometry, int precision)
void qudaGaugeFixingFFT(int precision, unsigned int gauge_dir, int Nsteps, int verbose_interval, double alpha, unsigned int autotune, double tolerance, unsigned int stopWtheta, void *milc_sitelink)
Gauge fixing with Steepest descent method with FFTs with support for single GPU only.
QudaExtLibType deflation_ext_lib
void qudaDslash(int external_precision, int quda_precision, QudaInvertArgs_t inv_args, const void *const milc_fatlink, const void *const milc_longlink, void *source, void *solution, int *num_iters)
double kappa
Definition: test_util.cpp:1647
QudaPrecision prec_ritz
enum QudaSolveType_s QudaSolveType
void qudaInit(QudaInitArgs_t input)
QudaMemoryType mem_type_ritz
void qudaLoadGaugeField(int external_precision, int quda_precision, QudaInvertArgs_t inv_args, const void *milc_link)
void qudaInvert(int external_precision, int quda_precision, double mass, QudaInvertArgs_t inv_args, double target_residual, double target_fermilab_residual, const void *const milc_fatlink, const void *const milc_longlink, void *source, void *solution, double *const final_resid, double *const final_rel_resid, int *num_iters)
void qudaEigCGInvert(int external_precision, int quda_precision, double mass, QudaInvertArgs_t inv_args, double target_residual, double target_fermilab_residual, const void *const fatlink, const void *const longlink, void *source, void *solution, QudaEigArgs_t eig_args, const int rhs_idx, const int last_rhs_flag, double *const final_residual, double *const final_fermilab_residual, int *num_iters)
void qudaCloverMultishiftInvert(int external_precision, int quda_precision, int num_offsets, double *const offset, double kappa, double clover_coeff, QudaInvertArgs_t inv_args, const double *target_residual, const void *milc_link, void *milc_clover, void *milc_clover_inv, void *source, void **solutionArray, double *const final_residual, int *num_iters)
void qudaEigCGCloverInvert(int external_precision, int quda_precision, double kappa, double clover_coeff, QudaInvertArgs_t inv_args, double target_residual, double target_fermilab_residual, const void *milc_link, void *milc_clover, void *milc_clover_inv, void *source, void *solution, QudaEigArgs_t eig_args, const int rhs_idx, const int last_rhs_flag, double *const final_residual, double *const final_fermilab_residual, int *num_iters)
void qudaGaugeFixingOVR(const int precision, const unsigned int gauge_dir, const int Nsteps, const int verbose_interval, const double relax_boost, const double tolerance, const unsigned int reunit_interval, const unsigned int stopWtheta, void *milc_sitelink)
Gauge fixing with overrelaxation with support for single and multi GPU.
void qudaSaveGaugeField(void *gauge, void *inGauge)
double qudaMomAction(int precision, void *momentum)
void qudaSetLayout(QudaLayout_t layout)
void qudaHisqForce(int precision, int num_terms, int num_naik_terms, double dt, double **coeff, void **quark_field, const double level2_coeff[6], const double fat7_coeff[6], const void *const w_link, const void *const v_link, const void *const u_link, void *const milc_momentum)
QudaSolutionType solution_type
Definition: test_util.cpp:1664
double tol
Definition: test_util.cpp:1656
void qudaFreePinned(void *ptr)
void qudaUpdateU(int precision, double eps, QudaMILCSiteArg_t *arg)
QIO_Layout layout
Definition: qio_field.cpp:6
void qudaFreeCloverField()
void qudaInvertMsrc(int external_precision, int quda_precision, double mass, QudaInvertArgs_t inv_args, double target_residual, double target_fermilab_residual, const void *const fatlink, const void *const longlink, void **sourceArray, void **solutionArray, double *const final_residual, double *const final_fermilab_residual, int *num_iters, int num_src)
int nvec[QUDA_MAX_MG_LEVEL]
Definition: test_util.cpp:1637
void * qudaCreateGaugeField(void *gauge, int geometry, int precision)
QudaInverterType solver_type
void qudaRephase(int prec, void *gauge, int flag, double i_mu)
enum QudaSolutionType_s QudaSolutionType
void qudaComputeOprod(int precision, int num_terms, int num_naik_terms, double **coeff, double scale, void **quark_field, void *oprod[3])
const int * machsize
void qudaFreeGaugeField()
enum QudaParity_s QudaParity
void qudaLoadCloverField(int external_precision, int quda_precision, QudaInvertArgs_t inv_args, void *milc_clover, void *milc_clover_inv, QudaSolutionType solution_type, QudaSolveType solve_type, double clover_coeff, int compute_trlog, double *trlog)
void qudaFinalize()
QudaFieldLocation location_ritz
void qudaCloverInvert(int external_precision, int quda_precision, double kappa, double clover_coeff, QudaInvertArgs_t inv_args, double target_residual, double target_fermilab_residual, const void *milc_link, void *milc_clover, void *milc_clover_inv, void *source, void *solution, double *const final_residual, double *const final_fermilab_residual, int *num_iters)
QudaSolveType solve_type
Definition: test_util.cpp:1663
void qudaCloverDerivative(void *out, void *gauge, void *oprod, int mu, int nu, double coeff, int precision, int parity, int conjugate)
void qudaLoadUnitarizedLink(int precision, QudaFatLinkArgs_t fatlink_args, const double path_coeff[6], void *inlink, void *fatlink, void *ulink)
enum QudaFieldLocation_s QudaFieldLocation
cpuColorSpinorField * out
void qudaSetMPICommHandle(void *mycomm)
double clover_coeff
Definition: test_util.cpp:1653
Main header file for the QUDA library.
void * qudaResidentExtendedGaugeField(void *gauge, int geometry, int precision)
void qudaCloverTrace(void *out, void *dummy, int mu, int nu)
const int * latsize
void qudaCloverForce(void *mom, double dt, void **x, void **p, double *coeff, double kappa, double ck, int nvec, double multiplicity, void *gauge, int precision, QudaInvertArgs_t inv_args)
void * qudaAllocatePinned(size_t bytes)
__host__ __device__ ValueType arg(const complex< ValueType > &z)
Returns the phase angle of z.
void * longlink
enum QudaVerbosity_s QudaVerbosity
void * fatlink
QudaVerbosity verbosity
void qudaLoadKSLink(int precision, QudaFatLinkArgs_t fatlink_args, const double act_path_coeff[6], void *inlink, void *fatlink, void *longlink)
QudaExtLibType solver_ext_lib
QudaParity parity
Definition: covdev_test.cpp:54
void qudaGaugeForce(int precision, int num_loop_types, double milc_loop_coeff[3], double eb3, QudaMILCSiteArg_t *arg)
QudaPrecision prec
Definition: test_util.cpp:1608
enum QudaInverterType_s QudaInverterType
void qudaMultishiftInvert(int external_precision, int precision, int num_offsets, double *const offset, QudaInvertArgs_t inv_args, const double *target_residual, const double *target_fermilab_residual, const void *const milc_fatlink, const void *const milc_longlink, void *source, void **solutionArray, double *const final_residual, double *const final_fermilab_residual, int *num_iters)
enum QudaMemoryType_s QudaMemoryType
void qudaAsqtadForce(int precision, const double act_path_coeff[6], const void *const one_link_src[4], const void *const naik_src[4], const void *const link, void *const milc_momentum)
unsigned long long bytes
Definition: blas_quda.cu:23
enum QudaExtLibType_s QudaExtLibType
QudaLayout_t layout
void qudaDestroyGaugeField(void *gauge)