su2hmc
Loading...
Searching...
No Matches
Clover
Collaboration diagram for Clover:

Topics

 Clover Force
 Clover Multiplication routines

Functions

void ByGenLeft (Complex_f a[nc], const unsigned short gen)
 Multiply leaf (or part of one) by generator from left.
void ByGenRight (Complex_f a[nc], const unsigned short gen)
 Multiply leaf (or part of one) by generator from right.
void Clover_SU2plaq (Complex_f *ut[nc], Complex_f Leaves[nc], unsigned int *iu, int i, int mu, int nu)
 Calculates the SU2 plaquette at site i in the \(\mu--\nu\) direction.
void Half_Leaves (Complex_f *hLeaves[2], Complex_f *ut[2], unsigned int *iu, unsigned int *id, const unsigned short mu, const unsigned short nu)
 Calculates the products of the first two links in a plaquette.
void Leaf (Complex_f Leaves[nc], Complex_f *ut[nc], unsigned int *iu, unsigned int *id, unsigned int i, const unsigned short mu, const unsigned short nu, const unsigned short leaf)
 Calculates a leaf for a clover term.
void Clover (Complex_f *clover[2], Complex_f *ut[2], unsigned int *iu, unsigned int *id)
 Calculates the clovers in all directions at all sites.
int Init_clover (Complex **sigval, Complex_f **sigval_f, unsigned short **sigin, float c_sw)
 Initialise values needed for the clover terms.
void Clover_free (Complex_f *clover[nc])
 Free's memory used for clover terms and leaves.
int cuClover (Complex_f *clover[nc], Complex_f *ut[nc], unsigned int *iu, unsigned int *id)
 CUDA wrapper for calculating the clovers in all directions at all sites \( F_{\mu\nu}(n)=\frac{-i}{8a^2}\left(Q_{\mu\nu}(n)-Q_{\nu\mu}(n)\right)\).
void Half_Leaf (Complex_f Leaves[nc], Complex_f *ut[nc], Complex_f a[nc], unsigned int *iu, unsigned int *id, const unsigned int i, const unsigned short mu, const unsigned short nu, const unsigned short leaf)
 Calculates the first half of the leaf for a clover term. We split it so that the force term can reuse the first half of the leaf.
template<typename T>
__device__ void Device::Half_Leaf (complex< T > Leaves[nc], complex< T > *u11t, complex< T > *u12t, complex< T > a[nc], unsigned int *iu, unsigned int *id, const unsigned int i, const unsigned short mu, const unsigned short nu, const unsigned short leaf)
 Calculates the first half of the leaf for a clover term. We split it so that the force term can reuse the first half of the leaf.
template<typename T>
__device__ void Device::Leaf (complex< T > *u11t, complex< T > *u12t, complex< T > Leaves[nc], unsigned int *iu, unsigned int *id, unsigned int i, const unsigned short mu, const unsigned short nu, const unsigned short leaf)
 Calculates a leaf for a clover term.
template<typename T>
__global__ void Kernels::Half_Leaves (complex< T > *hLeaves0, complex< T > *hLeaves1, complex< T > *u11t, complex< T > *u12t, unsigned int *iu, unsigned int *id, const unsigned short mu, const unsigned short nu)
 Calculates the products of the first two links in a plaquette.
template<typename T>
__global__ void Kernels::Full_Clover (complex< T > *clover1, complex< T > *clover2, complex< T > *u11t, complex< T > *u12t, unsigned int *iu, unsigned int *id, int mu, int nu)
 Calculates the clovers in all directions at all sites \( F_{\mu\nu}(n)=\frac{-i}{8a^2}\left(Q_{\mu\nu}(n)-Q_{\nu\mu}(n)\right)\).

Detailed Description

Function Documentation

◆ ByGenLeft()

void ByGenLeft ( Complex_f a[nc],
const unsigned short gen )

Multiply leaf (or part of one) by generator from left.

The leaves contributing to each force term need to be scaled by the generator, but the generator appears at different points in each leaf. This routine multiples by the generator from the left side.

Parameters
[in,out]aThe leaf or partial leaf
[in]genWhat generator are we multiplying by?
Postcondition
Product stored in a

\(i\sigma_x\)

\(i\sigma_y\)

\(i\sigma_z\)

Definition at line 9 of file clover.c.

9 {
10 Complex_f tmp = a[0];
11 switch(gen){
13 case(0):
14 a[0] = -cimagf(a[1])-crealf(a[1])*I;
15 a[1] = cimagf(tmp)+crealf(tmp)*I;
16 break;
18 case(1):
19 a[0] = -conjf(a[1]);
20 a[1] = conjf(tmp);
21 break;
23 case(2):
24 a[0] = -cimagf(a[0])+crealf(a[0])*I;
25 a[1] = -cimagf(a[1])+crealf(a[1])*I;
26 break;
27 }
28 return;
29}
Complex_f
#define Complex_f
Single precision complex number.
Definition sizes.h:62
I
#define I
Define I in double precision using C standard notation.
Definition thrust_complex.h:28

References Complex_f, I, and nc.

Here is the caller graph for this function:

◆ ByGenRight()

void ByGenRight ( Complex_f a[nc],
const unsigned short gen )

Multiply leaf (or part of one) by generator from right.

The leaves contributing to each force term need to be scaled by the generator, but the generator appears at different points in each leaf. This routine multiples by the generator from the right side.

Parameters
[in,out]aThe leaf or partial leaf
[in]genWhat generator are we multiplying by?
Postcondition
Product stored in a

\(i\sigma_x\)

