QUDA  1.0.0
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros
Classes | Macros | Typedefs | Functions
quda.h File Reference

Main header file for the QUDA library. More...

#include <enum_quda.h>
#include <stdio.h>
#include <quda_define.h>
#include <quda_constants.h>
Include dependency graph for quda.h:

Go to the source code of this file.

Classes

struct  QudaGaugeParam_s
 
struct  QudaInvertParam_s
 
struct  QudaEigParam_s
 
struct  QudaMultigridParam_s
 

Macros

#define double_complex   double _Complex
 

Typedefs

typedef struct QudaGaugeParam_s QudaGaugeParam
 
typedef struct QudaInvertParam_s QudaInvertParam
 
typedef struct QudaEigParam_s QudaEigParam
 
typedef struct QudaMultigridParam_s QudaMultigridParam
 
typedef int(* QudaCommsMap) (const int *coords, void *fdata)
 

Functions

void setVerbosityQuda (QudaVerbosity verbosity, const char prefix[], FILE *outfile)
 
void qudaSetCommHandle (void *mycomm)
 
void initCommsGridQuda (int nDim, const int *dims, QudaCommsMap func, void *fdata)
 
void initQudaDevice (int device)
 
void initQudaMemory ()
 
void initQuda (int device)
 
void endQuda (void)
 
void updateR ()
 update the radius for halos. More...
 
QudaGaugeParam newQudaGaugeParam (void)
 
QudaInvertParam newQudaInvertParam (void)
 
QudaMultigridParam newQudaMultigridParam (void)
 
QudaEigParam newQudaEigParam (void)
 
void printQudaGaugeParam (QudaGaugeParam *param)
 
void printQudaInvertParam (QudaInvertParam *param)
 
void printQudaMultigridParam (QudaMultigridParam *param)
 
void printQudaEigParam (QudaEigParam *param)
 
void loadGaugeQuda (void *h_gauge, QudaGaugeParam *param)
 
void freeGaugeQuda (void)
 
void saveGaugeQuda (void *h_gauge, QudaGaugeParam *param)
 
void loadCloverQuda (void *h_clover, void *h_clovinv, QudaInvertParam *inv_param)
 
void freeCloverQuda (void)
 
void lanczosQuda (int k0, int m, void *hp_Apsi, void *hp_r, void *hp_V, void *hp_alpha, void *hp_beta, QudaEigParam *eig_param)
 
void eigensolveQuda (void **h_evecs, double_complex *h_evals, QudaEigParam *param)
 
void invertQuda (void *h_x, void *h_b, QudaInvertParam *param)
 
void invertMultiSrcQuda (void **_hp_x, void **_hp_b, QudaInvertParam *param)
 
void invertMultiShiftQuda (void **_hp_x, void *_hp_b, QudaInvertParam *param)
 
void * newMultigridQuda (QudaMultigridParam *param)
 
void destroyMultigridQuda (void *mg_instance)
 Free resources allocated by the multigrid solver. More...
 
void updateMultigridQuda (void *mg_instance, QudaMultigridParam *param)
 Updates the multigrid preconditioner for the new gauge / clover field. More...
 
void dumpMultigridQuda (void *mg_instance, QudaMultigridParam *param)
 Dump the null-space vectors to disk. More...
 
void dslashQuda (void *h_out, void *h_in, QudaInvertParam *inv_param, QudaParity parity)
 
void dslashQuda_4dpc (void *h_out, void *h_in, QudaInvertParam *inv_param, QudaParity parity, int test_type)
 
void dslashQuda_mdwf (void *h_out, void *h_in, QudaInvertParam *inv_param, QudaParity parity, int test_type)
 
void cloverQuda (void *h_out, void *h_in, QudaInvertParam *inv_param, QudaParity *parity, int inverse)
 
void MatQuda (void *h_out, void *h_in, QudaInvertParam *inv_param)
 
void MatDagMatQuda (void *h_out, void *h_in, QudaInvertParam *inv_param)
 
void set_dim (int *)
 
void pack_ghost (void **cpuLink, void **cpuGhost, int nFace, QudaPrecision precision)
 
void computeKSLinkQuda (void *fatlink, void *longlink, void *ulink, void *inlink, double *path_coeff, QudaGaugeParam *param)
 
int computeGaugeForceQuda (void *mom, void *sitelink, int ***input_path_buf, int *path_length, double *loop_coeff, int num_paths, int max_length, double dt, QudaGaugeParam *qudaGaugeParam)
 
void updateGaugeFieldQuda (void *gauge, void *momentum, double dt, int conj_mom, int exact, QudaGaugeParam *param)
 
void staggeredPhaseQuda (void *gauge_h, QudaGaugeParam *param)
 
void projectSU3Quda (void *gauge_h, double tol, QudaGaugeParam *param)
 
double momActionQuda (void *momentum, QudaGaugeParam *param)
 
void * createGaugeFieldQuda (void *gauge, int geometry, QudaGaugeParam *param)
 
void saveGaugeFieldQuda (void *outGauge, void *inGauge, QudaGaugeParam *param)
 
void destroyGaugeFieldQuda (void *gauge)
 
void createCloverQuda (QudaInvertParam *param)
 
void computeCloverForceQuda (void *mom, double dt, void **x, void **p, double *coeff, double kappa2, double ck, int nvector, double multiplicity, void *gauge, QudaGaugeParam *gauge_param, QudaInvertParam *inv_param)
 
void computeStaggeredForceQuda (void *mom, double dt, double delta, void **x, void *gauge, QudaGaugeParam *gauge_param, QudaInvertParam *invert_param)
 
void computeHISQForceQuda (void *momentum, double dt, 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 **quark, int num, int num_naik, double **coeff, QudaGaugeParam *param)
 
void gaussGaugeQuda (unsigned long long seed, double sigma)
 Generate Gaussian distributed fields and store in the resident gauge field. We create a Gaussian-distributed su(n) field and exponentiate it, e.g., U = exp(sigma * H), where H is the distributed su(n) field and beta is the width of the distribution (beta = 0 results in a free field, and sigma = 1 has maximum disorder). More...
 
void plaqQuda (double plaq[3])
 
void copyExtendedResidentGaugeQuda (void *resident_gauge, QudaFieldLocation loc)
 
void performWuppertalnStep (void *h_out, void *h_in, QudaInvertParam *param, unsigned int nSteps, double alpha)
 
void performAPEnStep (unsigned int nSteps, double alpha)
 
void performSTOUTnStep (unsigned int nSteps, double rho)
 
void performOvrImpSTOUTnStep (unsigned int nSteps, double rho, double epsilon)
 
double qChargeQuda ()
 
void contractQuda (const void *x, const void *y, void *result, const QudaContractType cType, QudaInvertParam *param, const int *X)
 
double qChargeDensityQuda (void *qDensity)
 Calculates the topological charge from gaugeSmeared, if it exist, or from gaugePrecise if no smeared fields are present. More...
 
int computeGaugeFixingOVRQuda (void *gauge, const unsigned int gauge_dir, const unsigned int Nsteps, const unsigned int verbose_interval, const double relax_boost, const double tolerance, const unsigned int reunit_interval, const unsigned int stopWtheta, QudaGaugeParam *param, double *timeinfo)
 Gauge fixing with overrelaxation with support for single and multi GPU. More...
 
int computeGaugeFixingFFTQuda (void *gauge, const unsigned int gauge_dir, const unsigned int Nsteps, const unsigned int verbose_interval, const double alpha, const unsigned int autotune, const double tolerance, const unsigned int stopWtheta, QudaGaugeParam *param, double *timeinfo)
 Gauge fixing with Steepest descent method with FFTs with support for single GPU only. More...
 
void flushChronoQuda (int index)
 Flush the chronological history for the given index. More...
 
void openMagma ()
 
void closeMagma ()
 
void * newDeflationQuda (QudaEigParam *param)
 
void destroyDeflationQuda (void *df_instance)
 
void setMPICommHandleQuda (void *mycomm)
 

Detailed Description

Main header file for the QUDA library.

Note to QUDA developers: When adding new members to QudaGaugeParam and QudaInvertParam, be sure to update lib/check_params.h as well as the Fortran interface in lib/quda_fortran.F90.

Definition in file quda.h.

Macro Definition Documentation

◆ double_complex

#define double_complex   double _Complex

Definition at line 19 of file quda.h.

Typedef Documentation

◆ QudaCommsMap

typedef int(* QudaCommsMap) (const int *coords, void *fdata)

initCommsGridQuda() takes an optional "rank_from_coords" argument that should be a pointer to a user-defined function with this prototype.

Parameters
coordsNode coordinates
fdataAny auxiliary data needed by the function
Returns
MPI rank or QMP node ID cooresponding to the node coordinates
See also
initCommsGridQuda

Definition at line 696 of file quda.h.

◆ QudaEigParam

typedef struct QudaEigParam_s QudaEigParam

◆ QudaGaugeParam

Parameters having to do with the gauge field or the interpretation of the gauge field by various Dirac operators

◆ QudaInvertParam

Parameters relating to the solver and the choice of Dirac operator.

◆ QudaMultigridParam

Function Documentation

◆ closeMagma()

void closeMagma ( )

Definition at line 106 of file interface_quda.cpp.

References CloseMagma(), InitMagma, and printfQuda.

Referenced by destroyDeflationQuda().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ cloverQuda()

void cloverQuda ( void *  h_out,
void *  h_in,
QudaInvertParam inv_param,
QudaParity parity,
int  inverse 
)

Apply the clover operator or its inverse.

Parameters
h_outResult spinor field
h_inInput spinor field
paramContains all metadata regarding host and device storage
parityThe source and destination parity of the field
inverseWhether to apply the inverse of the clover term

◆ computeCloverForceQuda()

void computeCloverForceQuda ( void *  mom,
double  dt,
void **  x,
void **  p,
double *  coeff,
double  kappa2,
double  ck,
int  nvector,
double  multiplicity,
void *  gauge,
QudaGaugeParam gauge_param,
QudaInvertParam inv_param 
)

Compute the clover force contributions in each dimension mu given the array of solution fields, and compute the resulting momentum field.

Parameters
momForce matrix
dtIntegrating step size
xArray of solution vectors
pArray of intermediate vectors
coeffArray of residues for each contribution (multiplied by stepsize)
kappa2-kappa*kappa parameter
ck-clover_coefficient * kappa / 8
nvecNumber of vectors
multiplicityNumber fermions this bilinear reresents
gaugeGauge Field
gauge_paramGauge field meta data
inv_paramDirac and solver meta data

Definition at line 4684 of file interface_quda.cpp.

