|
int | Average_Plaquette (double *hg, double *avplaqs, double *avplaqt, Complex_f *ut[2], unsigned int *iu, float beta) |
| Calculates the gauge action using new (how new?) lookup table.
|
|
int | SU2plaq (Complex_f *ut[2], Complex_f Sigma[2], unsigned int *iu, int i, int mu, int nu) |
| Calculates the plaquette at site i in the \(\mu--\nu\) direction.
|
|
double | Polyakov (Complex_f *ut[2]) |
| Calculate the Polyakov loop (no prizes for guessing that one...)
|
|
Code for bosonic observables, Basically polyakov loop and Plaquette routines.
Code for bosonic observables.
Routines for polyakov loop, plaquettes and clovers.
- Author
- S J Hands (Original Fortran, March 2005)
-
P. Giudice (Hybrid Code, May 2013)
-
D. Lawlor (C version March 2021, CUDA/Mixed Precision/Clover Feb 2024 and beyond...)
Definition in file bosonic.c.
int Average_Plaquette |
( |
double * | hg, |
|
|
double * | avplaqs, |
|
|
double * | avplaqt, |
|
|
Complex_f * | ut[2], |
|
|
unsigned int * | iu, |
|
|
float | beta ) |
Calculates the gauge action using new (how new?) lookup table.
Follows a routine called qedplaq in some QED3 code
- Parameters
-
hg | Gauge component of Hamilton |
avplaqs | Average spacial Plaquette |
avplaqt | Average Temporal Plaquette |
ut | The trial fields |
iu | Upper halo indices |
beta | Inverse gauge coupling |
- See also
- Par_dsum
- Returns
- Zero on success, integer error code otherwise
Definition at line 20 of file bosonic.c.
20 {
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42 const char *funcname = "Average_Plaquette";
43
44
45
46
47
48
49#ifdef __NVCC__
50 __managed__ double hgs = 0; __managed__ double hgt = 0;
51 cuAverage_Plaquette(&hgs, &hgt, ut[0], ut[1], iu,dimGrid,dimBlock);
52#else
53 double hgs = 0; double hgt = 0;
54 for(
int mu=1;mu<
ndim;mu++)
55 for(int nu=0;nu<mu;nu++)
56
57
58#pragma omp parallel for simd reduction(+:hgs,hgt)
59 for(
int i=0;i<
kvol;i++){
62 switch(mu){
63
64 case(
ndim-1): hgt -= creal(Sigma[0]);
65 break;
66
67 default: hgs -= creal(Sigma[0]);
68 break;
69 }
70 }
71#endif
72#if(nproc>1)
74#endif
75 *avplaqs=-hgs/(3.0*
gvol); *avplaqt=-hgt/(
gvol*3.0);
76 *hg=(hgs+hgt)*beta;
77#ifdef _DEBUG
79 printf("hgs=%e hgt=%e hg=%e\n", hgs, hgt, *hg);
80#endif
81 return 0;
82}
int SU2plaq(Complex_f *ut[2], Complex_f Sigma[2], unsigned int *iu, int i, int mu, int nu)
Calculates the plaquette at site i in the direction.
int Par_dsum(double *dval)
Performs a reduction on a double dval to get a sum which is then distributed to all ranks.
#define kvol
Sublattice volume.
#define gvol
Lattice volume.
#define Complex_f
Single precision complex number.
References Complex_f, gvol, kvol, ndim, Par_dsum(), rank, and SU2plaq().
Calculate the Polyakov loop (no prizes for guessing that one...)
- Parameters
-
u11t,u12t | The gauge fields |
- See also
- Par_tmul, Par_dsum
- Returns
- Double corresponding to the polyakov loop
Calculate the Polyakov loop (no prizes for guessing that one...)
- Parameters
-
ut[0],ut[1] | The trial fields |
Calls:
Par_tmul(), Par_dsum()
- Returns
- Double corresponding to the polyakov loop
Definition at line 127 of file bosonic.c.
127 {
139 const char *funcname = "Polyakov";
140 double poly = 0;
142#ifdef __NVCC__
143 cuPolyakov(Sigma,ut,dimGrid,dimBlock);
144#else
147
148
149#ifdef USE_BLAS
150 cblas_ccopy(
kvol3, ut[0]+3,
ndim, Sigma[0], 1);
151 cblas_ccopy(
kvol3, ut[1]+3,
ndim, Sigma[1], 1);
152#else
153 for(
int i=0; i<
kvol3; i++){
154 Sigma[0][i]=ut[0][i*
ndim+3];
155 Sigma[1][i]=ut[1][i*
ndim+3];
156 }
157#endif
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172#pragma unroll
173 for(
int it=1;it<
ksizet;it++)
174#pragma omp parallel for simd
175 for(
int i=0;i<
kvol3;i++){
176
177
178 int indexu=it*
kvol3+i;
179 Complex_f a11=Sigma[0][i]*ut[0][indexu*
ndim+3]-Sigma[1][i]*conj(ut[1][indexu*
ndim+3]);
180
181 Sigma[1][i]=Sigma[0][i]*ut[1][indexu*
ndim+3]+Sigma[1][i]*conj(ut[0][indexu*
ndim+3]);
182 Sigma[0][i]=a11;
183 }
184 free(Sigma[1]);
185#endif
186
187
188
189
190
191
192#if (npt>1)
193#ifdef __NVCC_
194#error Par_tmul is not yet implimented in CUDA as Sigma[1] is device only memory
195#endif
196#ifdef _DEBUG
197 printf("Multiplying with MPI\n");
198#endif
199 Par_tmul(Sigma[0], Sigma[1]);
200
201#endif
202
203
204
205
206
207
208
210#pragma omp parallel for simd reduction(+:poly)
211 for(
int i=0;i<
kvol3;i++)
212 poly+=creal(Sigma[0][i]);
213#ifdef __NVCC__
214 cudaFree(Sigma[0]);
215#else
216 free(Sigma[0]);
217#endif
218
219#if(nproc>1)
221#endif
223 return poly;
224}
int * pcoord
The processor grid.
#define AVX
Alignment of arrays. 64 for AVX-512, 32 for AVX/AVX2. 16 for SSE. Since AVX is standard on modern x86...
#define gvol3
Lattice spatial volume.
#define ksizet
Sublattice t extent.
#define kvol3
Sublattice spatial volume.
References AVX, Complex_f, gvol3, ksizet, kvol3, ndim, Par_dsum(), pcoord, and rank.
int SU2plaq |
( |
Complex_f * | ut[2], |
|
|
Complex_f | Sigma[2], |
|
|
unsigned int * | iu, |
|
|
int | i, |
|
|
int | mu, |
|
|
int | nu ) |
|
inline |
Calculates the plaquette at site i in the \(\mu--\nu\) direction.
- Parameters
-
u2t | Trial fields |
Sigma | Plaquette components |
i | Lattice site |
iu | Upper halo indices |
mu,nu | Plaquette 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. |
- Returns
- double corresponding to the plaquette value
Calculates the trace of the plaquette at site i in the μ-ν direction
- Parameters
-
ut[0],ut[1] | Trial fields |
Sigma11,Sigma12 | Trial fields |
*iu | Upper halo indices |
i | site index |
mu,nu | Plaquette 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. |
Return:
double corresponding to the plaquette value
Definition at line 85 of file bosonic.c.
85 {
102 const char *funcname = "SU2plaq";
103 int uidm = iu[mu+
ndim*i];
104
105
106
107
108
109
110
111
112
113
114
115 Sigma[0]=ut[0][i*
ndim+mu]*ut[0][uidm*
ndim+nu]-ut[1][i*
ndim+mu]*conj(ut[1][uidm*
ndim+nu]);
116 Sigma[1]=ut[0][i*
ndim+mu]*ut[1][uidm*
ndim+nu]+ut[1][i*
ndim+mu]*conj(ut[0][uidm*
ndim+nu]);
117
118 int uidn = iu[nu+
ndim*i];
119 Complex_f a11=Sigma[0]*conj(ut[0][uidn*
ndim+mu])+Sigma[1]*conj(ut[1][uidn*
ndim+mu]);
121
122 Sigma[0]=a11*conj(ut[0][i*
ndim+nu])+a12*conj(ut[1][i*
ndim+nu]);
123 Sigma[1]=-a11*ut[1][i*
ndim+nu]+a12*ut[0][i*
ndim+mu];
124 return 0;
125}
References Complex_f, and ndim.