\(i\sigma_y\)

\(i\sigma_z\)

Definition at line 31 of file clover.c.

31 {
32 Complex_f tmp = a[0];
33 switch(gen){
35 case(0):
36 a[0] = -cimagf(a[1])+crealf(a[1])*I;
37 a[1] = -cimagf(tmp)+ crealf(tmp)*I;
38 break;
40 case(1):
41 a[0]=-a[1]; a[1]=tmp;
42 break;
44 case(2):
45 a[0] = -cimagf(a[0])+crealf(a[0])*I;
46 a[1] = cimagf(a[1])-crealf(a[1])*I;
47 break;
48 }
49 return;
50}

References Complex_f, I, and nc.

◆ Clover()

void Clover ( Complex_f * clover[2],
Complex_f * ut[2],
unsigned int * iu,
unsigned int * id )

Calculates the clovers in all directions at all sites.

\[ F_{\mu\nu}(n)=\frac{-i}{8a^2}\left(Q_{\mu\nu}(n)-Q_{\nu\mu}(n)\right)\]

Parameters
[out]cloverArray of clovers
[in]utGauge fields
[in]iu,idUpper and lower indices
Postcondition
Clover stored in clover

Definition at line 203 of file clover.c.

203 {
204 const char funcname[]="Full_Clover";
205#ifdef USE_GPU
206 cuClover(clover,ut,iu,id);
207#else
208 clover[0]=aligned_alloc(AVX,6*kvol*sizeof(Complex_f));
209 clover[1]=aligned_alloc(AVX,6*kvol*sizeof(Complex_f));
210 for(unsigned short mu=0;mu<ndim-1;mu++)
211 for(unsigned short nu=mu+1;nu<ndim;nu++)
212 if(mu!=nu){
213 //Clover index
214 unsigned short clov = (mu==0) ? nu-1 :mu+nu;
215#pragma omp parallel for
216 for(unsigned int i=0;i<kvol;i++){
217 clover[0][i+clov*kvol]=0;
218 clover[1][i+clov*kvol]=0;
219 Complex_f Leaves[nc];
220 for(unsigned short leaf=0;leaf<ndim;leaf++)
221 {
222 Leaf(Leaves,ut,iu,id,i,mu,nu,leaf);
223 clover[0][i+clov*kvol]+=Leaves[0]; clover[1][i+clov*kvol]+=Leaves[1];
224 }
227
232 clover[0][i+clov*kvol]=cimagf(clover[0][i+clov*kvol]); clover[0][i+clov*kvol]*=(1.0f/4.0f);
234 clover[1][i+clov*kvol]+=clover[1][i+clov*kvol]; clover[1][i+clov*kvol]*=(-I/8.0f);
235 }
236 }
237#endif
238 return;
239}
cuClover
int cuClover(Complex_f *clover[nc], Complex_f *ut[nc], unsigned int *iu, unsigned int *id)
CUDA wrapper for calculating the clovers in all directions at all sites .
Definition cuclover.cu:694
Leaf
void Leaf(Complex_f Leaves[nc], Complex_f *ut[nc], unsigned int *iu, unsigned int *id, unsigned int i, const unsigned short mu, const unsigned short nu, const unsigned short leaf)
Calculates a leaf for a clover term.
Definition clover.c:135
AVX
#define AVX
Alignment of arrays. 64 for AVX-512, 32 for AVX/AVX2. 16 for SSE. Since AVX is standard on modern x86...
Definition sizes.h:279
nc
#define nc
Colours.
Definition sizes.h:182
kvol
#define kvol
Sublattice volume.
Definition sizes.h:163
ndim
#define ndim
Dimensions.
Definition sizes.h:188

References AVX, Complex_f, cuClover(), I, kvol, Leaf(), nc, and ndim.

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

◆ Clover_free()

void Clover_free ( Complex_f * clover[nc])
inline

Free's memory used for clover terms and leaves.

Parameters
[in,out]cloverClovers
Postcondition
clover memory freed

Definition at line 750 of file clover.c.

750 {
751 for(unsigned short c=0;c<nc;c++){
752#ifdef USE_GPU
753#ifdef _DEBUG
754 cudaFree(clover[c]);
755#else
756 cudaFreeAsync(clover[c],streams[c]);
757#endif
758#else
759 free(clover[c]);
760#endif
761 }
762}
streams
cudaStream_t streams[ndirac *ndim *nadj]
An array of concurrent GPU streams to keep it busy.
Definition cusu2hmc.cu:29

References Complex_f, nc, and streams.

Here is the caller graph for this function:

◆ Clover_SU2plaq()

void Clover_SU2plaq ( Complex_f * ut[nc],
Complex_f Leaves[nc],
unsigned int * iu,
int i,
int mu,
int nu )

Calculates the SU2 plaquette at site i in the \(\mu--\nu\) direction.

Parameters
[in]utTrial fields
[out]LeavesTrial fields
[in]iuUpper halo indices
[in]isite index
[in]mu,nuPlaquette direction. Note that mu and nu can be negative to facilitate calculating plaquettes for Clover terms. No sanity checks are conducted on them in this routine.
Postcondition
Leaves overwritten by plaquette values

References Complex, Complex_f, and nc.

◆ cuClover()

int cuClover ( Complex_f * clover[nc],
Complex_f * ut[nc],
unsigned int * iu,
unsigned int * id )

CUDA wrapper for calculating the clovers in all directions at all sites \( F_{\mu\nu}(n)=\frac{-i}{8a^2}\left(Q_{\mu\nu}(n)-Q_{\nu\mu}(n)\right)\).