References checkCudaError, quda::cloverDerivative(), quda::computeCloverForce(), quda::computeCloverSigmaOprod(), quda::computeCloverSigmaTrace(), quda::copy(), cpuMom, quda::GaugeFieldParam::create, quda::ColorSpinorParam::create, quda::Dirac::create(), quda::ColorSpinorField::Create(), createExtendedGauge(), cudaForce, cudaMom, quda::Dirac::Dagger(), dirac, quda::Dirac::Dslash(), errorQuda, quda::ColorSpinorField::Even(), extendedGaugeResident, quda::ColorSpinorParam::fieldOrder, quda::gamma5(), quda::ColorSpinorParam::gammaBasis, QudaGaugeParam_s::gauge_order, quda::GaugeFieldParam::geometry, quda::GaugeFieldParam::link_type, quda::ColorSpinorParam::location, quda::Dirac::M(), quda::ColorSpinorParam::nColor, quda::LatticeFieldParam::nDim, quda::ColorSpinorParam::nSpin, quda::ColorSpinorField::Odd(), quda::GaugeFieldParam::order, quda::LatticeFieldParam::pad, param, quda::LatticeFieldParam::Precision(), profileCloverForce, QUDA_ASQTAD_MOM_LINKS, QUDA_CUDA_FIELD_LOCATION, QUDA_DAG_NO, QUDA_DAG_YES, QUDA_DEGRAND_ROSSI_GAMMA_BASIS, QUDA_DIRECT_PC_SOLVE, QUDA_DOUBLE_PRECISION, QUDA_EVEN_ODD_SITE_ORDER, QUDA_EVEN_PARITY, QUDA_FLOAT2_FIELD_ORDER, QUDA_FLOAT2_GAUGE_ORDER, QUDA_FULL_SITE_SUBSET, QUDA_GENERAL_LINKS, QUDA_NORMOP_PC_SOLVE, QUDA_NULL_FIELD_CREATE, QUDA_ODD_PARITY, QUDA_PARITY_SITE_SUBSET, quda::QUDA_PROFILE_COMPUTE, quda::QUDA_PROFILE_D2H, quda::QUDA_PROFILE_FREE, quda::QUDA_PROFILE_H2D, quda::QUDA_PROFILE_INIT, quda::QUDA_PROFILE_TOTAL, QUDA_RECONSTRUCT_10, QUDA_RECONSTRUCT_12, QUDA_RECONSTRUCT_NO, QUDA_REFERENCE_FIELD_CREATE, QUDA_SPACE_SPIN_COLOR_FIELD_ORDER, QUDA_TENSOR_GEOMETRY, QUDA_UKQCD_GAMMA_BASIS, QUDA_ZERO_FIELD_CREATE, R, quda::GaugeFieldParam::reconstruct, quda::cudaGaugeField::saveCPUField(), quda::setDiracParam(), quda::ColorSpinorParam::setPrecision(), quda::ColorSpinorParam::siteOrder, quda::LatticeFieldParam::siteSubset, QudaInvertParam_s::solve_type, tmp, quda::DiracParam::tmp1, quda::updateMomentum(), QudaInvertParam_s::use_resident_solution, quda::ColorSpinorParam::v, and quda::LatticeFieldParam::x.

Here is the call graph for this function:

◆ computeGaugeFixingFFTQuda()

int computeGaugeFixingFFTQuda ( void *  gauge,
const unsigned int  gauge_dir,
const unsigned int  Nsteps,
const unsigned int  verbose_interval,
const double  alpha,
const unsigned int  autotune,
const double  tolerance,
const unsigned int  stopWtheta,
QudaGaugeParam param,
double *  timeinfo 
)

Gauge fixing with Steepest descent method with FFTs with support for single GPU only.

Parameters
[in,out]gauge,gaugefield to be fixed
[in]gauge_dir,3for Coulomb gauge fixing, other for Landau gauge fixing
[in]Nsteps,maximumnumber of steps to perform gauge fixing
[in]verbose_interval,printgauge fixing info when iteration count is a multiple of this
[in]alpha,gaugefixing parameter of the method, most common value is 0.08
[in]autotune,1to autotune the method, i.e., if the Fg inverts its tendency we decrease the alpha value
[in]tolerance,torelancevalue to stop the method, if this value is zero then the method stops when iteration reachs the maximum number of steps defined by Nsteps
[in]stopWtheta,0for MILC criterium and 1 to use the theta value
[in]paramThe parameters of the external fields and the computation settings
[out]timeinfo

Definition at line 5716 of file interface_quda.cpp.

References checkCudaError, cpuGauge, quda::GaugeFieldParam::create, GaugeFixFFTQuda, quda::gaugefixingFFT(), gaugePrecise, gParam, quda::TimeProfile::Last(), quda::GaugeFieldParam::link_type, QudaGaugeParam_s::make_resident_gauge, quda::GaugeFieldParam::order, quda::LatticeFieldParam::Precision(), QUDA_DOUBLE_PRECISION, QUDA_FLOAT2_GAUGE_ORDER, QUDA_FLOAT4_GAUGE_ORDER, QUDA_NULL_FIELD_CREATE, quda::QUDA_PROFILE_COMPUTE, quda::QUDA_PROFILE_D2H, quda::QUDA_PROFILE_H2D, quda::QUDA_PROFILE_INIT, quda::QUDA_PROFILE_TOTAL, QUDA_RECONSTRUCT_NO, quda::GaugeFieldParam::reconstruct, QudaGaugeParam_s::reconstruct, and QudaGaugeParam_s::type.

Here is the call graph for this function:

◆ computeGaugeFixingOVRQuda()

int computeGaugeFixingOVRQuda ( void *  gauge,
const unsigned int  gauge_dir,
const unsigned int  Nsteps,
const unsigned int  verbose_interval,
const double  relax_boost,
const double  tolerance,
const unsigned int  reunit_interval,
const unsigned int  stopWtheta,
QudaGaugeParam param,
double *  timeinfo 
)

Gauge fixing with overrelaxation with support for single and multi GPU.

Parameters
[in,out]gauge,gaugefield to be fixed
[in]gauge_dir,3for Coulomb gauge fixing, other for Landau gauge fixing
[in]Nsteps,maximumnumber of steps to perform gauge fixing
[in]verbose_interval,printgauge fixing info when iteration count is a multiple of this
[in]relax_boost,gaugefixing parameter of the overrelaxation method, most common value is 1.5 or 1.7.
[in]tolerance,torelancevalue to stop the method, if this value is zero then the method stops when iteration reachs the maximum number of steps defined by Nsteps
[in]reunit_interval,reunitarizegauge field when iteration count is a multiple of this
[in]stopWtheta,0for MILC criterium and 1 to use the theta value
[in]paramThe parameters of the external fields and the computation settings
[out]timeinfo

