su2hmc
Loading...
Searching...
No Matches
integrate.h
Go to the documentation of this file.
1
5#ifndef INTEGRATE_H
6#define INTEGRATE_H
7#include <stdbool.h>
8#include <random.h>
9#include <sizes.h>
10
11#if (defined __cplusplus)
12extern "C"
13{
14#endif
25 int Gauge_Update(const double d, double *pp, Complex *ut[2],Complex_f *ut_f[2]);
35 int Momentum_Update(const double d,const double *dSdpi, double *pp);
65int Leapfrog(Complex *ut[2],Complex_f *ut_f[2],Complex *X0,Complex *X1, Complex *Phi,double *dk[2],float *dk_f[2],
66 double *dSdpi,double *pp, int *iu,int *id, Complex *gamval, Complex_f *gamval_f, int *gamin, Complex jqq,
67 float beta, float akappa, int stepl, float dt, double *ancg, int *itot, float proby);
95int OMF2(Complex *ut[2],Complex_f *ut_f[2],Complex *X0,Complex *X1, Complex *Phi,double *dk[2],float *dk_f[2],
96 double *dSdpi,double *pp, int *iu,int *id, Complex *gamval, Complex_f *gamval_f, int *gamin, Complex jqq,
97 float beta, float akappa, int stepl, float dt, double *ancg, int *itot, float proby);
125int OMF4(Complex *ut[2],Complex_f *ut_f[2],Complex *X0,Complex *X1, Complex *Phi,double *dk[2],float *dk_f[2],
126 double *dSdpi,double *pp, int *iu,int *id, Complex *gamval, Complex_f *gamval_f, int *gamin, Complex jqq,
127 float beta, float akappa, int stepl, float dt, double *ancg, int *itot, float proby);
128 //CUDA Calling functions
129#ifdef __NVCC__
130 void cuGauge_Update(const double d, double *pp, Complex *u11t, Complex *u12t, dim3 dimGrid, dim3 dimBlock);
131#endif
132
133#if (defined __cplusplus)
134}
135#endif
136//Actual CUDA
137#ifdef __CUDACC__
138//Update Gauge fields
139__global__ void cuGauge_Update(const double d, double *pp, Complex *u11t, Complex *u12t, int mu);
140#endif
141#endif
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.
Definition integrate.c:116
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 t...
Definition integrate.c:65
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.
Definition integrate.c:194
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.
Definition integrate.c:7
int Momentum_Update(const double d, const double *dSdpi, double *pp)
Wrapper for the momentum update during the integration step of the HMC.
Definition integrate.c:53
Header for random number configuration.
Defines the constants of the code and other parameters for loop dimensions. Each subroutine includes ...
#define Complex
Double precision complex number.
Definition sizes.h:58
#define Complex_f
Single precision complex number.
Definition sizes.h:56