Parameters
[out]cloverArray of clovers
[in]utGauge fields
[in]iu,idUpper and lower indices
Returns
Zero on success, integer error code otherwise
Postcondition
Clover stored in clover

Definition at line 694 of file cuclover.cu.

694 {
695 const char funcname[]="cuClover";
696#ifdef _DEBUG
697 cudaMallocManaged((void **)&clover[0],6*kvol*sizeof(Complex_f),cudaMemAttachGlobal);
698 cudaMallocManaged((void **)&clover[1],6*kvol*sizeof(Complex_f),cudaMemAttachGlobal);
699#else
700 cudaMallocAsync((void **)&clover[0],6*kvol*sizeof(Complex_f),streams[0]);
701 cudaMallocAsync((void **)&clover[1],6*kvol*sizeof(Complex_f),streams[1]);
702#endif
703 for(unsigned short mu=0;mu<ndim-1;mu++)
704 for(unsigned short nu=mu+1;nu<ndim;nu++)
705 if(mu!=nu){
706 //Clover index
707 unsigned short clov = (mu==0) ? nu-1 :mu+nu;
708 //Allocate clover memory
709 //Note that the clover is completely local, so doesn't need a halo for MPI
710 Full_Clover<<<dimGrid,dimBlock,0,streams[clov]>>>(clover[0]+clov*kvol,clover[1]+clov*kvol,\
711 ut[0],ut[1],iu,id,mu,nu);
712 }
713 cudaDeviceSynchronise();
714 return 0;
715}
Kernels::Full_Clover
__global__ void Full_Clover(complex< T > *clover1, complex< T > *clover2, complex< T > *u11t, complex< T > *u12t, unsigned int *iu, unsigned int *id, int mu, int nu)
Calculates the clovers in all directions at all sites .
Definition cuclover.cu:323
cudaDeviceSynchronise
#define cudaDeviceSynchronise()
Get rid of that bastardised yankee English.
Definition sizes.h:53

References Complex_f, cudaDeviceSynchronise, Kernels::Full_Clover(), kvol, nc, ndim, and streams.

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

◆ Full_Clover()

template<typename T>
__global__ void Kernels::Full_Clover ( complex< T > * clover1,
complex< T > * clover2,
complex< T > * u11t,
complex< T > * u12t,
unsigned int * iu,
unsigned int * id,
int mu,
int nu )

Calculates the clovers in all directions at all sites \( F_{\mu\nu}(n)=\frac{-i}{8a^2}\left(Q_{\mu\nu}(n)-Q_{\nu\mu}(n)\right)\).

Parameters
[out]clover1,clover2Array of clovers
[in]u11t,u12tGauge fields
[in]iu,idUpper and lower indices
[in]mu,nuClover direction
Postcondition
Contents of clover1 and clover2 overwritten

The clover is given by \(F_{\mu\nu}=\frac{-i}{8}\left(Q_{\mu\nu}-Q_{\nu\mu}\right)\). We do that manually below.

The \(\alpha\) component. Only the imaginary part survives. And since it is multiplied by \(-i\) it is real. Need to be extra cautious here though .imag() returns a real value. So we multiply by I_f manually The 8.0f becomes a 4.0f to account for the factor of two

The \(\beta\) component. Both real and imaginary components survive. It ends up getting doubled.

Definition at line 323 of file cuclover.cu.

324 {
325 const volatile int gsize = gridDim.x*gridDim.y*gridDim.z;
326 const volatile int bsize = blockDim.x*blockDim.y*blockDim.z;
327 const volatile int blockId = blockIdx.x+ blockIdx.y * gridDim.x+ gridDim.x * gridDim.y * blockIdx.z;
328 const volatile int bthreadId= (threadIdx.z * blockDim.y+ threadIdx.y)* blockDim.x+ threadIdx.x;
329 const volatile int gthreadId= blockId * bsize+bthreadId;
330
331 complex<T> Leaves[2];
332 for(unsigned int i=gthreadId;i<kvol;i+=gsize*bsize){
333 clover1[i]=0;clover2[i]=0;
334 for(unsigned short leaf=0;leaf<ndim;leaf++)
335 {
336 //Pointer arithemetic on the leaves.
337 Leaf(u11t,u12t,Leaves,iu,id,i,mu,nu,leaf);
338 clover1[i]+=Leaves[0]; clover2[i]+=Leaves[1];
339 }
342
346 clover1[i]=clover1[i].imag(); clover1[i]*=(1.0f/4.0f);
347 // clover1[i]=clover1[1].imag()/4.0f;
348
350 clover2[i]+=clover2[i]; clover2[i]*=(-I_f/8.0f);
351 }
352 return;
353 }
Device::Leaf
__device__ void Leaf(complex< T > *u11t, complex< T > *u12t, complex< T > Leaves[nc], unsigned int *iu, unsigned int *id, unsigned int i, const unsigned short mu, const unsigned short nu, const unsigned short leaf)
Calculates a leaf for a clover term.
Definition cuclover.cu:153
I_f
#define I_f
Define I in single precision.
Definition thrust_complex.h:30

References I_f, kvol, Device::Leaf(), and ndim.

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

◆ Half_Leaf() [1/2]

template<typename T>
__device__ void Device::Half_Leaf ( complex< T > Leaves[nc],
complex< T > * u11t,
complex< T > * u12t,
complex< T > a[nc],
unsigned int * iu,
unsigned int * id,
const unsigned int i,
const unsigned short mu,
const unsigned short nu,
const unsigned short leaf )

Calculates the first half of the leaf for a clover term. We split it so that the force term can reuse the first half of the leaf.