if (!param->use_resident_gauge) { // load fields onto the device

Definition at line 5634 of file interface_quda.cpp.

References checkCudaError, comm_size(), quda::copyExtendedGauge(), cpuGauge, quda::GaugeFieldParam::create, createExtendedGauge(), quda::gaugefixingOVR(), GaugeFixOVRQuda, gaugePrecise, gParam, quda::TimeProfile::Last(), quda::GaugeFieldParam::link_type, QudaGaugeParam_s::make_resident_gauge, quda::GaugeFieldParam::order, quda::LatticeFieldParam::Precision(), QUDA_CUDA_FIELD_LOCATION, QUDA_DOUBLE_PRECISION, QUDA_FLOAT2_GAUGE_ORDER, QUDA_FLOAT4_GAUGE_ORDER, QUDA_NULL_FIELD_CREATE, quda::QUDA_PROFILE_COMPUTE, quda::QUDA_PROFILE_D2H, quda::QUDA_PROFILE_H2D, quda::QUDA_PROFILE_INIT, quda::QUDA_PROFILE_TOTAL, QUDA_RECONSTRUCT_NO, R, quda::GaugeFieldParam::reconstruct, QudaGaugeParam_s::reconstruct, and QudaGaugeParam_s::type.

Here is the call graph for this function:

◆ computeGaugeForceQuda()

int computeGaugeForceQuda ( void *  mom,
void *  sitelink,
int ***  input_path_buf,
int *  path_length,
double *  loop_coeff,
int  num_paths,
int  max_length,
double  dt,
QudaGaugeParam qudaGaugeParam 
)

Compute the gauge force and update the mometum field

Parameters
momThe momentum field to be updated
sitelinkThe gauge field from which we compute the force
input_path_buf[dim][num_paths][path_length]
path_lengthOne less that the number of links in a loop (e.g., 3 for a staple)
loop_coeffCoefficients of the different loops in the Symanzik action
num_pathsHow many contributions from path_length different "staples"
max_lengthThe maximum number of non-zero of links in any path in the action
dtThe integration step size (for MILC this is dt*beta/3)
paramThe parameters of the external fields and the computation settings

Definition at line 4073 of file interface_quda.cpp.

References checkCudaError, cpuMom, quda::GaugeFieldParam::create, quda::GaugeField::Create(), createExtendedGauge(), QudaGaugeParam_s::cuda_prec, cudaGauge, cudaMom, errorQuda, extendedGaugeResident, quda::forceMonitor(), QudaGaugeParam_s::gauge_offset, quda::gaugeForce(), gaugePrecise, gParam, quda::GaugeFieldParam::link_type, quda::cudaGaugeField::loadCPUField(), QudaGaugeParam_s::make_resident_gauge, QudaGaugeParam_s::make_resident_mom, QudaGaugeParam_s::mom_offset, momResident, quda::GaugeFieldParam::order, QudaGaugeParam_s::overwrite_mom, profileGaugeForce, QUDA_ASQTAD_MOM_LINKS, QUDA_DOUBLE_PRECISION, QUDA_FLOAT2_GAUGE_ORDER, QUDA_FLOAT4_GAUGE_ORDER, QUDA_MILC_GAUGE_ORDER, QUDA_NULL_FIELD_CREATE, quda::QUDA_PROFILE_COMPUTE, quda::QUDA_PROFILE_D2H, quda::QUDA_PROFILE_FREE, quda::QUDA_PROFILE_H2D, quda::QUDA_PROFILE_INIT, quda::QUDA_PROFILE_TOTAL, QUDA_QDP_GAUGE_ORDER, QUDA_RECONSTRUCT_10, QUDA_RECONSTRUCT_NO, QUDA_TIFR_GAUGE_ORDER, QUDA_TIFR_PADDED_GAUGE_ORDER, QUDA_ZERO_FIELD_CREATE, R, quda::GaugeFieldParam::reconstruct, QudaGaugeParam_s::reconstruct, QudaGaugeParam_s::return_result_mom, quda::cudaGaugeField::saveCPUField(), quda::GaugeFieldParam::setPrecision(), quda::GaugeFieldParam::site_offset, quda::GaugeFieldParam::site_size, QudaGaugeParam_s::site_size, quda::updateMomentum(), QudaGaugeParam_s::use_resident_gauge, QudaGaugeParam_s::use_resident_mom, and quda::cudaGaugeField::zero().

Referenced by compute_gauge_force_quda_(), and gauge_force_test().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ computeHISQForceQuda()

void computeHISQForceQuda ( void *  momentum,
double  dt,
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 **  quark,
int  num,
int  num_naik,
double **  coeff,
QudaGaugeParam param 
)

Compute the fermion force for the HISQ quark action and integrate the momentum.

Parameters
momentumThe momentum field we are integrating
dtThe stepsize used to integrate the momentum
level2_coeffThe coefficients for the second level of smearing in the quark action.
fat7_coeffThe coefficients for the first level of smearing (fat7) in the quark action.
w_linkUnitarized link variables obtained by applying fat7 smearing and unitarization to the original links.
v_linkFat7 link variables.
u_linkSU(3) think link variables.
quarkThe input fermion field.
numThe number of quark fields
num_naikThe number of naik contributions
coeffThe coefficient multiplying the fermion fields in the outer product
param.The field parameters.

Definition at line 4433 of file interface_quda.cpp.

References quda::ax(), quda::GaugeField::Bytes(), comm_dim_partitioned(), quda::computeStaggeredOprod(), quda::cudaGaugeField::copy(), quda::copyExtendedGauge(), QudaGaugeParam_s::cpu_prec, cpuMom, cpuULink, quda::GaugeFieldParam::create, quda::ColorSpinorParam::create, cudaGauge, cudaMom, errorQuda, quda::cudaGaugeField::exchangeExtendedGhost(), quda::ColorSpinorParam::fieldOrder, quda::GaugeFieldParam::gauge, QudaGaugeParam_s::gauge_order, quda::cudaGaugeField::Gauge_p(), quda::LatticeFieldParam::ghostExchange, quda::fermion_force::hisqCompleteForce(), quda::fermion_force::hisqLongLinkForce(), quda::fermion_force::hisqStaplesForce(), quda::GaugeFieldParam::link_type, quda::cudaGaugeField::loadCPUField(), QudaGaugeParam_s::make_resident_mom, momResident, quda::ColorSpinorParam::nColor, quda::LatticeFieldParam::nDim, quda::GaugeFieldParam::nFace, quda::ColorSpinorParam::nSpin, num_failures_d, num_failures_h, quda::GaugeFieldParam::order, quda::LatticeFieldParam::pad, param, quda::LatticeFieldParam::Precision(), profileHISQForce, QUDA_ASQTAD_MOM_LINKS, QUDA_CUDA_FIELD_LOCATION, QUDA_EVEN_ODD_SITE_ORDER, QUDA_FLOAT2_FIELD_ORDER, QUDA_FLOAT2_GAUGE_ORDER, QUDA_FULL_SITE_SUBSET, QUDA_GENERAL_LINKS, QUDA_GHOST_EXCHANGE_EXTENDED, QUDA_GHOST_EXCHANGE_NO, QUDA_MILC_GAUGE_ORDER, QUDA_NULL_FIELD_CREATE, quda::QUDA_PROFILE_COMPUTE, quda::QUDA_PROFILE_FREE, quda::QUDA_PROFILE_H2D, quda::QUDA_PROFILE_INIT, quda::QUDA_PROFILE_TOTAL, QUDA_RECONSTRUCT_10, QUDA_RECONSTRUCT_NO, QUDA_REFERENCE_FIELD_CREATE, QUDA_SPACE_COLOR_SPIN_FIELD_ORDER, QUDA_ZERO_FIELD_CREATE, quda::LatticeFieldParam::r, R, quda::GaugeFieldParam::reconstruct, QudaGaugeParam_s::return_result_mom, quda::cudaGaugeField::saveCPUField(), quda::GaugeFieldParam::setPrecision(), quda::ColorSpinorParam::setPrecision(), quda::fermion_force::setUnitarizeForceConstants(), quda::ColorSpinorParam::siteOrder, quda::LatticeFieldParam::siteSubset, unitarize_eps, quda::fermion_force::unitarizeForce(), quda::updateMomentum(), QudaGaugeParam_s::use_resident_mom, quda::ColorSpinorParam::v, and quda::LatticeFieldParam::x.

Here is the call graph for this function:

◆ computeKSLinkQuda()

void computeKSLinkQuda ( void *  fatlink,
void *  longlink,
void *  ulink,
void *  inlink,
double *  path_coeff,
QudaGaugeParam param 
)

◆ computeStaggeredForceQuda()

void computeStaggeredForceQuda ( void *  mom,
double  dt,
double  delta,
void **  x,
void *  gauge,
QudaGaugeParam gauge_param,
QudaInvertParam invert_param 
)

Compute the naive staggered force. All fields must be in the same precision.

Parameters
momMomentum field
dtIntegrating step size
deltaAdditional scale factor when updating momentum (mom += delta * [force]_TA
gaugeGauge field (at present only supports resident gauge field)
xArray of single-parity solution vectors (at present only supports resident solutions)
gauge_paramGauge field meta data
invert_paramDirac and solver meta data

◆ contractQuda()

void contractQuda ( const void *  x,
const void *  y,
void *  result,
const QudaContractType  cType,
QudaInvertParam param,
const int *  X 
)

Public function to perform color contractions of the host spinors x and y.

Parameters
[in]xpointer to host data
[in]ypointer to host data
[out]resultpointer to the 16 spin projections per lattice site
[in]cTypeWhich type of contraction (open, degrand-rossi, etc)
[in]parammeta data for construction of ColorSpinorFields.
[in]Xspacetime data for construction of ColorSpinorFields.

Definition at line 5790 of file interface_quda.cpp.

References quda::contractQuda(), quda::ColorSpinorParam::create, quda::ColorSpinorField::Create(), quda::ColorSpinorParam::gammaBasis, QudaInvertParam_s::input_location, quda::ColorSpinorParam::location, pool_device_free, pool_device_malloc, quda::LatticeFieldParam::Precision(), profileContract, QUDA_CUDA_FIELD_LOCATION, QUDA_DEGRAND_ROSSI_GAMMA_BASIS, QUDA_NULL_FIELD_CREATE, quda::QUDA_PROFILE_COMPUTE, quda::QUDA_PROFILE_D2H, quda::QUDA_PROFILE_FREE, quda::QUDA_PROFILE_H2D, quda::QUDA_PROFILE_INIT, quda::QUDA_PROFILE_TOTAL, qudaMemcpy, quda::ColorSpinorParam::setPrecision(), and quda::ColorSpinorParam::v.

Here is the call graph for this function:

◆ copyExtendedResidentGaugeQuda()

void copyExtendedResidentGaugeQuda ( void *  resident_gauge,
QudaFieldLocation  loc 
)

Definition at line 5441 of file interface_quda.cpp.

References quda::copyExtendedGauge(), createExtendedGauge(), errorQuda, extendedGaugeResident, profilePlaq, and R.

Referenced by main().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ createCloverQuda()

void createCloverQuda ( QudaInvertParam param)

◆ createGaugeFieldQuda()

void* createGaugeFieldQuda ( void *  gauge,
int  geometry,
QudaGaugeParam param 
)

Allocate a gauge (matrix) field on the device and optionally download a host gauge field.

Parameters
gaugeThe host gauge field (optional - if set to 0 then the gauge field zeroed)
geometryThe geometry of the matrix field to create (1 - scalar, 4 - vector, 6 - tensor)
paramThe parameters of the external field and the field to be created
Returns
Pointer to the gauge field (cast as a void*)

Definition at line 4229 of file interface_quda.cpp.

References cpuGauge, quda::GaugeFieldParam::create, cudaGauge, errorQuda, quda::GaugeFieldParam::geometry, gParam, quda::cudaGaugeField::loadCPUField(), quda::GaugeFieldParam::order, QUDA_FLOAT2_GAUGE_ORDER, QUDA_GENERAL_LINKS, QUDA_SCALAR_GEOMETRY, QUDA_VECTOR_GEOMETRY, and QUDA_ZERO_FIELD_CREATE.

Here is the call graph for this function:

◆ destroyDeflationQuda()

void destroyDeflationQuda ( void *  df_instance)

Free resources allocated by the deflated solver

Definition at line 2823 of file interface_quda.cpp.

References closeMagma().

Referenced by main().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ destroyGaugeFieldQuda()

void destroyGaugeFieldQuda ( void *  gauge)

Reinterpret gauge as a pointer to cudaGaugeField and call destructor.

Parameters
gaugeGauge field to be freed

Definition at line 4265 of file interface_quda.cpp.

◆ destroyMultigridQuda()

void destroyMultigridQuda ( void *  mg_instance)

Free resources allocated by the multigrid solver.

Parameters
mg_instancePointer to instance of multigrid_solver
paramContains all metadata regarding host and device storage and solver parameters

Definition at line 2641 of file interface_quda.cpp.

References quda::multigrid_solver::mg.

Referenced by main().

Here is the caller graph for this function:

◆ dslashQuda()

void dslashQuda ( void *  h_out,
void *  h_in,
QudaInvertParam inv_param,
QudaParity  parity 
)

◆ dslashQuda_4dpc()

void dslashQuda_4dpc ( void *  h_out,
void *  h_in,
QudaInvertParam inv_param,
QudaParity  parity,
int  test_type 
)

Apply the Dslash operator (D_{eo} or D_{oe}) for 4D EO preconditioned DWF.

Parameters
h_outResult spinor field
h_inInput spinor field
paramContains all metadata regarding host and device storage
parityThe destination parity of the field
test_typeChoose a type of dslash operators

Definition at line 1945 of file interface_quda.cpp.

References quda::blas::ax(), quda::ColorSpinorParam::create, quda::ColorSpinorField::Create(), dirac, QudaInvertParam_s::dirac_order, quda::DiracDomainWall4D::Dslash4(), quda::DiracDomainWall4D::Dslash5(), quda::DiracDomainWall4DPC::Dslash5inv(), QudaInvertParam_s::dslash_type, errorQuda, gaugeLongPrecise, getVerbosity(), in, QudaInvertParam_s::input_location, QudaInvertParam_s::kappa, quda::blas::norm2(), out, QudaInvertParam_s::output_location, popVerbosity(), printfQuda, printQudaInvertParam(), pushVerbosity(), QUDA_ASQTAD_DSLASH, QUDA_CPS_WILSON_DIRAC_ORDER, QUDA_DEBUG_VERBOSE, QUDA_EVEN_PARITY, QUDA_NULL_FIELD_CREATE, QUDA_ODD_PARITY, quda::setDiracParam(), and QudaInvertParam_s::verbosity.

Referenced by dslashCUDA().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ dslashQuda_mdwf()

void dslashQuda_mdwf ( void *  h_out,
void *  h_in,
QudaInvertParam inv_param,
QudaParity  parity,
int  test_type 
)

Apply the Dslash operator (D_{eo} or D_{oe}) for Mobius DWF.

Parameters
h_outResult spinor field
h_inInput spinor field
paramContains all metadata regarding host and device storage
parityThe destination parity of the field
test_typeChoose a type of dslash operators

Definition at line 2015 of file interface_quda.cpp.

References quda::blas::ax(), quda::ColorSpinorParam::create, quda::ColorSpinorField::Create(), dirac, QudaInvertParam_s::dirac_order, quda::DiracMobius::Dslash4(), quda::DiracMobius::Dslash4pre(), quda::DiracMobius::Dslash5(), quda::DiracMobiusPC::Dslash5inv(), QudaInvertParam_s::dslash_type, errorQuda, gaugeLongPrecise, getVerbosity(), in, QudaInvertParam_s::input_location, quda::blas::norm2(), out, QudaInvertParam_s::output_location, popVerbosity(), printfQuda, printQudaInvertParam(), pushVerbosity(), QUDA_ASQTAD_DSLASH, QUDA_CPS_WILSON_DIRAC_ORDER, QUDA_DEBUG_VERBOSE, QUDA_EVEN_PARITY, QUDA_NULL_FIELD_CREATE, QUDA_ODD_PARITY, quda::setDiracParam(), and QudaInvertParam_s::verbosity.

Referenced by dslashCUDA().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ dumpMultigridQuda()

void dumpMultigridQuda ( void *  mg_instance,
QudaMultigridParam param 
)

Dump the null-space vectors to disk.

Parameters
[in]mg_instancePointer to the instance of multigrid_solver
[in]paramContains all metadata regarding host and device storage and solver parameters (QudaMultigridParam::vec_outfile sets the output filename prefix).

Definition at line 2726 of file interface_quda.cpp.

References checkGauge(), quda::MG::dumpNullVectors(), QudaMultigridParam_s::invert_param, quda::multigrid_solver::mg, popVerbosity(), profileInvert, profilerStart(), profilerStop(), pushVerbosity(), quda::QUDA_PROFILE_TOTAL, and QudaInvertParam_s::verbosity.

Referenced by main().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ eigensolveQuda()

void eigensolveQuda ( void **  h_evecs,
double_complex h_evals,
QudaEigParam param 
)

Perform the eigensolve. The problem matrix is defined by the invert param, the mode of solution is specified by the eig param. It is assumed that the gauge field has already been loaded via loadGaugeQuda().

Parameters
h_evecsArray of pointers to application eigenvectors
h_evalsHost side eigenvalues
paramContains all metadata regarding the type of solve.

Referenced by eigensolve_test(), and main().

Here is the caller graph for this function:

◆ endQuda()

void endQuda ( void  )

Finalize the library.

Definition at line 1461 of file interface_quda.cpp.

References quda::assertAllMemFree(), comm_finalize(), comms_initialized, quda::cublas::destroy(), quda::destroyDslashEvents(), quda::blas::end(), quda::pool::flush_device(), quda::pool::flush_pinned(), flushChronoQuda(), quda::flushForceMonitor(), freeCloverQuda(), freeGaugeQuda(), quda::LatticeField::freeGhostBuffer(), quda::cpuColorSpinorField::freeGhostBuffer(), getVerbosity(), host_free, initialized, momResident, quda::Nstream, num_failures_d, num_failures_h, quda::TimeProfile::Print(), quda::printAPIProfile(), printfQuda, quda::TimeProfile::PrintGlobal(), quda::printLaunchTimer(), quda::printPeakMemUsage(), profileAPE, profileClover, profileCloverForce, profileContract, profileCovDev, profileDslash, profileEigensolve, profileEnd, profileExtendedGauge, profileFatLink, profileGauge, profileGaugeForce, profileGaugeUpdate, profileHISQForce, profileInit, profileInit2End, profileInvert, profileMomAction, profileMulti, profilePhase, profilePlaq, profileProject, profileQCharge, profileStaggeredForce, profileSTOUT, QUDA_MAX_CHRONO, quda::QUDA_PROFILE_TOTAL, QUDA_SUMMARIZE, quda::saveProfile(), quda::saveTuneCache(), and streams.

Referenced by end(), end_quda_(), gauge_force_test(), hisq_force_end(), hisq_test(), llfat_test(), main(), plaq_test(), SU3test(), StaggeredDslashTest::TearDownTestCase(), DslashTest::TearDownTestCase(), and unitarize_link_test().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ flushChronoQuda()

void flushChronoQuda ( int  index)

Flush the chronological history for the given index.

Parameters
[in]indexIndex for which we are flushing

Definition at line 1448 of file interface_quda.cpp.

References chronoResident(), errorQuda, and QUDA_MAX_CHRONO.

Referenced by endQuda(), and flush_chrono_quda_().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ freeCloverQuda()

void freeCloverQuda ( void  )

Free QUDA's internal copy of the clover term and/or clover inverse.

Definition at line 1440 of file interface_quda.cpp.

References cloverPrecise, errorQuda, freeSloppyCloverQuda(), and initialized.

Referenced by endQuda(), free_clover_quda_(), and main().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ freeGaugeQuda()

void freeGaugeQuda ( void  )

Free QUDA's internal copy of the gauge field.

Definition at line 1259 of file interface_quda.cpp.

References errorQuda, extendedGaugeResident, freeSloppyGaugeQuda(), gaugeExtended, gaugeFatPrecise, gaugeLongExtended, gaugeLongPrecise, gaugePrecise, gaugeSmeared, and initialized.

Referenced by end(), endQuda(), free_gauge_quda_(), main(), plaq_test(), and SU3test().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ gaussGaugeQuda()

void gaussGaugeQuda ( unsigned long long  seed,
double  sigma 
)

Generate Gaussian distributed fields and store in the resident gauge field. We create a Gaussian-distributed su(n) field and exponentiate it, e.g., U = exp(sigma * H), where H is the distributed su(n) field and beta is the width of the distribution (beta = 0 results in a free field, and sigma = 1 has maximum disorder).

Parameters
seedThe seed used for the RNG
sigmaWidth of Gaussian distrubution

Definition at line 5386 of file interface_quda.cpp.

References errorQuda, quda::cudaGaugeField::exchangeExtendedGhost(), quda::gaugeGauss(), gaugePrecise, quda::LatticeFieldParam::ghostExchange, param, profileGauss, QUDA_GHOST_EXCHANGE_NO, quda::QUDA_PROFILE_COMPUTE, quda::QUDA_PROFILE_TOTAL, QUDA_RECONSTRUCT_12, R, quda::GaugeFieldParam::reconstruct, and redundant_comms.

Referenced by plaq_test().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ initCommsGridQuda()

void initCommsGridQuda ( int  nDim,
const int *  dims,
QudaCommsMap  func,
void *  fdata 
)

Declare the grid mapping ("logical topology" in QMP parlance) used for communications in a multi-GPU grid. This function should be called prior to initQuda(). The only case in which it's optional is when QMP is used for communication and the logical topology has already been declared by the application.

Parameters
nDimNumber of grid dimensions. "4" is the only supported value currently.
dimsArray of grid dimensions. dims[0]*dims[1]*dims[2]*dims[3] must equal the total number of MPI ranks or QMP nodes.
funcPointer to a user-supplied function that maps coordinates in the communication grid to MPI ranks (or QMP node IDs). If the pointer is NULL, the default mapping depends on whether QMP or MPI is being used for communication. With QMP, the existing logical topology is used if it's been declared. With MPI or as a fallback with QMP, the default ordering is lexicographical with the fourth ("t") index varying fastest.
fdataPointer to any data required by "func" (may be NULL)
See also
QudaCommsMap

Definition at line 401 of file interface_quda.cpp.

References comm_init(), comms_initialized, LexMapData::dims, errorQuda, lex_rank_from_coords(), LexMapData::ndim, and warningQuda.

Referenced by comm_set_gridsize_(), init_default_comms(), and initComms().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ initQuda()

void initQuda ( int  device)

Initialize the library. This function is actually a wrapper around calls to initQudaDevice() and initQudaMemory().

Parameters
deviceCUDA device number to use. In a multi-GPU build, this parameter may either be set explicitly on a per-process basis or set to -1 to enable a default allocation of devices to processes.

Definition at line 679 of file interface_quda.cpp.

References comms_initialized, init_default_comms(), initQudaDevice(), and initQudaMemory().

Referenced by eigensolve_test(), gauge_force_test(), hisq_force_init(), hisq_test(), init(), init_quda_(), invert_test(), llfat_test(), main(), plaq_test(), StaggeredDslashTest::SetUpTestCase(), DslashTest::SetUpTestCase(), SU3test(), and unitarize_link_test().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ initQudaDevice()

void initQudaDevice ( int  device)

Initialize the library. This is a low-level interface that is called by initQuda. Calling initQudaDevice requires that the user also call initQudaMemory before using QUDA.

Parameters
deviceCUDA device number to use. In a multi-GPU build, this parameter may either be set explicitly on a per-process basis or set to -1 to enable a default allocation of devices to processes.

Definition at line 483 of file interface_quda.cpp.

References checkCudaErrorNoSync, comm_gpuid(), comms_initialized, deviceProp, errorQuda, getVerbosity(), gitversion, initialized, length, printfQuda, profileInit, profileInit2End, QUDA_CPU_FIELD_LOCATION, QUDA_CUDA_FIELD_LOCATION, quda::QUDA_PROFILE_INIT, quda::QUDA_PROFILE_TOTAL, QUDA_SILENT, QUDA_SUMMARIZE, quda::quda_version, quda::reorder_location_set(), setNumaAffinityNVML(), and warningQuda.

Referenced by init_quda_device_(), and initQuda().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ initQudaMemory()

void initQudaMemory ( )

Initialize the library persistant memory allocations (both host and device). This is a low-level interface that is called by initQuda. Calling initQudaMemory requires that the user has previously called initQudaDevice.

Definition at line 638 of file interface_quda.cpp.

References checkCudaError, commDimPartitioned(), comms_initialized, quda::createDslashEvents(), quda::cublas::init(), quda::blas::init(), quda::pool::init(), init_default_comms(), quda::loadTuneCache(), mapped_malloc, quda::Nstream, num_failures_d, num_failures_h, profileInit, quda::QUDA_PROFILE_INIT, quda::QUDA_PROFILE_TOTAL, R, redundant_comms, and streams.

Referenced by init_quda_memory_(), and initQuda().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ invertMultiShiftQuda()

void invertMultiShiftQuda ( void **  _hp_x,
void *  _hp_b,
QudaInvertParam param 
)

Solve for multiple shifts (e.g., masses).

Parameters
_hp_xArray of solution spinor fields
_hp_bSource spinor fields
paramContains all metadata regarding host and device storage and solver parameters

Generic version of the multi-shift solver. Should work for most fermions. Note that offset[0] is not folded into the mass parameter.

For Wilson-type fermions, the solution_type must be MATDAG_MAT or MATPCDAG_MATPC, and solve_type must be NORMOP or NORMOP_PC. The solution and solve preconditioning have to match.

For Staggered-type fermions, the solution_type must be MATPC, and the solve type must be DIRECT_PC. This difference in convention is because preconditioned staggered operator is normal, unlike with Wilson-type fermions.

Definition at line 3579 of file interface_quda.cpp.

References QudaInvertParam_s::action, quda::blas::ax(), quda::blas::cDotProduct(), checkGauge(), QudaInvertParam_s::compute_action, QudaInvertParam_s::compute_true_res, quda::blas::copy(), quda::ColorSpinorParam::create, quda::ColorSpinorField::Create(), quda::createDirac(), QudaInvertParam_s::cuda_prec, QudaInvertParam_s::cuda_prec_refinement_sloppy, QudaInvertParam_s::cuda_prec_sloppy, quda::deflated_solver::d, quda::SolverParam::delta, dirac, QudaInvertParam_s::dslash_type, errorQuda, getVerbosity(), QudaInvertParam_s::gflops, initialized, QudaInvertParam_s::input_location, quda::SolverParam::iter, QudaInvertParam_s::iter, QudaInvertParam_s::iter_res_offset, quda::ColorSpinorParam::location, quda::deflated_solver::m, QudaInvertParam_s::make_resident_solution, QudaInvertParam_s::mass, quda::massRescale(), quda::blas::norm2(), QudaInvertParam_s::num_offset, QudaInvertParam_s::offset, QudaInvertParam_s::output_location, param, popVerbosity(), quda::pow(), quda::LatticeFieldParam::Precision(), printfQuda, profileMulti, profilerStart(), profilerStop(), pushVerbosity(), QUDA_ASQTAD_DSLASH, QUDA_COPY_FIELD_CREATE, QUDA_DIRECT_PC_SOLVE, QUDA_DIRECT_SOLVE, QUDA_HEAVY_QUARK_RESIDUAL, QUDA_MAT_SOLUTION, QUDA_MATPC_SOLUTION, QUDA_MATPCDAG_MATPC_SOLUTION, QUDA_MAX_MULTI_SHIFT, QUDA_NORMOP_PC_SOLVE, QUDA_NULL_FIELD_CREATE, quda::QUDA_PROFILE_D2H, quda::QUDA_PROFILE_EPILOGUE, quda::QUDA_PROFILE_FREE, quda::QUDA_PROFILE_H2D, quda::QUDA_PROFILE_INIT, quda::QUDA_PROFILE_PREAMBLE, quda::QUDA_PROFILE_TOTAL, QUDA_SOURCE_NORMALIZATION, QUDA_STAGGERED_DSLASH, QUDA_SUMMARIZE, QUDA_USE_INIT_GUESS_YES, QUDA_VERBOSE, QUDA_ZERO_FIELD_CREATE, QudaInvertParam_s::reliable_delta_refinement, QudaInvertParam_s::residual_type, QudaInvertParam_s::residue, quda::saveTuneCache(), QudaInvertParam_s::secs, quda::Dirac::setMass(), quda::DiracMatrix::shift, QudaInvertParam_s::solution_type, QudaInvertParam_s::solve_type, QudaInvertParam_s::solver_normalization, quda::sqrt(), tmp, quda::SolverParam::tol, tol_hq, quda::SolverParam::tol_hq, QudaInvertParam_s::tol_hq_offset, QudaInvertParam_s::tol_offset, quda::SolverParam::true_res, quda::SolverParam::true_res_hq, quda::SolverParam::true_res_hq_offset, QudaInvertParam_s::true_res_hq_offset, quda::SolverParam::true_res_offset, QudaInvertParam_s::true_res_offset, quda::unscaled_shifts, quda::SolverParam::updateInvertParam(), quda::SolverParam::use_init_guess, quda::ColorSpinorParam::v, QudaInvertParam_s::verbosity, X, and quda::LatticeField::X().

Referenced by invert_multishift_quda_(), invert_test(), and main().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ invertMultiSrcQuda()

void invertMultiSrcQuda ( void **  _hp_x,
void **  _hp_b,
QudaInvertParam param 
)

Perform the solve like but for multiples right hand sides.

Parameters
_hp_xArray of solution spinor fields
_hp_bArray of source spinor fields
paramContains all metadata regarding
paramContains all metadata regarding host and device storage and solver parameters

Generic version of the multi-shift solver. Should work for most fermions. Note that offset[0] is not folded into the mass parameter.

At present, the solution_type must be MATDAG_MAT or MATPCDAG_MATPC, and solve_type must be NORMOP or NORMOP_PC. The solution and solve preconditioning have to match.

Definition at line 3234 of file interface_quda.cpp.

References quda::blas::ax(), quda::Solver::blocksolve(), checkGauge(), quda::ColorSpinorField::Component(), quda::ColorSpinorParam::composite_dim, quda::blas::copy(), quda::ColorSpinorParam::create, quda::Solver::create(), quda::ColorSpinorField::Create(), quda::createDirac(), cudaGauge, quda::deflated_solver::d, dirac, errorQuda, getVerbosity(), QudaInvertParam_s::gflops, in, initialized, QudaInvertParam_s::input_location, QudaInvertParam_s::inv_type_precondition, quda::ColorSpinorParam::is_component, quda::ColorSpinorParam::is_composite, QudaInvertParam_s::iter, quda::ColorSpinorParam::location, quda::deflated_solver::m, QudaInvertParam_s::make_resident_solution, quda::massRescale(), quda::Dirac::Mdag(), quda::blas::norm2(), QudaInvertParam_s::num_src, out, QudaInvertParam_s::output_location, popVerbosity(), quda::Dirac::prepare(), printfQuda, printQudaInvertParam(), profileInvert, pushVerbosity(), QUDA_DEBUG_VERBOSE, QUDA_DIRECT_PC_SOLVE, QUDA_DIRECT_SOLVE, QUDA_MAT_SOLUTION, QUDA_MATDAG_MAT_SOLUTION, QUDA_MATPC_SOLUTION, QUDA_MATPCDAG_MATPC_SOLUTION, QUDA_MG_INVERTER, QUDA_NORMERR_PC_SOLVE, QUDA_NORMERR_SOLVE, QUDA_NORMOP_PC_SOLVE, QUDA_NULL_FIELD_CREATE, quda::QUDA_PROFILE_D2H, quda::QUDA_PROFILE_EPILOGUE, quda::QUDA_PROFILE_H2D, quda::QUDA_PROFILE_TOTAL, QUDA_SOURCE_NORMALIZATION, QUDA_USE_INIT_GUESS_YES, QUDA_VERBOSE, QUDA_ZERO_FIELD_CREATE, quda::Dirac::reconstruct(), quda::saveTuneCache(), QudaInvertParam_s::secs, QudaInvertParam_s::solution_type, QudaInvertParam_s::solve_type, QudaInvertParam_s::solver_normalization, quda::sqrt(), tmp, quda::SolverParam::updateInvertParam(), QudaInvertParam_s::use_init_guess, quda::ColorSpinorParam::v, QudaInvertParam_s::verbosity, X, and quda::LatticeField::X().

Referenced by invert_test(), and main().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ invertQuda()

void invertQuda ( void *  h_x,
void *  h_b,
QudaInvertParam param 
)

Perform the solve, according to the parameters set in param. It is assumed that the gauge field has already been loaded via loadGaugeQuda().

Parameters
h_xSolution spinor field
h_bSource spinor field
paramContains all metadata regarding host and device storage and solver parameters

Definition at line 2830 of file interface_quda.cpp.

References QudaInvertParam_s::action, quda::blas::ax(), quda::blas::cDotProduct(), checkGauge(), QudaInvertParam_s::chrono_index, QudaInvertParam_s::chrono_make_resident, QudaInvertParam_s::chrono_max_dim, QudaInvertParam_s::chrono_precision, QudaInvertParam_s::chrono_replace_last, QudaInvertParam_s::chrono_use_resident, chronoResident(), QudaInvertParam_s::compute_action, quda::blas::copy(), quda::ColorSpinorParam::create, quda::Solver::create(), quda::ColorSpinorField::Create(), quda::createDirac(), QudaInvertParam_s::cuda_prec, QudaInvertParam_s::cuda_prec_sloppy, cudaGauge, quda::deflated_solver::d, dirac, errorQuda, getVerbosity(), QudaInvertParam_s::gflops, in, initialized, QudaInvertParam_s::input_location, QudaInvertParam_s::inv_type_precondition, QudaInvertParam_s::iter, quda::ColorSpinorParam::location, quda::deflated_solver::m, QudaInvertParam_s::make_resident_solution, quda::massRescale(), quda::Dirac::Mdag(), quda::blas::norm2(), out, QudaInvertParam_s::output_location, popVerbosity(), quda::LatticeField::Precision(), quda::Dirac::prepare(), printfQuda, printQudaInvertParam(), profileInvert, profilerStart(), profilerStop(), pushVerbosity(), QUDA_COPY_FIELD_CREATE, QUDA_DEBUG_VERBOSE, QUDA_DIRECT_PC_SOLVE, QUDA_DIRECT_SOLVE, QUDA_MAT_SOLUTION, QUDA_MATDAG_MAT_SOLUTION, QUDA_MATPC_SOLUTION, QUDA_MATPCDAG_MATPC_SOLUTION, QUDA_MAX_CHRONO, QUDA_MG_INVERTER, QUDA_NORMERR_PC_SOLVE, QUDA_NORMERR_SOLVE, QUDA_NORMOP_PC_SOLVE, QUDA_NULL_FIELD_CREATE, quda::QUDA_PROFILE_CHRONO, quda::QUDA_PROFILE_D2H, quda::QUDA_PROFILE_EPILOGUE, quda::QUDA_PROFILE_FREE, quda::QUDA_PROFILE_H2D, quda::QUDA_PROFILE_PREAMBLE, quda::QUDA_PROFILE_TOTAL, QUDA_SOURCE_NORMALIZATION, QUDA_USE_INIT_GUESS_YES, QUDA_VERBOSE, quda::Dirac::reconstruct(), quda::saveTuneCache(), QudaInvertParam_s::secs, quda::ColorSpinorParam::setPrecision(), quda::ColorSpinorField::SiteSubset(), QudaInvertParam_s::solution_type, QudaInvertParam_s::solve_type, QudaInvertParam_s::solver_normalization, quda::sqrt(), tmp, tmp2, quda::SolverParam::updateInvertParam(), QudaInvertParam_s::use_init_guess, QudaInvertParam_s::use_resident_solution, quda::ColorSpinorParam::v, QudaInvertParam_s::verbosity, X, quda::LatticeField::X(), and quda::blas::zero().

Referenced by invert_quda_(), invert_test(), and main().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ lanczosQuda()

void lanczosQuda ( int  k0,
int  m,
void *  hp_Apsi,
void *  hp_r,
void *  hp_V,
void *  hp_alpha,
void *  hp_beta,
QudaEigParam eig_param 
)

Perform the solve, according to the parameters set in param. It is assumed that the gauge field has already been loaded via loadGaugeQuda().

Parameters
h_xSolution spinor field
h_bSource spinor field
paramContains all metadata regarding host and device storage and solver parameters

◆ loadCloverQuda()

void loadCloverQuda ( void *  h_clover,
void *  h_clovinv,
QudaInvertParam inv_param 
)

Load the clover term and/or the clover inverse from the host. Either h_clover or h_clovinv may be set to NULL.

Parameters
h_cloverBase pointer to host clover field
h_cloverinvBase pointer to host clover inverse field
inv_paramContains all metadata regarding host and device storage

Definition at line 985 of file interface_quda.cpp.

References quda::GaugeField::Anisotropy(), quda::CloverField::Bytes(), checkCudaError, QudaInvertParam_s::cl_pad, quda::CloverFieldParam::clover, QudaInvertParam_s::clover_coeff, QudaInvertParam_s::clover_cpu_prec, QudaInvertParam_s::clover_cuda_prec, QudaInvertParam_s::clover_cuda_prec_precondition, QudaInvertParam_s::clover_cuda_prec_refinement_sloppy, QudaInvertParam_s::clover_cuda_prec_sloppy, QudaInvertParam_s::clover_location, QudaInvertParam_s::clover_order, QudaInvertParam_s::clover_rho, quda::CloverFieldParam::cloverInv, quda::cloverInvert(), cloverPrecise, QudaInvertParam_s::compute_clover, QudaInvertParam_s::compute_clover_inverse, QudaInvertParam_s::compute_clover_trlog, quda::cudaCloverField::copy(), quda::CloverFieldParam::create, createCloverQuda(), quda::CloverFieldParam::csw, quda::CloverField::Csw(), quda::CloverFieldParam::direct, QudaInvertParam_s::dslash_type, errorQuda, freeSloppyCloverQuda(), getVerbosity(), in, initialized, invalidate_clover, quda::CloverFieldParam::inverse, quda::inverse(), quda::CloverFieldParam::invNorm, QudaInvertParam_s::kappa, loadSloppyCloverQuda(), QudaInvertParam_s::matpc_type, QudaInvertParam_s::mu, quda::CloverFieldParam::mu2, quda::LatticeFieldParam::nDim, quda::CloverFieldParam::norm, quda::CloverFieldParam::order, quda::LatticeFieldParam::pad, popVerbosity(), prec, printfQuda, printQudaInvertParam(), profileClover, pushVerbosity(), QUDA_CLOVER_WILSON_DSLASH, QUDA_CPU_FIELD_LOCATION, QUDA_DEBUG_VERBOSE, QUDA_DIRECT_PC_SOLVE, QUDA_FULL_SITE_SUBSET, QUDA_HALF_PRECISION, QUDA_MATPC_EVEN_EVEN_ASYMMETRIC, QUDA_MATPC_ODD_ODD_ASYMMETRIC, QUDA_MATPC_SOLUTION, QUDA_MATPCDAG_MATPC_SOLUTION, QUDA_NORMERR_PC_SOLVE, QUDA_NORMOP_PC_SOLVE, QUDA_NULL_FIELD_CREATE, quda::QUDA_PROFILE_COMPUTE, quda::QUDA_PROFILE_D2H, quda::QUDA_PROFILE_EPILOGUE, quda::QUDA_PROFILE_FREE, quda::QUDA_PROFILE_H2D, quda::QUDA_PROFILE_INIT, quda::QUDA_PROFILE_TOTAL, QUDA_REFERENCE_FIELD_CREATE, QUDA_TWISTED_CLOVER_DSLASH, QUDA_VERBOSE, qudaMemcpy, QudaInvertParam_s::return_clover, QudaInvertParam_s::return_clover_inverse, quda::CloverFieldParam::setPrecision(), quda::CloverField::setRho(), quda::LatticeFieldParam::siteSubset, QudaInvertParam_s::solution_type, QudaInvertParam_s::solve_type, quda::CloverField::TrLog(), QudaInvertParam_s::trlogA, quda::CloverFieldParam::twisted, quda::CloverField::V(), QudaInvertParam_s::verbosity, warningQuda, quda::LatticeFieldParam::x, and quda::LatticeField::X().

Referenced by init(), load_clover_quda_(), and main().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ loadGaugeQuda()

void loadGaugeQuda ( void *  h_gauge,
QudaGaugeParam param 
)

Load the gauge field from the host.

Parameters
h_gaugeBase pointer to host gauge field (regardless of dimensionality)
paramContains all metadata regarding host and device storage

Definition at line 729 of file interface_quda.cpp.

References quda::GaugeField::checksum(), commDimPartitioned(), quda::GaugeFieldParam::compute_fat_link_max, quda::cudaGaugeField::copy(), quda::GaugeFieldParam::create, createExtendedGauge(), QudaGaugeParam_s::cuda_prec, QudaGaugeParam_s::cuda_prec_precondition, QudaGaugeParam_s::cuda_prec_refinement_sloppy, QudaGaugeParam_s::cuda_prec_sloppy, errorQuda, quda::cudaGaugeField::exchangeGhost(), extendedGaugeResident, QudaGaugeParam_s::ga_pad, gauge_param, gaugeFatPrecise, gaugeFatPrecondition, gaugeFatRefinement, gaugeFatSloppy, gaugeLongPrecise, gaugeLongPrecondition, gaugeLongRefinement, gaugeLongSloppy, gaugePrecise, gaugePrecondition, gaugeRefinement, gaugeSloppy, gaugeSmeared, getVerbosity(), quda::LatticeFieldParam::ghostExchange, in, initialized, invalidate_clover, QudaGaugeParam_s::location, quda::GaugeFieldParam::order, quda::GaugeField::Order(), QudaGaugeParam_s::overlap, quda::LatticeFieldParam::pad, printfQuda, printQudaGaugeParam(), profileGauge, QUDA_ASQTAD_FAT_LINKS, QUDA_ASQTAD_LONG_LINKS, QUDA_BQCD_GAUGE_ORDER, QUDA_CPU_FIELD_LOCATION, QUDA_DEBUG_VERBOSE, QUDA_GHOST_EXCHANGE_NO, QUDA_GHOST_EXCHANGE_PAD, QUDA_NULL_FIELD_CREATE, quda::QUDA_PROFILE_COMPUTE, quda::QUDA_PROFILE_FREE, quda::QUDA_PROFILE_H2D, quda::QUDA_PROFILE_INIT, quda::QUDA_PROFILE_TOTAL, QUDA_SMEARED_LINKS, QUDA_VERBOSE, QUDA_WILSON_LINKS, R, quda::LatticeField::R(), quda::GaugeFieldParam::reconstruct, QudaGaugeParam_s::reconstruct, quda::GaugeField::Reconstruct(), QudaGaugeParam_s::reconstruct_precondition, QudaGaugeParam_s::reconstruct_refinement_sloppy, QudaGaugeParam_s::reconstruct_sloppy, quda::GaugeFieldParam::setPrecision(), QudaGaugeParam_s::type, and QudaGaugeParam_s::use_resident_gauge.

Referenced by computeStaggeredPlaquetteQDPOrder(), eigensolve_test(), init(), invert_test(), load_gauge_quda_(), main(), plaq_test(), and SU3test().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MatDagMatQuda()

void MatDagMatQuda ( void *  h_out,
void *  h_in,
QudaInvertParam inv_param 
)

◆ MatQuda()

void MatQuda ( void *  h_out,
void *  h_in,
QudaInvertParam inv_param 
)

◆ momActionQuda()

double momActionQuda ( void *  momentum,
QudaGaugeParam param 
)

Evaluate the momentum contribution to the Hybrid Monte Carlo action.

Parameters
momentumThe momentum field
paramThe parameters of the external fields and the computation settings
Returns
momentum action

Definition at line 5084 of file interface_quda.cpp.

References checkCudaError, quda::computeMomAction(), cpuMom, quda::GaugeFieldParam::create, cudaMom, errorQuda, gParam, quda::cudaGaugeField::loadCPUField(), QudaGaugeParam_s::make_resident_mom, momResident, quda::GaugeFieldParam::order, profileMomAction, QUDA_ASQTAD_MOM_LINKS, QUDA_FLOAT2_GAUGE_ORDER, QUDA_NULL_FIELD_CREATE, quda::QUDA_PROFILE_COMPUTE, quda::QUDA_PROFILE_FREE, quda::QUDA_PROFILE_H2D, quda::QUDA_PROFILE_INIT, quda::QUDA_PROFILE_TOTAL, QUDA_RECONSTRUCT_10, QUDA_RECONSTRUCT_NO, QUDA_TIFR_GAUGE_ORDER, QUDA_TIFR_PADDED_GAUGE_ORDER, quda::GaugeFieldParam::reconstruct, and QudaGaugeParam_s::use_resident_mom.

Referenced by kinetic_quda_().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ newDeflationQuda()

void* newDeflationQuda ( QudaEigParam param)

Create deflation solver resources.

Definition at line 2809 of file interface_quda.cpp.

References quda::deflated_solver::defl, quda::deflated_solver::deflated_solver(), quda::flushProfile(), openMagma(), profileInvert, quda::QUDA_PROFILE_TOTAL, and quda::saveProfile().

Referenced by main().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ newMultigridQuda()

void* newMultigridQuda ( QudaMultigridParam param)

Setup the multigrid solver, according to the parameters set in param. It is assumed that the gauge field has already been loaded via loadGaugeQuda().

Parameters
paramContains all metadata regarding host and device storage and solver parameters

Definition at line 2624 of file interface_quda.cpp.

References QudaMultigridParam_s::invert_param, quda::multigrid_solver::mg, quda::multigrid_solver::multigrid_solver(), popVerbosity(), profileInvert, profilerStart(), profilerStop(), pushVerbosity(), quda::QUDA_PROFILE_TOTAL, quda::saveTuneCache(), and QudaInvertParam_s::verbosity.

Referenced by main().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ newQudaEigParam()

QudaEigParam newQudaEigParam ( void  )

A new QudaEigParam should always be initialized immediately after it's defined (and prior to explicitly setting its members) using this function. Typical usage is as follows:

QudaEigParam eig_param = newQudaEigParam();

Referenced by eigensolve_test(), main(), and printQudaGaugeParam().

Here is the caller graph for this function:

◆ newQudaGaugeParam()

QudaGaugeParam newQudaGaugeParam ( void  )

A new QudaGaugeParam should always be initialized immediately after it's defined (and prior to explicitly setting its members) using this function. Typical usage is as follows:

QudaGaugeParam gauge_param = newQudaGaugeParam();

Referenced by computeStaggeredPlaquetteQDPOrder(), eigensolve_test(), gauge_force_test(), hisq_test(), init(), invert_test(), llfat_test(), main(), new_quda_gauge_param_(), plaq_test(), GaugeAlgTest::SetUp(), SU3test(), and unitarize_link_test().

Here is the caller graph for this function:

◆ newQudaInvertParam()

QudaInvertParam newQudaInvertParam ( void  )

A new QudaInvertParam should always be initialized immediately after it's defined (and prior to explicitly setting its members) using this function. Typical usage is as follows:

QudaInvertParam invert_param = newQudaInvertParam();

Referenced by eigensolve_test(), init(), invert_test(), main(), new_quda_invert_param_(), printQudaCloverParam(), and test().

Here is the caller graph for this function:

◆ newQudaMultigridParam()

QudaMultigridParam newQudaMultigridParam ( void  )

A new QudaMultigridParam should always be initialized immediately after it's defined (and prior to explicitly setting its members) using this function. Typical usage is as follows:

QudaMultigridParam mg_param = newQudaMultigridParam();

Referenced by main(), and printQudaInvertParam().

Here is the caller graph for this function:

◆ openMagma()

void openMagma ( )

Open/Close MAGMA library

Definition at line 95 of file interface_quda.cpp.

References InitMagma, OpenMagma(), and printfQuda.

Referenced by newDeflationQuda().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ pack_ghost()

void pack_ghost ( void **  cpuLink,
void **  cpuGhost,
int  nFace,
QudaPrecision  precision 
)

◆ performAPEnStep()

void performAPEnStep ( unsigned int  nSteps,
double  alpha 
)

Performs APE smearing on gaugePrecise and stores it in gaugeSmeared

Parameters
nStepsNumber of steps to apply.
alphaAlpha coefficient for APE smearing.

Definition at line 5528 of file interface_quda.cpp.

References quda::APEStep(), createExtendedGauge(), errorQuda, quda::cudaGaugeField::exchangeExtendedGhost(), gaugeSmeared, getVerbosity(), gParam, quda::plaquette(), printfQuda, profileAPE, quda::QUDA_PROFILE_TOTAL, QUDA_SUMMARIZE, R, and redundant_comms.

Referenced by SU3test().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ performOvrImpSTOUTnStep()

void performOvrImpSTOUTnStep ( unsigned int  nSteps,
double  rho,
double  epsilon 
)

Performs Over Imroved STOUT smearing on gaugePrecise and stores it in gaugeSmeared

Parameters
nStepsNumber of steps to apply.
rhoRho coefficient for STOUT smearing.
epsilonEpsilon coefficient for Over Improved STOUT smearing.

Definition at line 5598 of file interface_quda.cpp.

References createExtendedGauge(), errorQuda, quda::cudaGaugeField::exchangeExtendedGhost(), gaugeSmeared, getVerbosity(), gParam, quda::OvrImpSTOUTStep(), quda::plaquette(), printfQuda, profileOvrImpSTOUT, profileSTOUT, quda::QUDA_PROFILE_TOTAL, QUDA_SUMMARIZE, R, and redundant_comms.

Referenced by SU3test().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ performSTOUTnStep()

void performSTOUTnStep ( unsigned int  nSteps,
double  rho 
)

Performs STOUT smearing on gaugePrecise and stores it in gaugeSmeared

Parameters
nStepsNumber of steps to apply.
rhoRho coefficient for STOUT smearing.

Definition at line 5563 of file interface_quda.cpp.

References createExtendedGauge(), errorQuda, quda::cudaGaugeField::exchangeExtendedGhost(), gaugeSmeared, getVerbosity(), gParam, quda::plaquette(), printfQuda, profileSTOUT, quda::QUDA_PROFILE_TOTAL, QUDA_SUMMARIZE, R, redundant_comms, and quda::STOUTStep().

Referenced by SU3test().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ performWuppertalnStep()

void performWuppertalnStep ( void *  h_out,
void *  h_in,
QudaInvertParam param,
unsigned int  nSteps,
double  alpha 
)

Performs Wuppertal smearing on a given spinor using the gauge field gaugeSmeared, if it exist, or gaugePrecise if no smeared field is present.

Parameters
h_outResult spinor field
h_inInput spinor field
paramContains all metadata regarding host and device storage and operator which will be applied to the spinor
nStepsNumber of steps to apply.
alphaAlpha coefficient for Wuppertal smearing.

Definition at line 5457 of file interface_quda.cpp.

References quda::copyExtendedGauge(), quda::GaugeFieldParam::create, quda::ColorSpinorParam::create, quda::ColorSpinorField::Create(), errorQuda, quda::cudaGaugeField::exchangeGhost(), gaugePrecise, getVerbosity(), gParam, in, QudaInvertParam_s::input_location, quda::ColorSpinorParam::location, quda::norm(), quda::blas::norm2(), out, QudaInvertParam_s::output_location, parity, popVerbosity(), printfQuda, printQudaInvertParam(), profileWuppertal, pushVerbosity(), QUDA_CUDA_FIELD_LOCATION, QUDA_DEBUG_VERBOSE, QUDA_NULL_FIELD_CREATE, quda::QUDA_PROFILE_TOTAL, QUDA_VERBOSE, quda::ColorSpinorParam::v, QudaInvertParam_s::verbosity, quda::wuppertalStep(), and quda::LatticeField::X().

Here is the call graph for this function:

◆ plaqQuda()

void plaqQuda ( double  plaq[3])

Computes the total, spatial and temporal plaquette averages of the loaded gauge configuration.

Parameters
Arrayfor storing the averages (total, spatial, temporal)

Definition at line 5419 of file interface_quda.cpp.

References createExtendedGauge(), errorQuda, extendedGaugeResident, quda::plaquette(), profilePlaq, quda::QUDA_PROFILE_COMPUTE, quda::QUDA_PROFILE_TOTAL, and R.

Referenced by computeStaggeredPlaquetteQDPOrder(), main(), plaq_quda_(), plaq_test(), and SU3test().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ printQudaEigParam()

void printQudaEigParam ( QudaEigParam param)

Print the members of QudaEigParam.

Parameters
paramThe QudaEigParam whose elements we are to print.

Definition at line 143 of file check_params.h.

References eig_type, INVALID_DOUBLE, INVALID_INT, mem_type_ritz, param, printfQuda, QUDA_BOOLEAN_FALSE, QUDA_BOOLEAN_INVALID, QUDA_BOOLEAN_TRUE, QUDA_CUDA_FIELD_LOCATION, QUDA_EIG_INVALID, QUDA_EIG_TR_LANCZOS, QUDA_EIGEN_EXTLIB, QUDA_EXTLIB_INVALID, QUDA_INVALID_FIELD_LOCATION, QUDA_MEMORY_DEVICE, QUDA_MEMORY_INVALID, QUDA_SPECTRUM_LR_EIG, and tol.

Referenced by eigensolveQuda().

Here is the caller graph for this function:

◆ printQudaGaugeParam()

void printQudaGaugeParam ( QudaGaugeParam param)

◆ printQudaInvertParam()

void printQudaInvertParam ( QudaInvertParam param)

Print the members of QudaInvertParam.

Parameters
paramThe QudaInvertParam whose elements we are to print.

Whether to use a pipelined solver

< Number of offsets in the multi-shift solver

< Number of offsets in the multi-shift solver

< width of domain overlaps

Definition at line 277 of file check_params.h.

References ca_basis, QudaInvertParam_s::ca_basis, ca_lambda_max, ca_lambda_min, QudaInvertParam_s::chrono_precision, QudaInvertParam_s::compute_action, cpu_prec, cuda_prec, QudaInvertParam_s::cuda_prec, cuda_prec_precondition, QudaInvertParam_s::cuda_prec_precondition, cuda_prec_refinement_sloppy, QudaInvertParam_s::cuda_prec_refinement_sloppy, cuda_prec_ritz, cuda_prec_sloppy, QudaInvertParam_s::cuda_prec_sloppy, dagger, deflation_grid, dslash_type, QudaInvertParam_s::dslash_type, eigcg_max_restarts, eigenval_tol, gcrNkrylov, quda::get_pointer_location(), inc_tol, QudaInvertParam_s::input_location, inv_type, QudaInvertParam_s::inv_type, QudaInvertParam_s::inv_type_precondition, INVALID_DOUBLE, INVALID_INT, kappa, laplace3D, Ls, mass, matpc_type, max_restart_num, max_search_dim, mu, nev, newQudaMultigridParam(), QudaInvertParam_s::num_offset, omega, QudaInvertParam_s::output_location, pipeline, printfQuda, printQudaCloverParam(), QUDA_ADDITIVE_SCHWARZ, QUDA_ASQTAD_DSLASH, QUDA_BICGSTAB_INVERTER, QUDA_CA_CG_INVERTER, QUDA_CA_CGNE_INVERTER, QUDA_CA_CGNR_INVERTER, QUDA_CA_GCR_INVERTER, QUDA_CG_INVERTER, QUDA_CHEBYSHEV_BASIS, QUDA_CPU_FIELD_LOCATION, QUDA_DAG_INVALID, QUDA_DEFAULT_NORMALIZATION, QUDA_DOMAIN_WALL_4D_DSLASH, QUDA_DOMAIN_WALL_DSLASH, QUDA_EIGEN_EXTLIB, QUDA_EXTLIB_INVALID, QUDA_GCR_INVERTER, QUDA_HEAVY_QUARK_RESIDUAL, QUDA_INVALID_BASIS, QUDA_INVALID_DIRAC_ORDER, QUDA_INVALID_DSLASH, QUDA_INVALID_FIELD_LOCATION, QUDA_INVALID_GAMMA_BASIS, QUDA_INVALID_INVERTER, QUDA_INVALID_NORMALIZATION, QUDA_INVALID_PRECISION, QUDA_INVALID_RESIDUAL, QUDA_INVALID_SCHWARZ, QUDA_INVALID_SOLUTION, QUDA_INVALID_SOLVE, QUDA_INVALID_VERBOSITY, QUDA_L2_RELATIVE_RESIDUAL, QUDA_MATPC_INVALID, QUDA_MOBIUS_DWF_DSLASH, QUDA_MPBICGSTAB_INVERTER, QUDA_MPCG_INVERTER, QUDA_MR_INVERTER, QUDA_POWER_BASIS, QUDA_PRESERVE_SOURCE_INVALID, QUDA_SINGLE_PRECISION, QUDA_STAGGERED_DSLASH, QUDA_TWIST_INVALID, QUDA_TWISTED_MASS_DSLASH, QUDA_USE_INIT_GUESS_INVALID, QUDA_USE_INIT_GUESS_NO, reliable_delta, QudaInvertParam_s::reliable_delta, QudaInvertParam_s::reliable_delta_refinement, QudaInvertParam_s::residual_type, schwarz_type, solution_accumulator_pipeline, solution_type, solve_type, tol, tol_hq, tol_restart, twist_flavor, verbosity, and warningQuda.

Referenced by cloverQuda(), dslashQuda(), dslashQuda_4dpc(), dslashQuda_mdwf(), eigensolveQuda(), invertMultiSrcQuda(), invertQuda(), loadCloverQuda(), MatDagMatQuda(), MatQuda(), performWuppertalnStep(), and printQudaMultigridParam().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ printQudaMultigridParam()

void printQudaMultigridParam ( QudaMultigridParam param)

Print the members of QudaMultigridParam.

Parameters
paramThe QudaMultigridParam whose elements we are to print.

Definition at line 598 of file check_params.h.

References coarse_solver, coarse_solver_ca_basis, coarse_solver_ca_basis_size, coarse_solver_ca_lambda_max, coarse_solver_ca_lambda_min, coarse_solver_maxiter, coarse_solver_tol, QudaMultigridParam_s::compute_null_vector, errorQuda, generate_all_levels, QudaMultigridParam_s::generate_all_levels, geo_block_size, INVALID_DOUBLE, INVALID_INT, QudaMultigridParam_s::invert_param, mu_factor, n_block_ortho, QudaMultigridParam_s::n_level, QudaMultigridParam_s::n_vec, nu_post, nu_pre, num_setup_iter, omega, post_orthonormalize, pre_orthonormalize, printfQuda, printQudaInvertParam(), QUDA_BICGSTAB_INVERTER, QUDA_BOOLEAN_FALSE, QUDA_BOOLEAN_INVALID, QUDA_BOOLEAN_TRUE, QUDA_COMPUTE_NULL_VECTOR_INVALID, QUDA_COMPUTE_NULL_VECTOR_YES, QUDA_CUDA_FIELD_LOCATION, QUDA_INVALID_BASIS, QUDA_INVALID_FIELD_LOCATION, QUDA_INVALID_INVERTER, QUDA_INVALID_PRECISION, QUDA_INVALID_SCHWARZ, QUDA_INVALID_SETUP_TYPE, QUDA_INVALID_SOLUTION, QUDA_INVALID_SOLVE, QUDA_INVALID_VERBOSITY, QUDA_MAX_MG_LEVEL, QUDA_MG_CYCLE_INVALID, QUDA_NULL_VECTOR_SETUP, QUDA_POWER_BASIS, QUDA_SILENT, QUDA_SINGLE_PRECISION, setup_ca_basis, setup_ca_basis_size, setup_ca_lambda_max, setup_ca_lambda_min, setup_location, setup_maxiter, setup_maxiter_refresh, setup_tol, setup_type, smoother_solve_type, smoother_tol, and verbosity.

Referenced by quda::multigrid_solver::multigrid_solver().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ projectSU3Quda()

void projectSU3Quda ( void *  gauge_h,
double  tol,
QudaGaugeParam param 
)

◆ qChargeDensityQuda()

double qChargeDensityQuda ( void *  qDensity)

Calculates the topological charge from gaugeSmeared, if it exist, or from gaugePrecise if no smeared fields are present.

Parameters
[out]qDensityarray holding Q charge density

Definition at line 5880 of file interface_quda.cpp.

References quda::computeFmunu(), quda::computeQChargeDensity(), createExtendedGauge(), device_free, device_malloc, extendedGaugeResident, gaugeSmeared, quda::LatticeField::Precision(), profileQCharge, QUDA_FLOAT2_GAUGE_ORDER, QUDA_FULL_SITE_SUBSET, QUDA_GHOST_EXCHANGE_NO, quda::QUDA_PROFILE_COMPUTE, quda::QUDA_PROFILE_D2H, quda::QUDA_PROFILE_FREE, quda::QUDA_PROFILE_INIT, quda::QUDA_PROFILE_TOTAL, QUDA_RECONSTRUCT_NO, QUDA_TENSOR_GEOMETRY, qudaMemcpy, R, quda::LatticeFieldParam::siteSubset, quda::size, quda::LatticeField::Volume(), and quda::LatticeField::X().

Referenced by SU3test().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ qChargeQuda()

double qChargeQuda ( )

Calculates the topological charge from gaugeSmeared, if it exist, or from gaugePrecise if no smeared fields are present.

Definition at line 5846 of file interface_quda.cpp.

References quda::computeFmunu(), quda::computeQCharge(), createExtendedGauge(), extendedGaugeResident, gaugeSmeared, quda::LatticeField::Precision(), profileQCharge, QUDA_FLOAT2_GAUGE_ORDER, QUDA_FULL_SITE_SUBSET, QUDA_GHOST_EXCHANGE_NO, quda::QUDA_PROFILE_COMPUTE, quda::QUDA_PROFILE_INIT, quda::QUDA_PROFILE_TOTAL, QUDA_RECONSTRUCT_NO, QUDA_TENSOR_GEOMETRY, R, quda::LatticeFieldParam::siteSubset, and quda::LatticeField::X().

Referenced by main(), and SU3test().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ qudaSetCommHandle()

void qudaSetCommHandle ( void *  mycomm)
Parameters
mycommUser provided MPI communicator in place of MPI_COMM_WORLD

◆ saveGaugeFieldQuda()

void saveGaugeFieldQuda ( void *  outGauge,
void *  inGauge,
QudaGaugeParam param 
)

Copy the QUDA gauge (matrix) field on the device to the CPU

Parameters
outGaugePointer to the host gauge field
inGaugePointer to the device gauge field (QUDA device field)
paramThe parameters of the host and device fields

Definition at line 4252 of file interface_quda.cpp.

References cpuGauge, cudaGauge, quda::GaugeFieldParam::geometry, quda::GaugeField::Geometry(), gParam, QUDA_GENERAL_LINKS, and quda::cudaGaugeField::saveCPUField().

Referenced by main().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ saveGaugeQuda()

void saveGaugeQuda ( void *  h_gauge,
QudaGaugeParam param 
)

◆ set_dim()

void set_dim ( int *  )

◆ setMPICommHandleQuda()

void setMPICommHandleQuda ( void *  mycomm)

Definition at line 368 of file interface_quda.cpp.

◆ setVerbosityQuda()

void setVerbosityQuda ( QudaVerbosity  verbosity,
const char  prefix[],
FILE *  outfile 
)

Set parameters related to status reporting.

In typical usage, this function will be called once (or not at all) just before the call to initQuda(), but it's valid to call it any number of times at any point during execution. Prior to the first time it's called, the parameters take default values as indicated below.

Parameters
verbosityDefault verbosity, ranging from QUDA_SILENT to QUDA_DEBUG_VERBOSE. Within a solver, this parameter is overridden by the "verbosity" member of QudaInvertParam. The default value is QUDA_SUMMARIZE.
prefixString to prepend to all messages from QUDA. This defaults to the empty string (""), but you may wish to specify something like "QUDA: " to distinguish QUDA's output from that of your application.
outfileFile pointer (such as stdout, stderr, or a handle returned by fopen()) where messages should be printed. The default is stdout.

Definition at line 323 of file interface_quda.cpp.

References setOutputFile(), setOutputPrefix(), and setVerbosity().

Referenced by gauge_force_test(), and init().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ staggeredPhaseQuda()

void staggeredPhaseQuda ( void *  gauge_h,
QudaGaugeParam param 
)

◆ updateGaugeFieldQuda()

void updateGaugeFieldQuda ( void *  gauge,
void *  momentum,
double  dt,
int  conj_mom,
int  exact,
QudaGaugeParam param 
)

Evolve the gauge field by step size dt, using the momentum field I.e., Evalulate U(t+dt) = e(dt pi) U(t)

Parameters
gaugeThe gauge field to be updated
momentumThe momentum field
dtThe integration step size step
conj_momWhether to conjugate the momentum matrix
exactWhether to use an exact exponential or Taylor expand
paramThe parameters of the external fields and the computation settings

Definition at line 4869 of file interface_quda.cpp.

References checkCudaError, cpuGauge, cpuMom, quda::GaugeFieldParam::create, cudaMom, errorQuda, QudaGaugeParam_s::gauge_offset, gaugePrecise, quda::LatticeFieldParam::ghostExchange, gParam, quda::GaugeFieldParam::link_type, quda::cudaGaugeField::loadCPUField(), QudaGaugeParam_s::make_resident_gauge, QudaGaugeParam_s::make_resident_mom, QudaGaugeParam_s::mom_offset, momResident, quda::GaugeFieldParam::order, quda::LatticeFieldParam::pad, profileGaugeUpdate, QUDA_ASQTAD_MOM_LINKS, QUDA_FLOAT2_GAUGE_ORDER, QUDA_GHOST_EXCHANGE_NO, QUDA_NULL_FIELD_CREATE, quda::QUDA_PROFILE_COMPUTE, quda::QUDA_PROFILE_D2H, quda::QUDA_PROFILE_FREE, quda::QUDA_PROFILE_H2D, quda::QUDA_PROFILE_INIT, quda::QUDA_PROFILE_TOTAL, QUDA_RECONSTRUCT_10, QUDA_RECONSTRUCT_NO, QUDA_SU3_LINKS, QUDA_TIFR_GAUGE_ORDER, QUDA_TIFR_PADDED_GAUGE_ORDER, quda::GaugeFieldParam::reconstruct, QudaGaugeParam_s::reconstruct, QudaGaugeParam_s::return_result_gauge, quda::GaugeFieldParam::site_offset, quda::GaugeFieldParam::site_size, QudaGaugeParam_s::site_size, quda::updateGaugeField(), QudaGaugeParam_s::use_resident_gauge, and QudaGaugeParam_s::use_resident_mom.

Referenced by update_gauge_field_quda_().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ updateMultigridQuda()

void updateMultigridQuda ( void *  mg_instance,
QudaMultigridParam param 
)

Updates the multigrid preconditioner for the new gauge / clover field.

Parameters
mg_instancePointer to instance of multigrid_solver
paramContains all metadata regarding host and device storage and solver parameters

Definition at line 2645 of file interface_quda.cpp.

References checkGauge(), quda::Dirac::create(), getVerbosity(), QudaMultigridParam_s::invert_param, quda::multigrid_solver::mg, param, popVerbosity(), printfQuda, profileInvert, profilerStart(), profilerStop(), pushVerbosity(), QUDA_DIRECT_PC_SOLVE, QUDA_NORMOP_PC_SOLVE, quda::QUDA_PROFILE_PREAMBLE, quda::QUDA_PROFILE_TOTAL, QUDA_SUMMARIZE, quda::MG::reset(), quda::saveTuneCache(), quda::setDiracPreParam(), quda::setDiracSloppyParam(), setOutputPrefix(), QudaMultigridParam_s::smoother_solve_type, QudaInvertParam_s::solve_type, and QudaInvertParam_s::verbosity.

Referenced by main().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ updateR()

void updateR ( )

update the radius for halos.

This should only be needed for automated testing when different partitioning is applied within a single run.

Definition at line 674 of file interface_quda.cpp.

References commDimPartitioned(), R, and redundant_comms.

Referenced by quda::CG::blocksolve(), quda::CG::operator()(), quda::PreconCG::operator()(), quda::BiCGstab::operator()(), quda::MultiShiftCG::operator()(), quda::reliable(), quda::CACG::reliable(), StaggeredDslashTest::SetUp(), and DslashTest::SetUp().

Here is the call graph for this function:
Here is the caller graph for this function: