CHROMA
fuzwilp.cc
Go to the documentation of this file.
1 
2 // // version with added tmax (ACI)
3 /*! \file
4  * \brief Calculate ape-fuzzed Wilson loops
5  */
6 
7 #include "chromabase.h"
8 #include "meas/smear/ape_smear.h"
9 #include "meas/glue/fuzwilp.h"
10 
11 namespace Chroma {
12 //! Calculate ape-fuzzed Wilson loops
13 /*!
14  * \ingroup glue
15  *
16  * Computes time-like APE_fuzzed Wilson loops, including non-planar loops,
17  *
18  * This version makes APE-smeared links with no blocking as required
19  * for potential clculations
20  *
21  * Warning: This version is VERY Slow as it has non-recursive shifting
22  * of some link products
23  * Warning: this works only for Nc = 2 and 3 ! (Projection of
24  * smeared/blocked links)
25  *
26  * Warning: this version assumes the space-like directions (perpendicular
27  * to j_decay) to have equal length.
28  *
29  * \param u gauge field ( Read )
30  * \param j_decay 'time' direction for 'fuzzed' Wilson loops ( Read )
31  * \param tmax maximum time-extent loops ( Read )
32  * \param n_smear number of applying smearing to the gauge links ( Read )
33  * \param sm_fact "smearing" factor = weight of old link w. r. to staples ( Read )
34  * \param BlkAccu accuracy in fuzzy link projection ( Read )
35  * \param BlkMax maximum number of iterations in fuzzy link projection ( Read )
36  */
37 
38 void fuzwilp( const multi1d<LatticeColorMatrix>& u,
39  int j_decay, int tmax, int n_smear,
40  const Real& sm_fact, const Real& BlkAccu, int BlkMax,
41  XMLWriter& xml, const std::string& xml_group)
42 
43 {
44  START_CODE();
45 
46 
47  int lengthr;
48  int lengtht;
49  int lengthrs;
50 
51  multi1d<int> nrow(Nd);
52  nrow = Layout::lattSize();
53 
54 
55  lengtht = nrow[j_decay] / 2;
56  lengthr = nrow[0]/2;
57  lengthrs = (lengthr) * (lengthr+1) / 2;
58 
59  if (tmax < lengtht)
60  {
61  lengtht = tmax;
62  QDPIO::cout << " using lengtht = " << lengtht << std::endl;
63 
64  }
65  multi1d<LatticeColorMatrix> u_smear(Nd);
66  multi1d<LatticeColorMatrix> u_tmp(Nd);
67  multi2d<LatticeColorMatrix> u_prod( (Nd-1), lengthr);
68  LatticeColorMatrix tmp_tog;
69  LatticeColorMatrix up_t;
70  LatticeColorMatrix u_corn;
71  LatticeColorMatrix tmp_1;
72  LatticeColorMatrix tmp_2;
73  LatticeReal wl_trace;
74 
75  multi2d<Double> fuz_wlp1(lengthr, lengtht);
76  multi2d<Double> fuz_wlp2(lengthrs, lengtht);
77  Double ddummy;
78  Real rdummy;
79  Real vol;
80  int bl_level;
81  int mu;
82  int nu;
83  int mum;
84  int nun;
85  int t;
86  int r;
87  int s;
88  int i;
89  int n;
90 
91 /* Check that all direction perpendicular to j_decay are equal. */
92  for(mu = 1;mu < ( Nd); ++mu )
93  if( mu != j_decay && nrow[mu] != nrow[0] )
94  {
95  QDPIO::cout << " j_decay = " << j_decay << "and lengthr =" << nrow[0] << std::endl;
96 
97  QDP_error_exit("Wrong lattice size for Wilson loops: ", mu, nrow[mu]);
98  }
99 
100  /* First construct the smeared links for mu != j_decay */
101 
102  /* Copy u's to u_smear */
103  u_smear = u;
104 
105  bl_level = 0;
106 
107  /* Smear the space-like links n_smear times */
108  for(i = 1;i <= ( n_smear); ++i )
109  {
110  for(mu = 0;mu < ( Nd); ++mu )
111  {
112  if( mu != j_decay )
113  APE_Smear (u_smear, u_tmp[mu], mu, bl_level, sm_fact, BlkAccu, BlkMax, j_decay);
114  }
115  u_smear = u_tmp;
116  }
117 
118  // Compute the average link
119  Double link = 0.0;
120 
121  for(int mu=0; mu < Nd; ++mu)
122  link += sum(real(trace(u_smear[mu])));
123 
124  link /= double(Layout::vol()*Nd*Nc);
125  QDPIO::cout << "Average link after smearing: " << link << std::endl;
126 
127 
128 /* Construct products of smeared links (for mu != j_decay) */
129 
130 // u_prod[mum][r] r=0,1,2,....
131 // mum =mu direction with mu.ne.j_decay, i.e. 0,1,2
132 // ----> mum
133 // Us Us Us
134 // --- --- --- ....
135 //
136  mum = -1;
137  for(mu = 0;mu < ( Nd); ++mu )
138  if( mu != j_decay )
139  {
140  mum = mum + 1;
141  u_prod[mum][0] = u_smear[mu];
142  for(r=1;r < lengthr; ++r)
143  {
144  tmp_tog = shift(u_prod[mum][r-1],FORWARD, mu);
145  tmp_1 = u_smear[mu] * tmp_tog ;
146 //
147 // Up[r] = Us.[Us.Us.Us....Us]
148 // [Up[r-1]shifted +1 in mu
149 //
150  u_prod[mum][r] = tmp_1;
151  }
152  }
153 
154  fuz_wlp1 = 0; /* initialize the planar Wilson loops */
155  fuz_wlp2 = 0; /* initialize the nplanar Wilson loops */
156 
157  /* Compute 'fuzzied' Wilson loops (with normal links in 'time' direction!) */
158 
159 //
160 // For t x r planar case will use:
161 // U3
162 // --------->------
163 // | |
164 // | | ^
165 // U4 ^ ^ U2 t |
166 // | | (j_decay) |---->
167 // |-------->------| r (mu)
168 // U1
169 // so that
170 // W = Tr[ U1 * U2 * U3+ * U4+ ]
171 //
172 
173  for(t = 0;t < ( lengtht); ++t )
174  {
175 
176  /* Compute product of un-fuzzed links in t (j_decay) direction */
177  QDPIO::cout << "t= " << t << std::endl;
178  if( t == 0 )
179  {
180  up_t = u[j_decay];
181  }
182  else
183  {
184  tmp_tog = shift(up_t, FORWARD, j_decay);
185  up_t = u[j_decay] * tmp_tog;
186 //
187 // | U
188 // | U last Up shifted up one in t direction
189 // | U
190 // | U
191 // ----
192 // | U
193 // use new up_t for U4
194 //
195  }
196 
197  mum = -1;
198  for(mu = 0;mu < ( Nd); ++mu )
199  {
200  if( mu != j_decay )
201  {
202  mum = mum + 1;
203 // QDPIO::cout << "mu= " << mu << std::endl;
204  for(r = 0;r < ( lengthr); ++r )
205  {
206  /* Gather 'time-like' link segment from r-direction */
207 
208  tmp_tog = shift(up_t, FORWARD, mu);
209  n = r;
210  while( n > 0 )
211  {
212  tmp_1 = shift(tmp_tog, FORWARD, mu);
213  tmp_tog = tmp_1;
214  n = n - 1;
215  }
216  tmp_1 = tmp_tog;
217 //
218 // . |
219 // . ^ tmp_1 (Up shifted r+1 in mu direction)
220 // . |
221 // . |
222 // use tmp_1 for U2
223 
224  /* Gather 'space-like' link segment from t-direction */
225 
226  tmp_tog = shift(u_prod[mum][r], FORWARD, j_decay);
227  n = t;
228  while( n > 0 )
229  {
230  tmp_2 = shift( tmp_tog, FORWARD, j_decay);
231  tmp_tog = tmp_2;
232  n = n - 1;
233  }
234 //
235 //
236 // tmp_tog (u_prod shifted t+1 in t direction)
237 // -------->--------
238 //
239 //
240 //
241 // .................
242 // use tmp_tog for U3
243 //
244 //
245  /* Now complete the planar Wilson loop */
246 
247  tmp_2 = tmp_1 * adj(tmp_tog);
248 //
249 // U2 * U3+
250 //
251  tmp_1 = tmp_2 * adj(up_t);
252 //
253 // U2 * U3+ * U4+
254 //
255  tmp_2 = u_prod[mum][r] * tmp_1;
256 //
257 // U1 * U2 * U3+ * U4+
258 //
259  wl_trace = real(trace(tmp_2));
260  fuz_wlp1[r][t] += sum(wl_trace);
261 
262 // QDPIO::cout << "t= " << t << " r= " << r << std::endl;
263 // QDPIO::cout << "fuz_wlp1= " << 2*fuz_wlp1[r][t]/
264 // double((Layout::vol())*(Nd-1)*(Nd-2)*Nc) << std::endl;
265 
266  /* Now do non-planar loops */
267 
268 //
269 // We are currently inside t, mu, and r loops
270 // _____________
271 // / |
272 // / |
273 // / |
274 // | |
275 // | |
276 // ^ t (j_decay) |
277 // | ____________|
278 // | / s (nu)
279 // | /r (mu)
280 // |/
281 //
282  nun = -1;
283  for(nu = 0;nu < ( Nd); ++nu )
284  {
285  if ( nu != j_decay && nu == mu )
286  /* advance nun, since next "if" is not satisfied! */
287  nun = nun + 1;
288 
289  if ( nu != j_decay && nu != mu )
290  {
291  nun = nun + 1;
292 // QDPIO::cout << "nu= " << nu << std::endl;
293  for(s = 0;s <= ( r); ++s )
294  {
295  /* Construct the 'forward' space-like segments */
296 
297  /* Gather 'nu link segment' from r-direction */
298 
299  tmp_tog = shift(u_prod[nun][s], FORWARD, mu);
300 
301  n = r;
302  while( n > 0 )
303  {
304  tmp_1 = shift( tmp_tog, FORWARD, mu);
305  tmp_tog = tmp_1;
306  n = n - 1;
307  }
308  u_corn = u_prod[mum][r] * tmp_tog;
309 //
310 // ___________
311 // / s (nu)
312 // / r (mu) 1st u_corn contribution
313 // /
314 //
315 
316  /* Gather 'mu link segment' from s-direction */
317 
318  tmp_tog = shift(u_prod[mum][r], FORWARD, nu);
319  n = s;
320  while( n > 0 )
321  {
322  tmp_1 = shift(tmp_tog, FORWARD, nu);
323  tmp_tog = tmp_1;
324  n = n -1;
325  }
326  u_corn += u_prod[nun][s] * tmp_tog;
327 //
328 // /
329 // / r (mu) 2nd u_corn contribution
330 // ____________/
331 // s (nu)
332 
333  /* Now collect to construct the Wilson loop */
334 
335  /* Gather 'time-like' link segment from r-direction first */
336 
337  tmp_tog = shift(up_t, FORWARD, mu);
338  n = r;
339  while( n > 0 )
340  {
341  tmp_1 = shift(tmp_tog, FORWARD, mu);
342  tmp_tog = tmp_1;
343  n = n - 1;
344  }
345  tmp_1 = tmp_tog;
346 
347  /* Gather 'time-like' link segment from s-direction next */
348 
349  tmp_tog = shift(tmp_1, FORWARD, nu);
350  n = s;
351  while( n > 0 )
352  {
353  tmp_2 = shift(tmp_tog, FORWARD, nu);
354  tmp_tog = tmp_2;
355  n = n - 1;
356  }
357  tmp_1 = tmp_tog;
358 //
359 // . |
360 // . |
361 // . . |
362 // . . ^ tmp_1 is now Up shifted by s,t in mu,nu
363 // . . |
364 // . . |
365 // . . . . . |
366 // . .
367 // .
368 //
369 
370  /* Gather 'space-like' link segment from t-direction */
371 
372  tmp_tog = shift(u_corn, FORWARD, j_decay);
373  n = t;
374  while( n > 0 )
375  {
376  tmp_2 = shift(tmp_tog, FORWARD, j_decay);
377  tmp_tog = tmp_2;
378  n = n - 1;
379  }
380 //
381 // _____________
382 // / .
383 // / .
384 // / .
385 // . tmp_tog is u_corn shifted by t+1
386 // . .
387 // . .
388 // . .............
389 // . . s (nu)
390 // . . r (mu)
391 // ..
392 //
393  /* Now complete the 'forward' non-planar Wilson loop */
394 
395  tmp_2 = tmp_1 * adj(tmp_tog);
396  tmp_1 = tmp_2 * adj(up_t);
397  tmp_2 = u_corn * tmp_1;
398 
399 
400 // U3 (tmp_tog)
401 // _____________
402 // / |
403 // / |
404 // / |
405 // | | U2 (tmp_1)
406 // | |
407 // U4 (Up) | |
408 // | ____________|
409 // | /
410 // | / U1 (u_corn)
411 // |/
412 //
413 // so W = Tr [ U1 * U2 * U3+ * U4+ ]
414 //
415  wl_trace = real(trace(tmp_2));
416  n = r * (r+1) / 2 + s;
417  fuz_wlp2[n][t] += sum(wl_trace);
418 
419  /* Construct the 'backward' space-like segments */
420 
421  /* Gather 'nu link segment' from r-direction */
422 
423  tmp_tog = shift(u_prod[nun][s], FORWARD, mu);
424 
425  n = r;
426  while( n > 0 )
427  {
428  tmp_1 = shift(tmp_tog, FORWARD, mu);
429  tmp_tog = tmp_1;
430  n = n - 1;
431  }
432  tmp_1 = tmp_tog;
433 
434  /* Now fetch this from backward s-direction */
435 
436  tmp_tog = shift(tmp_1, BACKWARD, nu);
437  n = s;
438  while( n > 0 )
439  {
440  tmp_2 = shift(tmp_tog, BACKWARD, nu);
441  tmp_tog = tmp_2;
442  n = n - 1;
443  }
444  u_corn = u_prod[mum][r] * adj(tmp_tog);
445 
446  /* Gather 'mu link segment' from s-direction */
447 
448  tmp_tog = shift(u_prod[mum][r], BACKWARD, nu);
449  n = s;
450  while( n > 0 )
451  {
452  tmp_1 = shift(tmp_tog, BACKWARD, nu);
453  tmp_tog = tmp_1;
454  n = n - 1;
455  }
456  tmp_1 = tmp_tog;
457 
458  /* Gather 'nu link segment' from backward s-direction */
459 
460  tmp_tog = shift(u_prod[nun][s], BACKWARD, nu);
461  n = s;
462  while( n > 0 )
463  {
464  tmp_2 = shift(tmp_tog, BACKWARD, nu);
465  tmp_tog = tmp_2;
466  n = n - 1;
467  }
468  u_corn += adj(tmp_tog) * tmp_1;
469 
470 
471  /* Now collect to construct the Wilson loop */
472 
473  /* Gather 'time-like' link segment from r-direction first */
474 
475  tmp_tog = shift(up_t, FORWARD, mu);
476  n = r;
477  while( n > 0 )
478  {
479  tmp_1 = shift(tmp_tog, FORWARD, mu);
480  tmp_tog = tmp_1;
481  n = n - 1;
482  }
483  tmp_1 = tmp_tog;
484 
485  /* Gather 'time-like' link segment from backward s-direction next */
486 
487  tmp_tog = shift(tmp_1, BACKWARD, nu);
488  n = s;
489  while( n > 0 )
490  {
491  tmp_2 = shift(tmp_tog, BACKWARD, nu);
492  tmp_tog = tmp_2;
493  n = n - 1;
494  }
495  tmp_1 = tmp_tog;
496 
497  /* Gather 'space-like' link segment from t-direction */
498 
499  tmp_tog = shift(u_corn, FORWARD, j_decay);
500  n = t;
501  while( n > 0 )
502  {
503  tmp_2 = shift(tmp_tog, FORWARD, j_decay);
504  tmp_tog = tmp_2;
505  n = n - 1 ;
506  }
507 
508  /* Now complete the 'backward' non-planar Wilson loop */
509 
510  tmp_2 = tmp_1 * adj(tmp_tog);
511  tmp_1 = tmp_2 * adj(up_t);
512  tmp_2 = u_corn * tmp_1;
513  wl_trace = real(trace(tmp_2));
514  n = r * (r+1) / 2 + s;
515  fuz_wlp2[n][t] += sum(wl_trace);
516 
517  } /* end s loop */
518  } /* end nu != j_decay & nu != mu */
519  } /* end nu loop */
520  } /* end r loop */
521  } /* end mu != j_decay */
522  } /* end mu loop */
523  } /* end t loop */
524 
525  ddummy = 1.0 / double (Layout::vol()*Nc*(Nd-1)) ;
526 
527  push(xml,"fuz_wlp1"); // XML tag for fuz_wlp1
528  write(xml, "lengthr", lengthr);
529 
530  multi1d<Real> wloopr(lengtht);
531  for(r = 0; r < lengthr; ++r)
532  {
533  for(t = 0; t < lengtht; ++t)
534  {
535  fuz_wlp1[r][t] = fuz_wlp1[r][t] * ddummy;
536  wloopr[t] = fuz_wlp1[r][t];
537  }
538  write(xml, "r", r);
539  write(xml, "wloopr", wloopr); // write out fuz_wlp1
540  } // end for r
541  pop(xml); // XML end tag for fuz_wlp1
542 
543  ddummy = 1.0 / double(Layout::vol()*8*Nc*(Nd-1)*(Nd-2) ) ;
544 
545  for(t = 0;t < ( lengtht); ++t )
546  {
547  n = -1;
548  for(r = 0;r < ( lengthr); ++r )
549  for(s = 0;s <= ( r); ++s )
550  {
551  n = n + 1;
552  fuz_wlp2[n][t] = fuz_wlp2[n][t] * ddummy;
553  }
554  }
555  multi1d<Real> wlooprs(lengtht);
556  push(xml,"fuz_wlp2"); // XML tag for fuz_wlp2
557  write(xml, "lengthr", lengthr);
558 
559  n = -1;
560  for(r = 0; r < lengthr; ++r)
561  {
562  write(xml, "r", r);
563  for(s = 0;s <= ( r); ++s )
564  {
565  n = n + 1;
566  write(xml, "s", s);
567  for(t = 0; t < lengtht; ++t)
568  wlooprs[t] = fuz_wlp2[n][t];
569  write(xml, "wlooprs", wlooprs); // write out fuz_wlp2
570  } // end for s
571  } // end for r
572  pop(xml); // XML end tag for fuz_wlp2
573  QDPIO::cout << "fuz_wlp1 and fuz_wlp2 written to .xml file " << std::endl;
574 }
575 
576 }
Primary include file for CHROMA library code.
int mu
Definition: cool.cc:24
int nu
Definition: cool.cc:25
Calculate ape-fuzzed Wilson loops.
void write(XMLWriter &xml, const std::string &path, const AsqtadFermActParams &param)
Writer parameters.
void fuzwilp(const multi1d< LatticeColorMatrix > &u, int j_decay, int tmax, int n_smear, const Real &sm_fact, const Real &BlkAccu, int BlkMax, XMLWriter &xml, const std::string &xml_group)
Calculate ape-fuzzed Wilson loops.
Definition: fuzwilp.cc:38
void APE_Smear(const multi1d< LatticeColorMatrix > &u, LatticeColorMatrix &u_smear, int mu, int bl_level, const Real &sm_fact, const Real &BlkAccu, int BlkMax, int j_decay)
Construct APE smeared links from:
Definition: ape_smear.cc:44
unsigned n
Definition: ldumul_w.cc:36
LatticeColorMatrix tmp_2
Definition: meslate.cc:51
int j_decay
Definition: meslate.cc:22
int t
Definition: meslate.cc:37
LatticeColorMatrix tmp_1
Definition: meslate.cc:50
Nd
Definition: meslate.cc:74
Asqtad Staggered-Dirac operator.
Definition: klein_gord.cc:10
QDP_error_exit("too many BiCG iterations", n_count, rsd_sq, cp, c, re_rvr, im_rvr, re_a, im_a, re_b, im_b)
static multi1d< LatticeColorMatrix > u
push(xml_out,"Condensates")
int i
Definition: pbg5p_w.cc:55
pop(xml_out)
START_CODE()
multi1d< LatticeFermion > s(Ncb)
FloatingPoint< double > Double
Definition: gtest.h:7351
::std::string string
Definition: gtest.h:1979
Double link
Definition: pade_trln_w.cc:146
#define FORWARD
Definition: primitives.h:82
#define BACKWARD
Definition: primitives.h:83
Double sum
Definition: qtopcor.cc:37
Real rdummy
Definition: qtopcor.cc:38