Parameters
[out]LeavesLeaf
[in]u11t,u12tGauge fields
[in]aBuffer array
[in]iu,idUpper and lower site indices
[in]iLattice index of the clover in question
[in]mu,nuDirection in which we're evaluating the leaf
[in]leafWhich leaf of the clover is being calculated

Both positive is just a standard plaquette

\(U_\mu(x)U_\nu(x+\hat{\mu})\)

Leaf in the forward nu and backwards mu direction

\(U_\nu(x)U^\dagger_\mu(x-\hat{\mu}+\hat{\nu})\)

Leaf in the backwards nu and forwards mu direction

\(U^\dagger_\nu(x-\hat{\nu})U_\mu(x-\hat{\nu})\)

Leaf in the backwards mu and backwards nu direction

\(U_\mu^\dagger(x-\hat{\mu})U_\nu^\dagger(x-\hat{\mu}-\hat{\nu})\)

Definition at line 91 of file cuclover.cu.

92 {
93
94 unsigned int ind; unsigned int double_ind=id[nu*kvol+ind];
95 switch(leaf){
96 case(0):
98 a[0]=u11t[i+kvolHalo*mu]; a[1]=u12t[i+kvolHalo*mu];
99 ind = iu[mu*kvol+i];
100
102 Leaves[0]=a[0]*u11t[ind+kvolHalo*nu]-a[1]*conj(u12t[ind+kvolHalo*nu]);
103 Leaves[1]=a[0]*u12t[ind+kvolHalo*nu]+a[1]*conj(u11t[ind+kvolHalo*nu]);
104 break;
105 case(1):
107 ind = id[mu*kvol+i];
108 a[0]=u11t[i+kvolHalo*nu]; a[1]=u12t[i+kvolHalo*nu];
109 //Awkward index...
110 double_ind=iu[nu*kvol+ind];
112 Leaves[0]=a[0]*conj(u11t[double_ind+kvolHalo*mu])+a[1]*conj(u12t[double_ind+kvolHalo*mu]);
113 Leaves[1]=-a[0]*u12t[double_ind+kvolHalo*mu]+a[1]*u11t[double_ind+kvolHalo*mu];
114 break;
115 case(2):
117 ind = id[nu*kvol+i];
118 //Daggered. So Conj what goes into a[0] and negate what goes into a[1]
119 a[0]=conj(u11t[ind+kvolHalo*nu]); a[1]=-u12t[ind+kvolHalo*nu];
120
122 Leaves[0]=a[0]*u11t[ind+kvolHalo*mu]-a[1]*conj(u12t[ind+kvolHalo*mu]);
123 //Don't forget negatiion of second term was handled earlier!
124 Leaves[1]=a[0]*u12t[ind+kvolHalo*mu]+a[1]*conj(u11t[ind+kvolHalo*mu]);
125 break;
126 case(3):
128 ind = id[i+kvol*mu];
129 //Daggered. So Conj what goes into a[0] and negate what goes into a[1]
130 a[0]=conj(u11t[ind+kvolHalo*mu]); a[1]=-u12t[ind+kvolHalo*mu];
131 //Another awkward index
132 double_ind=id[nu*kvol+ind];
133
135 Leaves[0]=a[0]*conj(u11t[double_ind+kvolHalo*nu])+a[1]*conj(u12t[double_ind+kvolHalo*nu]);
136 Leaves[1]=-a[0]*u12t[double_ind+kvolHalo*nu]+a[1]*u11t[double_ind+kvolHalo*nu];
137 break;
138 }
139 return;
140 }
conj
__device__ __forceinline__ T conj(const T &z)
Complex Conjugation.
Definition cusu2hmc.cu:33
kvolHalo
#define kvolHalo
Subvolume + halo size.
Definition sizes.h:234

References conj(), kvol, kvolHalo, and nc.

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

◆ Half_Leaf() [2/2]

void Half_Leaf ( Complex_f Leaves[nc],
Complex_f * ut[nc],
Complex_f a[nc],
unsigned int * iu,
unsigned int * id,
const unsigned int i,
const unsigned short mu,
const unsigned short nu,
const unsigned short leaf )

Calculates the first half of the leaf for a clover term. We split it so that the force term can reuse the first half of the leaf.

Parameters
[out]LeavesLeaf
[in]utGauge fields
[in]aBuffer array
[in]iu,idUpper and lower site indices
[in]iLattice index of the clover in question
[in]mu,nuDirection in which we're evaluating the leaf
[in]leafWhich leaf of the clover is being calculated
Postcondition
Leaves replaced with output

Both positive is just a standard plaquette

\(U_\mu(x)U_\nu(x+\hat{\mu})\)

Leaf in the forward nu and backwards mu direction

\(U_\nu(x)U^\dagger_\mu(x-\hat{\mu}+\hat{\nu})\)

Leaf in the backwards nu and forwards mu direction

\(U^\dagger_\nu(x-\hat{\nu})U_\mu(x-\hat{\nu})\)

Leaf in the backwards mu and backwards nu direction

\(U_\mu^\dagger(x-\hat{\mu})U_\nu^\dagger(x-\hat{\mu}-\hat{\nu})\)

Definition at line 70 of file clover.c.

