su2hmc
|
Integrators for the HMC. More...
Go to the source code of this file.
Functions | |
int | Gauge_Update (const double d, double *pp, Complex *ut[2], Complex_f *ut_f[2]) |
Gauge update for the integration step of the HMC. | |
int | Momentum_Update (const double d, const double *dSdpi, double *pp) |
Wrapper for the momentum update during the integration step of the HMC. | |
int | Leapfrog (Complex *ut[2], Complex_f *ut_f[2], Complex *X0, Complex *X1, Complex *Phi, double *dk[2], float *dk_f[2], double *dSdpi, double *pp, int *iu, int *id, Complex *gamval, Complex_f *gamval_f, int *gamin, Complex jqq, float beta, float akappa, int stepl, float dt, double *ancg, int *itot, float proby) |
Leapfrog integrator. Each trajectory step takes the form of p->p+dt/2,u->u+dt,p->p+dt/2 In practice this is implemented for the entire trajectory as p->p+dt/2,u->u+dt,p->p+dt,u->u+dt,p->p+dt,...p->p+dt/2,u->u+dt,p->p+dt/2. | |
int | OMF2 (Complex *ut[2], Complex_f *ut_f[2], Complex *X0, Complex *X1, Complex *Phi, double *dk[2], float *dk_f[2], double *dSdpi, double *pp, int *iu, int *id, Complex *gamval, Complex_f *gamval_f, int *gamin, Complex jqq, float beta, float akappa, int stepl, float dt, double *ancg, int *itot, float proby) |
OMF second order five step integrator. | |
int | OMF4 (Complex *ut[2], Complex_f *ut_f[2], Complex *X0, Complex *X1, Complex *Phi, double *dk[2], float *dk_f[2], double *dSdpi, double *pp, int *iu, int *id, Complex *gamval, Complex_f *gamval_f, int *gamin, Complex jqq, float beta, float akappa, int stepl, float dt, double *ancg, int *itot, float proby) |
OMF fourth order eleven step integrator. | |
Integrators for the HMC.
Definition in file integrate.h.
Gauge update for the integration step of the HMC.
dt | Gauge step size |
pp | Momentum field |
ut | Double precision gauge fields |
ut_f | Single precision gauge fields |
Definition at line 7 of file integrate.c.
References Complex, Complex_f, kvol, nadj, ndim, Reunitarise(), and Trial_Exchange().
int Leapfrog | ( | Complex * | ut[2], |
Complex_f * | ut_f[2], | ||
Complex * | X0, | ||
Complex * | X1, | ||
Complex * | Phi, | ||
double * | dk[2], | ||
float * | dk_f[2], | ||
double * | dSdpi, | ||
double * | pp, | ||
int * | iu, | ||
int * | id, | ||
Complex * | gamval, | ||
Complex_f * | gamval_f, | ||
int * | gamin, | ||
Complex | jqq, | ||
float | beta, | ||
float | akappa, | ||
int | stepl, | ||
float | dt, | ||
double * | ancg, | ||
int * | itot, | ||
float | proby ) |
Leapfrog integrator. Each trajectory step takes the form of p->p+dt/2,u->u+dt,p->p+dt/2 In practice this is implemented for the entire trajectory as p->p+dt/2,u->u+dt,p->p+dt,u->u+dt,p->p+dt,...p->p+dt/2,u->u+dt,p->p+dt/2.
ut | Double precision colour fields |
ut_f | Single precision colour fields |
X0 | Up/down partitioned pseudofermion field |
X1 | Holder for the partitioned fermion field, then the conjugate gradient output |
Phi | Pseudofermion field |
dk | \(\left(1+\gamma_0\right)e^{-\mu}\) and \(\left(1-\gamma_0\right)e^\mu\) |
dk_f | \(\left(1+\gamma_0\right)e^{-\mu}\) and \(\left(1-\gamma_0\right)e^\mu\) float |
dSdpi | The force |
pp | Momentum field |
iu,id | Lattice indices |
gamin | Gamma indices |
gamval | Double precision gamma matrices rescaled by kappa |
gamval_f | Single precision gamma matrices rescaled by kappa |
jqq | Diquark source |
akappa | Hopping parameter |
beta | Inverse gauge coupling |
stepl | Steps per trajectory |
dt | Step size |
ancg | Counter for average conjugate gradient iterations |
itot | Total average conjugate gradient iterations |
proby | Termination probability for random trajectory length |
Definition at line 65 of file integrate.c.
References Complex, Complex_f, Force(), Gauge_Update(), Momentum_Update(), rank, and rescgg.
|
inline |
Wrapper for the momentum update during the integration step of the HMC.
dt | Gauge step size |
pp | Momentum field |
dSdpi | Force field |
Definition at line 53 of file integrate.c.
References kmom.
int OMF2 | ( | Complex * | ut[2], |
Complex_f * | ut_f[2], | ||
Complex * | X0, | ||
Complex * | X1, | ||
Complex * | Phi, | ||
double * | dk[2], | ||
float * | dk_f[2], | ||
double * | dSdpi, | ||
double * | pp, | ||
int * | iu, | ||
int * | id, | ||
Complex * | gamval, | ||
Complex_f * | gamval_f, | ||
int * | gamin, | ||
Complex | jqq, | ||
float | beta, | ||
float | akappa, | ||
int | stepl, | ||
float | dt, | ||
double * | ancg, | ||
int * | itot, | ||
float | proby ) |
OMF second order five step integrator.
ut | Double precision colour fields |
ut_f | Single precision colour fields |
X0 | Up/down partitioned pseudofermion field |
X1 | Holder for the partitioned fermion field, then the conjugate gradient output |
Phi | Pseudofermion field |
dk | \(\left(1+\gamma_0\right)e^{-\mu}\) and \(\left(1-\gamma_0\right)e^\mu\) |
dk_f | \(\left(1+\gamma_0\right)e^{-\mu}\) and \(\left(1-\gamma_0\right)e^\mu\) float |
dSdpi | The force |
pp | Momentum field |
iu,id | Lattice indices |
gamin | Gamma indices |
gamval | Double precision gamma matrices rescaled by kappa |
gamval_f | Single precision gamma matrices rescaled by kappa |
jqq | Diquark source |
akappa | Hopping parameter |
beta | Inverse gauge coupling |
stepl | Steps per trajectory |
dt | Step size |
ancg | Counter for average conjugate gradient iterations |
itot | Total average conjugate gradient iterations |
proby | Termination probability for random trajectory length |
Definition at line 116 of file integrate.c.
References Complex, Complex_f, Force(), Gauge_Update(), Momentum_Update(), rank, and rescgg.
int OMF4 | ( | Complex * | ut[2], |
Complex_f * | ut_f[2], | ||
Complex * | X0, | ||
Complex * | X1, | ||
Complex * | Phi, | ||
double * | dk[2], | ||
float * | dk_f[2], | ||
double * | dSdpi, | ||
double * | pp, | ||
int * | iu, | ||
int * | id, | ||
Complex * | gamval, | ||
Complex_f * | gamval_f, | ||
int * | gamin, | ||
Complex | jqq, | ||
float | beta, | ||
float | akappa, | ||
int | stepl, | ||
float | dt, | ||
double * | ancg, | ||
int * | itot, | ||
float | proby ) |
OMF fourth order eleven step integrator.
ut | Double precision colour fields |
ut_f | Single precision colour fields |
X0 | Up/down partitioned pseudofermion field |
X1 | Holder for the partitioned fermion field, then the conjugate gradient output |
Phi | Pseudofermion field |
dk | \(\left(1+\gamma_0\right)e^{-\mu}\) and \(\left(1-\gamma_0\right)e^\mu\) |
dk_f | \(\left(1+\gamma_0\right)e^{-\mu}\) and \(\left(1-\gamma_0\right)e^\mu\) float |
dSdpi | The force |
pp | Momentum field |
iu,id | Lattice indices |
gamin | Gamma indices |
gamval | Double precision gamma matrices rescaled by kappa |
gamval_f | Single precision gamma matrices rescaled by kappa |
jqq | Diquark source |
akappa | Hopping parameter |
beta | Inverse gauge coupling |
stepl | Steps per trajectory |
dt | Step size |
ancg | Counter for average conjugate gradient iterations |
itot | Total average conjugate gradient iterations |
proby | Termination probability for random trajectory length |
Momentum updates
Outer updates depend on r1. We have two of these, doubled for between full steps
Middle updates depend on r3
Inner updates depend on r1 and r3
Gauge updates. These depend on r2 and r4
Outer gauge update depends on r2
Middle gauge update depends on r4
Inner gauge update depends on r2 and r4
Definition at line 194 of file integrate.c.
References Complex, Complex_f, Force(), Gauge_Update(), Momentum_Update(), rank, and rescgg.