QUDA
v1.1.0
A library for QCD on GPUs
|
Go to the source code of this file.
Classes | |
struct | QudaMILCSiteArg_t |
struct | QudaInvertArgs_t |
struct | QudaEigArgs_t |
struct | QudaLayout_t |
struct | QudaInitArgs_t |
struct | QudaHisqParams_t |
struct | QudaFatLinkArgs_t |
Functions | |
void | qudaSetMPICommHandle (void *mycomm) |
void | qudaInit (QudaInitArgs_t input) |
void | qudaSetLayout (QudaLayout_t layout) |
void | qudaFinalize () |
void * | qudaAllocatePinned (size_t bytes) |
void | qudaFreePinned (void *ptr) |
void * | qudaAllocateManaged (size_t bytes) |
void | qudaFreeManaged (void *ptr) |
void | qudaHisqParamsInit (QudaHisqParams_t hisq_params) |
void | qudaLoadKSLink (int precision, QudaFatLinkArgs_t fatlink_args, const double act_path_coeff[6], void *inlink, void *fatlink, void *longlink) |
void | qudaLoadUnitarizedLink (int precision, QudaFatLinkArgs_t fatlink_args, const double path_coeff[6], void *inlink, void *fatlink, void *ulink) |
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) |
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) |
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 * | qudaMultigridCreate (int external_precision, int quda_precision, double mass, QudaInvertArgs_t inv_args, const void *const milc_fatlink, const void *const milc_longlink, const char *const mg_param_file) |
void | qudaInvertMG (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 *mg_pack_ptr, int mg_rebuild_type, void *source, void *solution, double *const final_residual, double *const final_fermilab_residual, int *num_iters) |
void | qudaMultigridDestroy (void *mg_pack_ptr) |
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) |
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) |
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 | 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) |
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 | qudaLoadGaugeField (int external_precision, int quda_precision, QudaInvertArgs_t inv_args, const void *milc_link) |
void | qudaFreeGaugeField () |
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 | qudaFreeCloverField () |
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 | 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) |
void | qudaGaugeForce (int precision, int num_loop_types, double milc_loop_coeff[3], double eb3, QudaMILCSiteArg_t *arg) |
void | qudaGaugeForcePhased (int precision, int num_loop_types, double milc_loop_coeff[3], double eb3, QudaMILCSiteArg_t *arg, int phase_in) |
void | qudaUpdateU (int precision, double eps, QudaMILCSiteArg_t *arg) |
void | qudaUpdateUPhased (int precision, double eps, QudaMILCSiteArg_t *arg, int phase_in) |
void | qudaUpdateUPhasedPipeline (int precision, double eps, QudaMILCSiteArg_t *arg, int phase_in, int want_gaugepipe) |
void | qudaMomLoad (int precision, QudaMILCSiteArg_t *arg) |
void | qudaMomSave (int precision, QudaMILCSiteArg_t *arg) |
double | qudaMomAction (int precision, QudaMILCSiteArg_t *arg) |
void | qudaRephase (int prec, void *gauge, int flag, double i_mu) |
void | qudaUnitarizeSU3 (int prec, double tol, QudaMILCSiteArg_t *arg) |
void | qudaUnitarizeSU3Phased (int prec, double tol, QudaMILCSiteArg_t *arg, int phase_in) |
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 | qudaCloverTrace (void *out, void *dummy, int mu, int nu) |
void | qudaCloverDerivative (void *out, void *gauge, void *oprod, int mu, int nu, double coeff, int precision, int parity, int conjugate) |
void * | qudaCreateExtendedGaugeField (void *gauge, int geometry, int precision) |
void * | qudaResidentExtendedGaugeField (void *gauge, int geometry, int precision) |
void * | qudaCreateGaugeField (void *gauge, int geometry, int precision) |
void | qudaSaveGaugeField (void *gauge, void *inGauge) |
void | qudaDestroyGaugeField (void *gauge) |
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. More... | |
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. More... | |
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) |
void | qudaComputeOprod (int precision, int num_terms, int num_naik_terms, double **coeff, double scale, void **quark_field, void *oprod[3]) |
The header file defines the milc interface to enable easy interfacing between QUDA and the MILC software packed.
Definition in file quda_milc_interface.h.
void* qudaAllocateManaged | ( | size_t | bytes | ) |
Allocate managed memory to reduce CPU-GPU transfers
bytes | The size of the requested allocation |
Definition at line 165 of file milc_interface.cpp.
void* qudaAllocatePinned | ( | size_t | bytes | ) |
Allocate pinned memory suitable for CPU-GPU transfers
bytes | The size of the requested allocation |
Definition at line 161 of file milc_interface.cpp.
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 | ||
) |
Note this interface function has been removed. This stub remains for compatibility only.
Definition at line 313 of file milc_interface.cpp.
void qudaCloverDerivative | ( | void * | out, |
void * | gauge, | ||
void * | oprod, | ||
int | mu, | ||
int | nu, | ||
double | coeff, | ||
int | precision, | ||
int | parity, | ||
int | conjugate | ||
) |
Compute the derivative of the clover term (part of clover force computation). All the pointers here are for QUDA native device objects. The precisions of all fields must match.
out | Clover derivative field (QUDA device field, geometry = 1) |
gauge | Gauge field (extended QUDA device field, gemoetry = 4) |
oprod | Matrix field (outer product) which is multiplied by the derivative |
mu | mu direction |
nu | nu direction |
coeff | Coefficient of the clover derviative (including stepsize and clover coefficient) |
precision | Precision of the fields (2 = double, 1 = single) |
parity | Parity for which we are computing |
conjugate | Whether to make the oprod field anti-hermitian prior to multiplication |
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 | ||
) |
Compute the clover force contributions in each dimension mu given the array solution fields, and compute the resulting momentum field.
mom | Momentum matrix |
dt | Integrating step size |
x | Array of solution vectors |
p | Array of intermediate vectors |
coeff | Array of residues for each contribution |
kappa | kappa parameter |
ck | -clover_coefficient * kappa / 8 |
nvec | Number of vectors |
multiplicity | Number of fermions represented by this bilinear |
gauge | Gauge Field |
precision | Precision of the fields |
inv_args | Struct setting some solver metadata |
Definition at line 2303 of file milc_interface.cpp.
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 | ||
) |
Solve Ax=b using a Wilson-Clover operator. All fields are fields passed and returned are host (CPU) field in MILC order. This function creates the gauge and clover field from the host fields. Reliable updates are used with a reliable_delta parameter of 0.1.
external_precision | Precision of host fields passed to QUDA (2 - double, 1 - single) |
quda_precision | Precision for QUDA to use (2 - double, 1 - single) |
kappa | Kappa value |
clover_coeff | Clover coefficient |
inv_args | Struct setting some solver metadata |
target_residual | Target residual |
milc_link | Gauge field on the host |
milc_clover | Clover field on the host |
milc_clover_inv | Inverse clover on the host |
clover_coeff | Clover coefficient |
source | Right-hand side source field |
solution | Solution spinor field |
final_residual | True residual returned by the solver |
final_residual | True Fermilab residual returned by the solver |
num_iters | Number of iterations taken |
Definition at line 2492 of file milc_interface.cpp.
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 | ||
) |
Solve for multiple shifts (e.g., masses) using a Wilson-Clover operator with multi-shift CG. All fields are fields passed and returned are host (CPU) field in MILC order. This function requires that persistent gauge and clover fields have been created prior. When a pure double-precision solver is requested no reliable updates are used, else reliable updates are used with a reliable_delta parameter of 0.1.
external_precision | Precision of host fields passed to QUDA (2 - double, 1 - single) |
quda_precision | Precision for QUDA to use (2 - double, 1 - single) |
num_offsets | Number of shifts to solve for |
offset | Array of shift offset values |
kappa | Kappa value |
clover_coeff | Clover coefficient |
inv_args | Struct setting some solver metadata |
target_residual | Array of target residuals per shift |
milc_link | Ignored |
milc_clover | Ignored |
milc_clover_inv | Ignored |
clover_coeff | Clover coefficient |
source | Right-hand side source field |
solutionArray | Array of solution spinor fields |
final_residual | Array of true residuals |
num_iters | Number of iterations taken |
Definition at line 2628 of file milc_interface.cpp.
void qudaCloverTrace | ( | void * | out, |
void * | dummy, | ||
int | mu, | ||
int | nu | ||
) |
Compute the sigma trace field (part of clover force computation). All the pointers here are for QUDA native device objects. The precisions of all fields must match. This function requires that there is a persistent clover field.
out | Sigma trace field (QUDA device field, geometry = 1) |
dummy | (not used) |
mu | mu direction |
nu | nu direction |
void qudaComputeOprod | ( | int | precision, |
int | num_terms, | ||
int | num_naik_terms, | ||
double ** | coeff, | ||
double | scale, | ||
void ** | quark_field, | ||
void * | oprod[3] | ||
) |
Note this interface function has been removed. This stub remains for compatibility only.
Definition at line 322 of file milc_interface.cpp.
void* qudaCreateExtendedGaugeField | ( | void * | gauge, |
int | geometry, | ||
int | precision | ||
) |
Take a gauge field on the host, load it onto the device and extend it. Return a pointer to the extended gauge field object.
gauge | The CPU gauge field (optional - if set to 0 then the gauge field zeroed) |
geometry | The geometry of the matrix field to create (1 - scaler, 4 - vector, 6 - tensor) |
precision | The precision of the fields (2 - double, 1 - single) |
void* qudaCreateGaugeField | ( | void * | gauge, |
int | geometry, | ||
int | precision | ||
) |
Allocate a gauge (matrix) field on the device and optionally download a host gauge field.
gauge | The host gauge field (optional - if set to 0 then the gauge field zeroed) |
geometry | The geometry of the matrix field to create (1 - scaler, 4 - vector, 6 - tensor) |
precision | The precision of the field to be created (2 - double, 1 - single) |
Definition at line 2271 of file milc_interface.cpp.
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 | ||
) |
Solve Ax=b using an improved staggered operator with a domain-decomposition preconditioner. All fields are fields passed and returned are host (CPU) field in MILC order. This function requires that persistent gauge and clover fields have been created prior. This interface is experimental.
external_precision | Precision of host fields passed to QUDA (2 - double, 1 - single) |
precision | Precision for QUDA to use (2 - double, 1 - single) |
mass | Fermion mass parameter |
inv_args | Struct setting some solver metadata |
target_residual | Target residual |
target_relative_residual | Target Fermilab residual |
domain_overlap | Array specifying the overlap of the domains in each dimension |
fatlink | Fat-link field on the host |
longlink | Long-link field on the host |
source | Right-hand side source field |
solution | Solution spinor field |
final_residual | True residual |
final_relative_residual | True Fermilab residual |
num_iters | Number of iterations taken |
void qudaDestroyGaugeField | ( | void * | gauge | ) |
Reinterpret gauge as a pointer to cudaGaugeField and call destructor.
gauge | Gauge field to be freed |
Definition at line 2292 of file milc_interface.cpp.
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 | ||
) |
Apply the improved staggered operator to a field. All fields passed and returned are host (CPU) field in MILC order.
external_precision | Precision of host fields passed to QUDA (2 - double, 1 - single) |
quda_precision | Precision for QUDA to use (2 - double, 1 - single) |
inv_args | Struct setting some solver metadata |
milc_fatlink | Fat-link field on the host |
milc_longlink | Long-link field on the host |
source | Right-hand side source field |
solution | Solution spinor field |
Definition at line 1066 of file milc_interface.cpp.
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 | ||
) |
Solve for a system with many RHS using using a Wilson-Clover operator. The solving procedure consists of two computation phases : 1) incremental pahse : call eigCG solver to accumulate low eigenmodes 2) deflation phase : use computed eigenmodes to deflate a regular CG All fields are fields passed and returned are host (CPU) field in MILC order. This function requires that persistent gauge and clover fields have been created prior.
external_precision | Precision of host fields passed to QUDA (2 - double, 1 - single) |
quda_precision | Precision for QUDA to use (2 - double, 1 - single) |
kappa | Kappa value |
clover_coeff | Clover coefficient |
inv_args | Struct setting some solver metadata |
target_residual | Target residual |
milc_link | Gauge field on the host |
milc_clover | Clover field on the host |
milc_clover_inv | Inverse clover on the host |
clover_coeff | Clover coefficient |
source | Right-hand side source field |
solution | Solution spinor field |
eig_args | contains info about deflation space |
rhs_idx | bookkeep current rhs |
last_rhs_flag | is this the last rhs to solve? |
final_residual | Array of true residuals |
final_relative_residual | Array of true Fermilab residuals |
num_iters | Number of iterations taken |
Definition at line 2548 of file milc_interface.cpp.
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 | ||
) |
Solve for a system with many RHS using an improved staggered operator. The solving procedure consists of two computation phases : 1) incremental pahse : call eigCG solver to accumulate low eigenmodes 2) deflation phase : use computed eigenmodes to deflate a regular CG All fields are fields passed and returned are host (CPU) field in MILC order. This function requires that persistent gauge and clover fields have been created prior.
external_precision | Precision of host fields passed to QUDA (2 - double, 1 - single) |
precision | Precision for QUDA to use (2 - double, 1 - single) |
num_offsets | Number of shifts to solve for |
offset | Array of shift offset values |
inv_args | Struct setting some solver metadata |
target_residual | Array of target residuals per shift |
target_relative_residual | Array of target Fermilab residuals per shift |
milc_fatlink | Fat-link field on the host |
milc_longlink | Long-link field on the host |
source | Right-hand side source field |
solution | Array of solution spinor fields |
eig_args | contains info about deflation space |
rhs_idx | bookkeep current rhs |
last_rhs_flag | is this the last rhs to solve? |
final_residual | Array of true residuals |
final_relative_residual | Array of true Fermilab residuals |
num_iters | Number of iterations taken |
Definition at line 1188 of file milc_interface.cpp.
void qudaFinalize | ( | ) |
Destroy the QUDA context.
Definition at line 106 of file milc_interface.cpp.
void qudaFreeCloverField | ( | ) |
Free the clover field allocated in QUDA.
Definition at line 2480 of file milc_interface.cpp.
void qudaFreeGaugeField | ( | ) |
Free the gauge field allocated in QUDA.
Definition at line 2435 of file milc_interface.cpp.
void qudaFreeManaged | ( | void * | ptr | ) |
Free managed memory
ptr | Pointer to memory to be free |
Definition at line 167 of file milc_interface.cpp.
void qudaFreePinned | ( | void * | ptr | ) |
Free pinned memory
ptr | Pointer to memory to be free |
Definition at line 163 of file milc_interface.cpp.
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.
[in] | precision,1 | for single precision else for double precision |
[in] | gauge_dir,3 | for Coulomb gauge fixing, other for Landau gauge fixing |
[in] | Nsteps,maximum | number of steps to perform gauge fixing |
[in] | verbose_interval,print | gauge fixing info when iteration count is a multiple of this |
[in] | alpha,gauge | fixing parameter of the method, most common value is 0.08 |
[in] | autotune,1 | to autotune the method, i.e., if the Fg inverts its tendency we decrease the alpha value |
[in] | tolerance,torelance | value 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,0 | for MILC criterium and 1 to use the theta value |
[in,out] | milc_sitelink,MILC | gauge field to be fixed |
Definition at line 2730 of file milc_interface.cpp.
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.
[in] | precision,1 | for single precision else for double precision |
[in] | gauge_dir,3 | for Coulomb gauge fixing, other for Landau gauge fixing |
[in] | Nsteps,maximum | number of steps to perform gauge fixing |
[in] | verbose_interval,print | gauge fixing info when iteration count is a multiple of this |
[in] | relax_boost,gauge | fixing parameter of the overrelaxation method, most common value is 1.5 or 1.7. |
[in] | tolerance,torelance | value 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,reunitarize | gauge field when iteration count is a multiple of this |
[in] | stopWtheta,0 | for MILC criterium and 1 to use the theta value |
[in,out] | milc_sitelink,MILC | gauge field to be fixed |
Definition at line 2711 of file milc_interface.cpp.
void qudaGaugeForce | ( | int | precision, |
int | num_loop_types, | ||
double | milc_loop_coeff[3], | ||
double | eb3, | ||
QudaMILCSiteArg_t * | arg | ||
) |
Compute the gauge force and update the momentum field. All fields here are CPU fields in MILC order, and their precisions should match.
precision | The precision of the field (2 - double, 1 - single) |
num_loop_types | 1, 2 or 3 |
milc_loop_coeff | Coefficients of the different loops in the Symanzik action |
eb3 | The integration step size (for MILC this is dt*beta/3) |
arg | Metadata for MILC's internal site struct array |
Definition at line 655 of file milc_interface.cpp.
void qudaGaugeForcePhased | ( | int | precision, |
int | num_loop_types, | ||
double | milc_loop_coeff[3], | ||
double | eb3, | ||
QudaMILCSiteArg_t * | arg, | ||
int | phase_in | ||
) |
Compute the gauge force and update the momentum field. All fields here are CPU fields in MILC order, and their precisions should match.
precision | The precision of the field (2 - double, 1 - single) |
num_loop_types | 1, 2 or 3 |
milc_loop_coeff | Coefficients of the different loops in the Symanzik action |
eb3 | The integration step size (for MILC this is dt*beta/3) |
arg | Metadata for MILC's internal site struct array |
phase_in | whether staggered phases are applied |
Definition at line 550 of file milc_interface.cpp.
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 | ||
) |
Compute the fermion force for the HISQ quark action. All fields are host fields in MILC order, and the precision of these fields must match.
precision | The precision of the fields |
num_terms | The number of quark fields |
num_naik_terms | The number of naik contributions |
dt | Integrating step size |
coeff | The coefficients multiplying the fermion fields in the outer product |
quark_field | The input fermion field. |
level2_coeff | The coefficients for the second level of smearing in the quark action. |
fat7_coeff | The coefficients for the first level of smearing (fat7) in the quark action. |
w_link | Unitarized link variables obtained by applying fat7 smearing and unitarization to the original links. |
v_link | Fat7 link variables. |
u_link | SU(3) think link variables. |
milc_momentum | The momentum contribution from the quark action. |
Definition at line 283 of file milc_interface.cpp.
void qudaHisqParamsInit | ( | QudaHisqParams_t | hisq_params | ) |
Set the algorithms to use for HISQ fermion calculations, e.g., SVD parameters for reunitarization.
hisq_params | Meta data desribing the algorithms to use for HISQ fermions |
Definition at line 169 of file milc_interface.cpp.
void qudaInit | ( | QudaInitArgs_t | input | ) |
Initialize the QUDA context.
input | Meta data for the QUDA context |
Definition at line 96 of file milc_interface.cpp.
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 | ||
) |
Solve Ax=b for an improved staggered operator. All fields are fields passed and returned are host (CPU) field in MILC order. This function requires that persistent gauge and clover fields have been created prior. This interface is experimental.
external_precision | Precision of host fields passed to QUDA (2 - double, 1 - single) |
quda_precision | Precision for QUDA to use (2 - double, 1 - single) |
mass | Fermion mass parameter |
inv_args | Struct setting some solver metadata |
target_residual | Target residual |
target_relative_residual | Target Fermilab residual |
milc_fatlink | Fat-link field on the host |
milc_longlink | Long-link field on the host |
source | Right-hand side source field |
solution | Solution spinor field |
final_residual | True residual |
final_relative_residual | True Fermilab residual |
num_iters | Number of iterations taken |
Definition at line 985 of file milc_interface.cpp.
void qudaInvertMG | ( | 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 * | mg_pack_ptr, | ||
int | mg_rebuild_type, | ||
void * | source, | ||
void * | solution, | ||
double *const | final_residual, | ||
double *const | final_fermilab_residual, | ||
int * | num_iters | ||
) |
Solve Ax=b for an improved staggered operator using MG. All fields are fields passed and returned are host (CPU) field in MILC order. This function requires that persistent gauge and clover fields have been created prior. It also requires a multigrid parameter built from qudaSetupMultigrid This interface is experimental.
external_precision | Precision of host fields passed to QUDA (2 - double, 1 - single) |
quda_precision | Precision for QUDA to use (2 - double, 1 - single) |
mass | Fermion mass parameter |
inv_args | Struct setting some solver metadata |
target_residual | Target residual |
target_relative_residual | Target Fermilab residual |
milc_fatlink | Fat-link field on the host |
milc_longlink | Long-link field on the host |
mg_pack_ptr | MG preconditioner structure created by qudaSetupMultigrid |
mg_rebuild_type | whether to do a full (1) or thin (0) MG rebuild |
source | Right-hand side source field |
solution | Solution spinor field |
final_residual | True residual |
final_relative_residual | True Fermilab residual |
num_iters | Number of iterations taken |
Definition at line 2134 of file milc_interface.cpp.
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 | ||
) |
Solve Ax=b for an improved staggered operator with many right hand sides. All fields are fields passed and returned are host (CPU) field in MILC order. This function requires that persistent gauge and clover fields have been created prior. This interface is experimental.
external_precision | Precision of host fields passed to QUDA (2 - double, 1 - single) |
quda_precision | Precision for QUDA to use (2 - double, 1 - single) |
mass | Fermion mass parameter |
inv_args | Struct setting some solver metadata |
target_residual | Target residual |
target_relative_residual | Target Fermilab residual |
milc_fatlink | Fat-link field on the host |
milc_longlink | Long-link field on the host |
source | array of right-hand side source fields |
solution | array of solution spinor fields |
final_residual | True residual |
final_relative_residual | True Fermilab residual |
num_iters | Number of iterations taken |
num_src | Number of source fields |
Definition at line 1116 of file milc_interface.cpp.
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 | ||
) |
Load the clover field and its inverse from the host. If null pointers are passed, the clover field and / or its inverse will be computed dynamically from the resident gauge field.
external_precision | Precision of host fields passed to QUDA (2 - double, 1 - single) |
quda_precision | Precision for QUDA to use (2 - double, 1 - single) |
inv_args | Meta data |
milc_clover | Pointer to host clover field. If 0 then the clover field is computed dynamically within QUDA. |
milc_clover_inv | Pointer to host inverse clover field. If 0 then the inverse if computed dynamically within QUDA. |
solution_type | The type of solution required (mat, matpc) |
solve_type | The solve type to use (normal/direct/preconditioning) |
clover_coeff | Clover coefficient |
compute_trlog | Whether to compute the trlog of the clover field when inverting |
Array | for storing the trlog (length two, one for each parity) |
void qudaLoadGaugeField | ( | int | external_precision, |
int | quda_precision, | ||
QudaInvertArgs_t | inv_args, | ||
const void * | milc_link | ||
) |
Load the gauge field from the host.
external_precision | Precision of host fields passed to QUDA (2 - double, 1 - single) |
quda_precision | Precision for QUDA to use (2 - double, 1 - single) |
inv_args | Meta data |
milc_link | Base pointer to host gauge field (regardless of dimensionality) |
Definition at line 2422 of file milc_interface.cpp.
void qudaLoadKSLink | ( | int | precision, |
QudaFatLinkArgs_t | fatlink_args, | ||
const double | act_path_coeff[6], | ||
void * | inlink, | ||
void * | fatlink, | ||
void * | longlink | ||
) |
Compute the fat and long links using the input gauge field. All fields passed here are host fields, that must be preallocated. The precision of all fields must match.
precision | The precision of the fields |
fatlink_args | Meta data for the algorithms to deploy |
act_path_coeff | Array of coefficients for each path in the action |
inlink | Host gauge field used for input |
fatlink | Host fat-link field that is computed |
longlink | Host long-link field that is computed |
Definition at line 239 of file milc_interface.cpp.
void qudaLoadUnitarizedLink | ( | int | precision, |
QudaFatLinkArgs_t | fatlink_args, | ||
const double | path_coeff[6], | ||
void * | inlink, | ||
void * | fatlink, | ||
void * | ulink | ||
) |
Compute the fat links and unitzarize using the input gauge field. All fields passed here are host fields, that must be preallocated. The precision of all fields must match.
precision | The precision of the fields |
fatlink_args | Meta data for the algorithms to deploy |
path_coeff | Array of coefficients for each path in the action |
inlink | Host gauge field used for input |
fatlink | Host fat-link field that is computed |
ulink | Host unitarized field that is computed |
Definition at line 263 of file milc_interface.cpp.
double qudaMomAction | ( | int | precision, |
QudaMILCSiteArg_t * | arg | ||
) |
Evaluate the momentum contribution to the Hybrid Monte Carlo action. MILC convention is applied, subtracting 4.0 from each momentum matrix to increase stability.
precision | Precision of the field (2 - double, 1 - single) |
arg | Metadata for MILC's internal site struct array |
Definition at line 467 of file milc_interface.cpp.
void qudaMomLoad | ( | int | precision, |
QudaMILCSiteArg_t * | arg | ||
) |
Download the momentum from MILC and place into QUDA's resident momentum field. The source momentum field can either be as part of a MILC site struct (QUDA_MILC_SITE_GAUGE_ORDER) or as a separate field (QUDA_MILC_GAUGE_ORDER).
precision | Precision of the field (2 - double, 1 - single) |
arg | Metadata for MILC's internal site struct array |
Definition at line 426 of file milc_interface.cpp.
void qudaMomSave | ( | int | precision, |
QudaMILCSiteArg_t * | arg | ||
) |
Upload the momentum to MILC from QUDA's resident momentum field. The destination momentum field can either be as part of a MILC site struct (QUDA_MILC_SITE_GAUGE_ORDER) or as a separate field (QUDA_MILC_GAUGE_ORDER).
precision | Precision of the field (2 - double, 1 - single) |
arg | Metadata for MILC's internal site struct array |
Definition at line 447 of file milc_interface.cpp.
void* qudaMultigridCreate | ( | int | external_precision, |
int | quda_precision, | ||
double | mass, | ||
QudaInvertArgs_t | inv_args, | ||
const void *const | milc_fatlink, | ||
const void *const | milc_longlink, | ||
const char *const | mg_param_file | ||
) |
Prepare a staggered/HISQ multigrid solve with given fat and long links. All fields passed are host (CPU) fields in MILC order. This function requires persistent gauge fields. This interface is experimental.
external_precision | Precision of host fields passed to QUDA (2 - double, 1 - single) |
quda_precision | Precision for QUDA to use (2 - double, 1 - single) |
mass | Fermion mass parameter |
inv_args | Struct setting some solver metadata; required for tadpole, naik coeff |
milc_fatlink | Fat-link field on the host |
milc_longlink | Long-link field on the host |
mg_param_file | Path to an input text file describing the MG solve, to be documented on QUDA wiki |
Definition at line 2073 of file milc_interface.cpp.
void qudaMultigridDestroy | ( | void * | mg_pack_ptr | ) |
Clean up a staggered/HISQ multigrid object, freeing all internal fields and otherwise allocated memory.
mg_pack_ptr | Void pointer mapping to the multigrid structure returned by qudaSetupMultigrid |
Definition at line 2255 of file milc_interface.cpp.
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 | ||
) |
Solve for multiple shifts (e.g., masses) using an improved staggered operator. All fields are fields passed and returned are host (CPU) field in MILC order. This function requires that persistent gauge and clover fields have been created prior. When a pure double-precision solver is requested no reliable updates are used, else reliable updates are used with a reliable_delta parameter of 0.1.
external_precision | Precision of host fields passed to QUDA (2 - double, 1 - single) |
precision | Precision for QUDA to use (2 - double, 1 - single) |
num_offsets | Number of shifts to solve for |
offset | Array of shift offset values |
inv_args | Struct setting some solver metadata |
target_residual | Array of target residuals per shift |
target_relative_residual | Array of target Fermilab residuals per shift |
milc_fatlink | Fat-link field on the host |
milc_longlink | Long-link field on the host |
source | Right-hand side source field |
solutionArray | Array of solution spinor fields |
final_residual | Array of true residuals |
final_relative_residual | Array of true Fermilab residuals |
num_iters | Number of iterations taken |
void qudaRephase | ( | int | prec, |
void * | gauge, | ||
int | flag, | ||
double | i_mu | ||
) |
Apply the staggered phase factors to the gauge field. If the imaginary chemical potential is non-zero then the phase factor exp(imu/T) will be applied to the links in the temporal direction.
prec | Precision of the gauge field |
gauge_h | The gauge field |
flag | Whether to apply to remove the staggered phase |
i_mu | Imaginary chemical potential |
Definition at line 376 of file milc_interface.cpp.
void* qudaResidentExtendedGaugeField | ( | void * | gauge, |
int | geometry, | ||
int | precision | ||
) |
Take the QUDA resident gauge field and extend it. Return a pointer to the extended gauge field object.
gauge | The CPU gauge field (optional - if set to 0 then the gauge field zeroed) |
geometry | The geometry of the matrix field to create (1 - scaler, 4 - vector, 6 - tensor) |
precision | The precision of the fields (2 - double, 1 - single) |
void qudaSaveGaugeField | ( | void * | gauge, |
void * | inGauge | ||
) |
Copy the QUDA gauge (matrix) field on the device to the CPU
outGauge | Pointer to the host gauge field |
inGauge | Pointer to the device gauge field (QUDA device field) |
Definition at line 2282 of file milc_interface.cpp.
void qudaSetLayout | ( | QudaLayout_t | layout | ) |
Set set the local dimensions and machine topology for QUDA to use
layout | Struct defining local dimensions and machine topology |
Definition at line 129 of file milc_interface.cpp.
void qudaSetMPICommHandle | ( | void * | mycomm | ) |
Optional: Set the MPI Comm Handle if it is not MPI_COMM_WORLD
input | Pointer to an MPI_Comm handle, static cast as a void * |
Definition at line 94 of file milc_interface.cpp.
void qudaUnitarizeSU3 | ( | int | prec, |
double | tol, | ||
QudaMILCSiteArg_t * | arg | ||
) |
Project the input field on the SU(3) group. If the target tolerance is not met, this routine will give a runtime error.
prec | Precision of the gauge field |
tol | The tolerance to which we iterate |
arg | Metadata for MILC's internal site struct array |
Definition at line 423 of file milc_interface.cpp.
void qudaUnitarizeSU3Phased | ( | int | prec, |
double | tol, | ||
QudaMILCSiteArg_t * | arg, | ||
int | phase_in | ||
) |
Project the input field on the SU(3) group. If the target tolerance is not met, this routine will give a runtime error.
prec | Precision of the gauge field |
tol | The tolerance to which we iterate |
arg | Metadata for MILC's internal site struct array |
phase_in | whether staggered phases are applied |
Definition at line 392 of file milc_interface.cpp.
void qudaUpdateU | ( | int | precision, |
double | eps, | ||
QudaMILCSiteArg_t * | arg | ||
) |
Evolve the gauge field by step size dt, using the momentum field I.e., Evalulate U(t+dt) = e(dt pi) U(t). All fields are CPU fields in MILC order.
precision | Precision of the field (2 - double, 1 - single) |
dt | The integration step size step |
arg | Metadata for MILC's internal site struct array |
Definition at line 374 of file milc_interface.cpp.
void qudaUpdateUPhased | ( | int | precision, |
double | eps, | ||
QudaMILCSiteArg_t * | arg, | ||
int | phase_in | ||
) |
Evolve the gauge field by step size dt, using the momentum field I.e., Evalulate U(t+dt) = e(dt pi) U(t). All fields are CPU fields in MILC order.
precision | Precision of the field (2 - double, 1 - single) |
dt | The integration step size step |
arg | Metadata for MILC's internal site struct array |
phase_in | whether staggered phases are applied |
Definition at line 369 of file milc_interface.cpp.
void qudaUpdateUPhasedPipeline | ( | int | precision, |
double | eps, | ||
QudaMILCSiteArg_t * | arg, | ||
int | phase_in, | ||
int | want_gaugepipe | ||
) |
Evolve the gauge field by step size dt, using the momentum field I.e., Evalulate U(t+dt) = e(dt pi) U(t). All fields are CPU fields in MILC order.
precision | Precision of the field (2 - double, 1 - single) |
dt | The integration step size step |
arg | Metadata for MILC's internal site struct array |
phase_in | whether staggered phases are applied |
want_gaugepipe | whether to enabled QUDA gaugepipe for HMC |
Definition at line 328 of file milc_interface.cpp.