71 {
72 unsigned int uidm;
73 switch(leaf){
74 case(0):
76 a[0]=ut[0][i+kvolHalo*mu]; a[1]=ut[1][i+kvolHalo*mu];
77 uidm = iu[mu*kvol+i];
78
80 Leaves[0]=a[0]*ut[0][uidm+kvolHalo*nu]-a[1]*conjf(ut[1][uidm+kvolHalo*nu]);
81 Leaves[1]=a[0]*ut[1][uidm+kvolHalo*nu]+a[1]*conjf(ut[0][uidm+kvolHalo*nu]);
82 break;
83 case(1):
85 //Should really read didm, but I've already declared this
86 uidm = id[mu*kvol+i];
87 a[0]=ut[0][i+kvolHalo*nu]; a[1]=ut[1][i+kvolHalo*nu];
88 //Awkward index...
89 const unsigned int uin_didm=iu[nu*kvol+uidm];
91 Leaves[0]=a[0]*conjf(ut[0][uin_didm+kvolHalo*mu])+a[1]*conjf(ut[1][uin_didm+kvolHalo*mu]);
92 Leaves[1]=-a[0]*ut[1][uin_didm+kvolHalo*mu]+a[1]*ut[0][uin_didm+kvolHalo*mu];
93 break;
94 case(2):
96 //Should really read didn, but I've already declared this
97 uidm = id[nu*kvol+i];
98 //Daggered. So Conj what goes into a[0] and negate what goes into a[1]
99 a[0]=conjf(ut[0][uidm+kvolHalo*nu]); a[1]=-ut[1][uidm+kvolHalo*nu];
100
102 //Don't forget negation of second term was handled earlier!
103 Leaves[0]=a[0]*ut[0][uidm+kvolHalo*mu]-a[1]*conjf(ut[1][uidm+kvolHalo*mu]);
104 Leaves[1]=a[0]*ut[1][uidm+kvolHalo*mu]+a[1]*conjf(ut[0][uidm+kvolHalo*mu]);
105 break;
106 case(3):
108 //Should really read didm, but I've already declared this
109 uidm = id[i+kvol*mu];
110 //Daggered. So Conj what goes into a[0] and negate what goes into a[1]
111 a[0]=conjf(ut[0][uidm+kvolHalo*mu]); a[1]=-ut[1][uidm+kvolHalo*mu];
112 //Another awkward index
113 const unsigned int din_didm=id[nu*kvol+uidm];
114
116 Leaves[0]=a[0]*conjf(ut[0][din_didm+kvolHalo*nu])+a[1]*conjf(ut[1][din_didm+kvolHalo*nu]);
117 Leaves[1]=-a[0]*ut[1][din_didm+kvolHalo*nu]+a[1]*ut[0][din_didm+kvolHalo*nu];
118 break;
119 }
120 return;
121}

References Complex_f, kvol, kvolHalo, and nc.

Here is the caller graph for this function:

◆ Half_Leaves() [1/2]

void Half_Leaves ( Complex_f * hLeaves[2],
Complex_f * ut[2],
unsigned int * iu,
unsigned int * id,
const unsigned short mu,
const unsigned short nu )

Calculates the products of the first two links in a plaquette.

Parameters
[out]hLeavesProduct of first two links
[in]utGauge fields
[in]iu,idUpper and lower indices
[in]mu,nuClover direction
Postcondition
Product of first two links stored in hLeaves

Definition at line 122 of file clover.c.

123 {
124
125#pragma omp parallel for simd collapse(2)
126 for(unsigned short leaf=0;leaf<ndim;leaf++)
127 for(unsigned int i=0;i<kvol;i++){
128 Complex_f Leaves[nc], a[nc];
129 Half_Leaf(Leaves,ut,a,iu,id,i,mu,nu,leaf);
130 hLeaves[0][i+kvol*leaf]=Leaves[0]; hLeaves[1][i+kvol*leaf]=Leaves[1];
131 }
132 return;
133}
Half_Leaf
void Half_Leaf(Complex_f Leaves[nc], Complex_f *ut[nc], Complex_f a[nc], unsigned int *iu, unsigned int *id, const unsigned int i, const unsigned short mu, const unsigned short nu, const unsigned short leaf)
Calculates the first half of the leaf for a clover term. We split it so that the force term can reuse...
Definition clover.c:70

References Complex_f, Half_Leaf(), kvol, nc, and ndim.

Here is the call graph for this function:

◆ Half_Leaves() [2/2]

template<typename T>
__global__ void Kernels::Half_Leaves ( complex< T > * hLeaves0,
complex< T > * hLeaves1,
complex< T > * u11t,
complex< T > * u12t,
unsigned int * iu,
unsigned int * id,
const unsigned short mu,
const unsigned short nu )

Calculates the products of the first two links in a plaquette.

Parameters
[out]hLeaves0,hLeaves1Product of first two links
[in]u11t,u12tGauge fields
[in]iu,idUpper and lower indices
[in]mu,nuClover direction
Postcondition
Contents of hLeaves0 and hLeaves1 overwritten

Definition at line 294 of file cuclover.cu.

295 {
296 const volatile int gsize = gridDim.x*gridDim.y*gridDim.z;
297 const volatile int bsize = blockDim.x*blockDim.y*blockDim.z;
298 const volatile int blockId = blockIdx.x+ blockIdx.y * gridDim.x+ gridDim.x * gridDim.y * blockIdx.z;
299 const volatile int bthreadId= (threadIdx.z * blockDim.y+ threadIdx.y)* blockDim.x+ threadIdx.x;
300 const volatile int gthreadId= blockId * bsize+bthreadId;
301
302 complex<T> Leaves[nc], a[nc];
303 for(unsigned short leaf=0;leaf<ndim;leaf++)
304 for(unsigned int i=gthreadId;i<kvol;i+=gsize*bsize){
305 Half_Leaf(Leaves,u11t,u12t,a,iu,id,i,mu,nu,leaf);
306 hLeaves0[i+kvol*leaf]=Leaves[0]; hLeaves1[i+kvol*leaf]=Leaves[1];
307 }
308 return;
309 }
Device::Half_Leaf
__device__ void Half_Leaf(complex< T > Leaves[nc], complex< T > *u11t, complex< T > *u12t, complex< T > a[nc], unsigned int *iu, unsigned int *id, const unsigned int i, const unsigned short mu, const unsigned short nu, const unsigned short leaf)
Calculates the first half of the leaf for a clover term. We split it so that the force term can reuse...
Definition cuclover.cu:91

References Device::Half_Leaf(), kvol, nc, and ndim.

Here is the call graph for this function:

◆ Init_clover()

int Init_clover ( Complex ** sigval,
Complex_f ** sigval_f,
unsigned short ** sigin,
float c_sw )

Initialise values needed for the clover terms.

Parameters
[out]sigval,sigval_f\( \sigma_{\mu\nu}=\frac{1}{2i}[\gamma_\mu,\gamma_\nu]\) in double and single precision scaled by \(c_{sw}\)
[out]siginWhich column does row idirac of \(\sigma_{\mu\nu}\) act on
[in]c_swClover coefficient
Returns
Zero on success, integer error code otherwise
Postcondition
sigval and sigval_f initialised with matrix entries. sigin initialised with index of non-zero entries

Definition at line 705 of file clover.c.

705 {
706 const char funcname[] = "Init_clover";
707 unsigned short __attribute__((aligned(AVX))) sigin_t[6][4] = {{0,1,2,3},{1,0,3,2},{1,0,3,2},{1,0,3,2},{1,0,3,2},{0,1,2,3}};
708 //The sigma matrices are the commutators of the gamma matrices. These are antisymmetric when you swap the indices
709 //0 is sigma_0,1
710 //1 is sigma_0,2
711 //2 is sigma_0,3
712 //3 is sigma_1,2
713 //4 is sigma_1,3
714 //5 is sigma_2,3
715 Complex __attribute__((aligned(AVX))) sigval_t[6][4] = {{1,-1,1,-1},{I,-I,I,-I},{-1,-1,1,1},{1,1,1,1},{I,-I,-I,I},{-1,1,1,-1}};
716 //Complex __attribute__((aligned(AVX))) sigval_t[6][4] = {{1,1,1,1},{1,1,1,1},{1,1,1,1},{1,1,1,1},{1,1,1,1},{1,1,1,1},{1,1,1,1}};
717 //We mutiply by 1/2 and c_sw here since sigval is never used without them.
718#if defined USE_BLAS
719 cblas_zdscal(6*4, 0.5*c_sw, sigval_t, 1);
720#else
721#pragma omp parallel for simd collapse(2) aligned(sigval,sigval_f:AVX)
722 for(int i=0;i<6;i++)
723 for(int j=0;j<4;j++)
724 sigval_t[i][j]*=c_sw*0.5;
725#endif
726
727#ifdef USE_GPU
728 int device = -1;
729 cudaGetDevice(&device);
730
731 cudaMalloc((void **)sigin,6*4*sizeof(short));
732 cudaMalloc((void **)sigval,6*4*sizeof(Complex));
733 cudaMalloc((void **)sigval_f,6*4*sizeof(Complex_f));
734
735 cudaMemcpy(*sigin,sigin_t,6*4*sizeof(short),cudaMemcpyDefault);
736 cudaMemcpy(*sigval,sigval_t,6*4*sizeof(Complex),cudaMemcpyDefault);
737
738 cuComplex_convert(*sigval_f,*sigval,24,true,dimBlockOne,dimGridOne);
739#else
740 *sigin = (unsigned short *)malloc(6*4*sizeof(short));
741 *sigval=(Complex *)malloc(6*4*sizeof(Complex));
742 *sigval_f=(Complex_f *)malloc(6*4*sizeof(Complex_f));;
743 memcpy(*sigval,sigval_t,6*4*sizeof(Complex));
744 memcpy(*sigin,sigin_t,6*4*sizeof(short));
745 for(int i=0;i<6*4;i++)
746 *(*sigval_f+i)=(Complex_f)*(*sigval+i);
747#endif
748 return 0;
749}
cuComplex_convert
void cuComplex_convert(Complex_f *a, Complex *b, const unsigned int len, const bool dtof, dim3 dimBlock, dim3 dimGrid)
takes an array of complex float and double precision numbers and converts the precision
Definition cusu2hmc.cu:284
Complex
#define Complex
Double precision complex number.
Definition sizes.h:64
dimBlockOne
dim3 dimBlockOne
block size of one
Definition cusu2hmc.cu:22
dimGridOne
dim3 dimGridOne
Grid size of one.
Definition cusu2hmc.cu:23

References AVX, Complex, Complex_f, cuComplex_convert(), dimBlockOne, dimGridOne, and I.

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

◆ Leaf() [1/2]

template<typename T>
__device__ void Device::Leaf ( complex< T > * u11t,
complex< T > * u12t,
complex< T > Leaves[nc],
unsigned int * iu,
unsigned int * id,
unsigned int i,
const unsigned short mu,
const unsigned short nu,
const unsigned short leaf )

Calculates a leaf for a clover term.

Parameters
[out]LeavesArray of leaves
[in]u11t,u12tGauge fields
[in]iu,idUpper and lower site indices
[in]iLattice index of the clover in question
[in]mu,nuDirection in which we're evaluating the leaf
[in]leafWhich leaf of the clover is being calculated

\(U_\mu(x)U_\nu(x+\hat{\mu})U^\dagger_\mu(x+\hat{\nu})\)

\(U_\mu(x)U_\nu(x+\hat{\mu})U^\dagger_\mu(x+\hat{\nu})U^\dagger_\nu(x)\)

Leaf in the forwards nu and backwards nu direction

\(U_\nu(x)U^\dagger_\mu(x-\hat{\mu}+\hat{\nu})U^\dagger_\nu(x-\hat{\mu})\)

\(U_\nu(x)U^\dagger_\mu(x-\hat{\mu}+\hat{\nu})U^\dagger_\nu(x-\hat{\mu})U_\mu(x-\hat{\mu})\)

Leaf in the forwards mu and backwards nu direction

\(U^\dagger_\nu(x-\hat{\nu})U_\mu(x-\hat{\nu})U_\nu(x-\hat{\nu}+\hat{\mu})\)

\(U^\dagger_\nu(x-\hat{\nu})U_\mu(x-\hat{\nu})U_\nu(x-\hat{\nu}+\hat{\mu})U^\dagger_\mu(x)\)

Leaf in the backwards mu and backwards nu direction

\(U_\mu^\dagger(x-\hat{\mu})U_\nu^\dagger(x-\hat{\mu}-\hat{\nu})U_\mu(n-\hat{\nu}-\hat{\mu})\)

\(U_\mu^\dagger(x-\hat{\mu})U_\nu^\dagger(x-\hat{\mu}-\hat{\nu})U_\mu(n-\hat{\nu}-\hat{\mu})U_\nu(n-\hat{\nu})\)

Definition at line 153 of file cuclover.cu.

155 {
156 complex<T> a[nc];
157 Half_Leaf(Leaves,u11t,u12t,a,iu,id,i,mu,nu,leaf);
158 unsigned int ind,double_ind;
159 switch(leaf){
160 case(0)://Need braces for strict C++ standard compliance as a variable is declared in the case
161 {
162 ind = iu[nu*kvol+i];
164 a[0]=Leaves[0]*conj(u11t[ind+kvolHalo*mu])+Leaves[1]*conj(u12t[ind+kvolHalo*mu]);
165 a[1]=-Leaves[0]*u12t[ind+kvolHalo*mu]+Leaves[1]*u11t[ind+kvolHalo*mu];
166
168 Leaves[0]=a[0]*conj(u11t[i+kvolHalo*nu])+a[1]*conj(u12t[i+kvolHalo*nu]);
169 Leaves[1]=-a[0]*u12t[i+kvolHalo*nu]+a[1]*u11t[i+kvolHalo*nu];
170
171 }
172 break;
173 case(1):
175 ind = id[mu*kvol+i];
176
178 a[0]=Leaves[0]*conj(u11t[ind+kvolHalo*nu])+Leaves[1]*conj(u12t[ind+kvolHalo*nu]);
179 a[1]=-Leaves[0]*u12t[ind+kvolHalo*nu]+Leaves[1]*u11t[ind+kvolHalo*nu];
180
182 Leaves[0]=a[0]*u11t[ind+kvolHalo*mu]-a[1]*conj(u12t[ind+kvolHalo*mu]);
183 Leaves[1]=a[0]*u12t[ind+kvolHalo*mu]+a[1]*conj(u11t[ind+kvolHalo*mu]);
184 break;
185 case(2):
187 ind = id[nu*kvol+i];
188 double_ind=iu[mu*kvol+ind];
190 a[0]=Leaves[0]*u11t[double_ind+kvolHalo*nu]-Leaves[1]*conj(u12t[double_ind+kvolHalo*nu]);
191 a[1]=Leaves[0]*u12t[double_ind+kvolHalo*nu]+Leaves[1]*conj(u11t[double_ind+kvolHalo*nu]);
192
194 Leaves[0]=a[0]*conj(u11t[i+kvolHalo*mu])+a[1]*conj(u12t[i+kvolHalo*mu]);
195 Leaves[1]=-a[0]*u12t[i+kvolHalo*mu]+a[1]*u11t[i+kvolHalo*mu];
196
197 break;
198 case(3)://Need braces for strict C++ standard compliance as a variable is declared in the case
200 ind = id[nu*kvol+i];
201 double_ind=id[mu*kvol+ind];
202
204 a[0]=Leaves[0]*u11t[double_ind+kvolHalo*mu]-Leaves[1]*conj(u12t[double_ind+kvolHalo*mu]);
205 a[1]=Leaves[0]*u12t[double_ind+kvolHalo*mu]+Leaves[1]*conj(u11t[double_ind+kvolHalo*mu]);
206
208 Leaves[0]=a[0]*u11t[ind+kvolHalo*nu]-a[1]*conj(u12t[ind+kvolHalo*nu]);
209 Leaves[1]=a[0]*u12t[ind+kvolHalo*nu]+a[1]*conj(u11t[ind+kvolHalo*nu]);
210 break;
211 }
212 return;
213 }

References conj(), Half_Leaf(), kvol, kvolHalo, and nc.

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

◆ Leaf() [2/2]

void Leaf ( Complex_f Leaves[nc],
Complex_f * ut[nc],
unsigned int * iu,
unsigned int * id,
unsigned int i,
const unsigned short mu,
const unsigned short nu,
const unsigned short leaf )

Calculates a leaf for a clover term.

Parameters
[out]LeavesArray of leaves
[in]utGauge fields
[in]iu,idUpper and lower site indices
[in]iLattice index of the clover in question
[in]mu,nuDirection in which we're evaluating the leaf
[in]leafWhich leaf of the clover is being calculated
Postcondition
Clover leaf stored in Leaves

Plaquette

\(U_\mu(x)U_\nu(x+\hat{\mu})U^\dagger_\mu(x+\hat{\nu})\)

\(U_\mu(x)U_\nu(x+\hat{\mu})U^\dagger_\mu(x+\hat{\nu})U^\dagger_\nu(x)\)

Leaf in the forwards nu and backwards mu direction

\(U_\nu(x)U^\dagger_\mu(x-\hat{\mu}+\hat{\nu})U^\dagger_\nu(x-\hat{\mu})\)

\(U_\nu(x)U^\dagger_\mu(x-\hat{\mu}+\hat{\nu})U^\dagger_\nu(x-\hat{\mu})U_\mu(x-\hat{\mu})\)

Leaf in the forwards mu and backwards nu direction

\(U^\dagger_\nu(x-\hat{\nu})U_\mu(x-\hat{\nu})U_\nu(x-\hat{\nu}+\hat{\mu})\)

\(U^\dagger_\nu(x-\hat{\nu})U_\mu(x-\hat{\nu})U_\nu(x-\hat{\nu}+\hat{\mu})U^\dagger_\mu(x)\)

Leaf in the backwards mu and backwards nu direction

\(U_\mu^\dagger(x-\hat{\mu})U_\nu^\dagger(x-\hat{\mu}-\hat{\nu})U_\mu(n-\hat{\nu}-\hat{\mu})\)

\(U_\mu^\dagger(x-\hat{\mu})U_\nu^\dagger(x-\hat{\mu}-\hat{\nu})U_\mu(n-\hat{\nu}-\hat{\mu})U_\nu(n-\hat{\nu})\)

Definition at line 135 of file clover.c.

136 {
137 Complex_f a[nc];
138 Half_Leaf(Leaves,ut,a,iu,id,i,mu,nu,leaf);
139 unsigned int didm,didn,uidm;
140 switch(leaf){
141 case(0):
143 unsigned int uidn = iu[nu*kvol+i];
145 a[0]=Leaves[0]*conjf(ut[0][uidn+kvolHalo*mu])+Leaves[1]*conjf(ut[1][uidn+kvolHalo*mu]);
146 a[1]=-Leaves[0]*ut[1][uidn+kvolHalo*mu]+Leaves[1]*ut[0][uidn+kvolHalo*mu];
147
149 Leaves[0]=a[0]*conjf(ut[0][i+kvolHalo*nu])+a[1]*conjf(ut[1][i+kvolHalo*nu]);
150 Leaves[1]=-a[0]*ut[1][i+kvolHalo*nu]+a[1]*ut[0][i+kvolHalo*nu];
151
152 //DEBUG
153 // Leaves[0]=0; Leaves[1]=0;
154 break;
155 case(1):
157 didm = id[mu*kvol+i];
158
160 a[0]=Leaves[0]*conjf(ut[0][didm+kvolHalo*nu])+Leaves[1]*conjf(ut[1][didm+kvolHalo*nu]);
161 a[1]=-Leaves[0]*ut[1][didm+kvolHalo*nu]+Leaves[1]*ut[0][didm+kvolHalo*nu];
162
164 Leaves[0]=a[0]*ut[0][didm+kvolHalo*mu]-a[1]*conjf(ut[1][didm+kvolHalo*mu]);
165 Leaves[1]=a[0]*ut[1][didm+kvolHalo*mu]+a[1]*conjf(ut[0][didm+kvolHalo*mu]);
166 //DEBUG
167 // Leaves[0]=0; Leaves[1]=0;
168 break;
169 case(2):
171 didn = id[nu*kvol+i];
172 unsigned int uim_didn=iu[mu*kvol+didn];
174 a[0]=Leaves[0]*ut[0][uim_didn+kvolHalo*nu]-Leaves[1]*conjf(ut[1][uim_didn+kvolHalo*nu]);
175 a[1]=Leaves[0]*ut[1][uim_didn+kvolHalo*nu]+Leaves[1]*conjf(ut[0][uim_didn+kvolHalo*nu]);
176
178 Leaves[0]=a[0]*conjf(ut[0][i+kvolHalo*mu])+a[1]*conjf(ut[1][i+kvolHalo*mu]);
179 Leaves[1]=-a[0]*ut[1][i+kvolHalo*mu]+a[1]*ut[0][i+kvolHalo*mu];
180
181 //DEBUG
182 // Leaves[0]=0; Leaves[1]=0;
183 break;
184 case(3):
186 didn = id[nu*kvol+i];
187 unsigned int din_didm=id[mu*kvol+didn];
188
190 a[0]=Leaves[0]*ut[0][din_didm+kvolHalo*mu]-Leaves[1]*conjf(ut[1][din_didm+kvolHalo*mu]);
191 a[1]=Leaves[0]*ut[1][din_didm+kvolHalo*mu]+Leaves[1]*conjf(ut[0][din_didm+kvolHalo*mu]);
192
194 Leaves[0]=a[0]*ut[0][didn+kvolHalo*nu]-a[1]*conjf(ut[1][didn+kvolHalo*nu]);
195 Leaves[1]=a[0]*ut[1][didn+kvolHalo*nu]+a[1]*conjf(ut[0][didn+kvolHalo*nu]);
196
197 //DEBUG
198 // Leaves[0]=0; Leaves[1]=0;
199 break;
200 }
201 return;
202}

References Complex_f, Half_Leaf(), kvol, kvolHalo, and nc.

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