CHROMA
inline_prop_distillation_w.cc
Go to the documentation of this file.
1 /*! \file
2  * \brief Compute propagators from distillation
3  *
4  * Propagator calculation in distillation
5  */
6 
7 #include "fermact.h"
10 #include "meas/glue/mesplq.h"
11 #include "qdp_map_obj.h"
12 #include "qdp_map_obj_disk.h"
13 #include "qdp_map_obj_disk_multiple.h"
14 #include "qdp_disk_map_slice.h"
17 #include "util/ferm/transf.h"
18 #include "util/ferm/spin_rep.h"
19 #include "util/ferm/diractodr.h"
21 #include "util/ft/time_slice_set.h"
22 #include "util/info/proginfo.h"
26 
28 
29 namespace Chroma
30 {
31  //----------------------------------------------------------------------------
32  namespace InlinePropDistillationEnv
33  {
34  //! Propagator input
35  void read(XMLReader& xml, const std::string& path, Params::NamedObject_t& input)
36  {
37  XMLReader inputtop(xml, path);
38 
39  read(inputtop, "gauge_id", input.gauge_id);
40  read(inputtop, "colorvec_files", input.colorvec_files);
41  read(inputtop, "soln_file", input.soln_file);
42  }
43 
44  //! Propagator output
45  void write(XMLWriter& xml, const std::string& path, const Params::NamedObject_t& input)
46  {
47  push(xml, path);
48 
49  write(xml, "gauge_id", input.gauge_id);
50  write(xml, "colorvec_files", input.colorvec_files);
51  write(xml, "soln_file", input.soln_file);
52 
53  pop(xml);
54  }
55 
56 
57  //! Propagator input
58  void read(XMLReader& xml, const std::string& path, Params::Param_t::Contract_t& input)
59  {
60  XMLReader inputtop(xml, path);
61 
62  read(inputtop, "num_vecs", input.num_vecs);
63  read(inputtop, "t_sources", input.t_sources);
64  read(inputtop, "decay_dir", input.decay_dir);
65  read(inputtop, "Nt_forward", input.Nt_forward);
66  read(inputtop, "Nt_backward", input.Nt_backward);
67  read(inputtop, "mass_label", input.mass_label);
68  }
69 
70  //! Propagator output
71  void write(XMLWriter& xml, const std::string& path, const Params::Param_t::Contract_t& input)
72  {
73  push(xml, path);
74 
75  write(xml, "num_vecs", input.num_vecs);
76  write(xml, "t_sources", input.t_sources);
77  write(xml, "decay_dir", input.decay_dir);
78  write(xml, "Nt_forward", input.Nt_forward);
79  write(xml, "Nt_backward", input.Nt_backward);
80  write(xml, "mass_label", input.mass_label);
81 
82  pop(xml);
83  }
84 
85 
86  //! Propagator input
87  void read(XMLReader& xml, const std::string& path, Params::Param_t& input)
88  {
89  XMLReader inputtop(xml, path);
90 
91  read(inputtop, "Propagator", input.prop);
92  read(inputtop, "Contractions", input.contract);
93  }
94 
95  //! Propagator output
96  void write(XMLWriter& xml, const std::string& path, const Params::Param_t& input)
97  {
98  push(xml, path);
99 
100  write(xml, "Propagator", input.prop);
101  write(xml, "Contractions", input.contract);
102 
103  pop(xml);
104  }
105 
106 
107  //! Propagator input
108  void read(XMLReader& xml, const std::string& path, Params& input)
109  {
110  Params tmp(xml, path);
111  input = tmp;
112  }
113 
114  //! Propagator output
115  void write(XMLWriter& xml, const std::string& path, const Params& input)
116  {
117  push(xml, path);
118 
119  write(xml, "Param", input.param);
120  write(xml, "NamedObject", input.named_obj);
121 
122  pop(xml);
123  }
124  } // namespace InlinePropDistillationEnv
125 
126 
127  //----------------------------------------------------------------------------
128  //----------------------------------------------------------------------------
129  namespace InlinePropDistillationEnv
130  {
131 
132 #if 1
133  // Anonymous namespace
134  namespace
135  {
136  // Convenience type
137  //typedef QDP::MapObjectDisk<KeyTimeSliceColorVec_t, TimeSliceIO<LatticeColorVectorF> > MOD_t;
138  typedef QDP::MapObjectDiskMultiple<KeyTimeSliceColorVec_t, TimeSliceIO<LatticeColorVectorF> > MODS_t;
139 
140  //----------------------------------------------------------------------------
141  //! Get source key
142  KeyTimeSliceColorVec_t getSrcKey(int t_source, int colorvec_src)
143  {
145 
146  key.t_slice = t_source;
147  key.colorvec = colorvec_src;
148 
149  return key;
150  }
151 
152 
153  //----------------------------------------------------------------------------
154  //! Read a source std::vector
155  LatticeColorVector getSrc(MODS_t& source_obj, int t_source, int colorvec_src)
156  {
157  QDPIO::cout << __func__ << ": on t_source= " << t_source << " colorvec_src= " << colorvec_src << std::endl;
158 
159  // Get the source std::vector
160  KeyTimeSliceColorVec_t src_key = getSrcKey(t_source, colorvec_src);
161  LatticeColorVectorF vec_srce = zero;
162 
163  TimeSliceIO<LatticeColorVectorF> time_slice_io(vec_srce, t_source);
164  source_obj.get(src_key, time_slice_io);
165 
166  return vec_srce;
167  }
168 
169 
170  //----------------------------------------------------------------------------
171  //! Get active time-slices
172  std::vector<bool> getActiveTSlices(int t_source, int Nt_forward, int Nt_backward)
173  {
174  // Initialize the active time slices
175  const int decay_dir = Nd-1;
176  const int Lt = Layout::lattSize()[decay_dir];
177 
178  std::vector<bool> active_t_slices(Lt);
179  for(int t=0; t < Lt; ++t)
180  {
181  active_t_slices[t] = false;
182  }
183 
184  // Forward
185  for(int dt=0; dt < Nt_forward; ++dt)
186  {
187  int t = t_source + dt;
188  active_t_slices[t % Lt] = true;
189  }
190 
191  // Backward
192  for(int dt=0; dt < Nt_backward; ++dt)
193  {
194  int t = t_source - dt;
195  while (t < 0) {t += Lt;}
196 
197  active_t_slices[t % Lt] = true;
198  }
199 
200  return active_t_slices;
201  }
202 
203 
204  //----------------------------------------------------------------------------
205  //! Get source keys
206  std::list<KeyTimeSliceColorVec_t> getSrcKeys(int t_source, int colorvec_src)
207  {
208  std::list<KeyTimeSliceColorVec_t> keys;
209 
210  keys.push_back(getSrcKey(t_source,colorvec_src));
211 
212  return keys;
213  }
214 
215 
216  //----------------------------------------------------------------------------
217  //! Get sink keys
218  std::list<KeyPropDistillation_t> getSnkKeys(int t_source, int colorvec_src, int Nt_forward, int Nt_backward, const std::string mass)
219  {
220  std::list<KeyPropDistillation_t> keys;
221 
222  std::vector<bool> active_t_slices = getActiveTSlices(t_source, Nt_forward, Nt_backward);
223 
224  const int decay_dir = Nd-1;
225  const int Lt = Layout::lattSize()[decay_dir];
226 
227  for(int spin_source=0; spin_source < Ns; ++spin_source)
228  {
229  for(int spin_sink=0; spin_sink < Ns; ++spin_sink)
230  {
231  for(int t=0; t < Lt; ++t)
232  {
233  if (! active_t_slices[t]) {continue;}
234 
236 
237  key.t_source = t_source;
238  key.t_slice = t;
239  key.colorvec_src = colorvec_src;
240  key.spin_src = spin_source;
241  key.spin_snk = spin_sink;
242  key.mass = mass;
243 
244  // QDPIO::cout << key << std::flush;
245 
246  keys.push_back(key);
247  } // for t
248  } // for spin_sink
249  } // for spin_source
250 
251  return keys;
252  }
253 
254  } // end anonymous
255 #else
256  // Anonymous namespace
257  namespace
258  {
259  // Convenience type
260  typedef QDP::MapObjectDisk<KeyPropDistillation_t, TimeSliceIO<LatticeColorVectorF> > MOD_t;
261 
262  //----------------------------------------------------------------------------
263  //! Get source key
264  KeyPropDistillation_t getSrcKey(int t_source, int colorvec_src)
265  {
267 
268  key.t_source = t_source;
269  key.t_slice = t_source;
270  key.colorvec_src = colorvec_src;
271  key.spin_src = -1;
272  key.spin_snk = -1;
273  // key.mass NOTE: no mass key in source
274 
275  return key;
276  }
277 
278 
279  //----------------------------------------------------------------------------
280  //! Read a source std::vector
281  LatticeColorVector getSrc(MOD_t& source_obj, int t_source, int colorvec_src)
282  {
283  QDPIO::cout << __func__ << ": on t_source= " << t_source << " colorvec_src= " << colorvec_src << std::endl;
284 
285  // Get the source std::vector
286  KeyPropDistillation_t src_key = getSrcKey(t_source, colorvec_src);
287  LatticeColorVectorF vec_srce = zero;
288 
289  TimeSliceIO<LatticeColorVectorF> time_slice_io(vec_srce, t_source);
290  source_obj.get(src_key, time_slice_io);
291 
292  return vec_srce;
293  }
294 
295 
296  //----------------------------------------------------------------------------
297  //! Get active time-slices
298  std::vector<bool> getActiveTSlices(int t_source, int Nt_forward, int Nt_backward)
299  {
300  // Initialize the active time slices
301  const int decay_dir = Nd-1;
302  const int Lt = Layout::lattSize()[decay_dir];
303 
304  std::vector<bool> active_t_slices(Lt);
305  for(int t=0; t < Lt; ++t)
306  {
307  active_t_slices[t] = false;
308  }
309 
310  // Forward
311  for(int dt=0; dt < Nt_forward; ++dt)
312  {
313  int t = t_source + dt;
314  active_t_slices[t % Lt] = true;
315  }
316 
317  // Backward
318  for(int dt=0; dt < Nt_backward; ++dt)
319  {
320  int t = t_source - dt;
321  while (t < 0) {t += Lt;}
322 
323  active_t_slices[t % Lt] = true;
324  }
325 
326  return active_t_slices;
327  }
328 
329 
330  //----------------------------------------------------------------------------
331  //! Get source keys
332  std::list<KeyPropDistillation_t> getSrcKeys(int t_source, int colorvec_src)
333  {
334  std::list<KeyPropDistillation_t> keys;
335 
336  keys.push_back(getSrcKey(t_source,colorvec_src));
337 
338  return keys;
339  }
340 
341 
342  //----------------------------------------------------------------------------
343  //! Get sink keys
344  std::list<KeyPropDistillation_t> getSnkKeys(int t_source, int colorvec_src, int Nt_forward, int Nt_backward, const std::string mass)
345  {
346  std::list<KeyPropDistillation_t> keys;
347 
348  std::vector<bool> active_t_slices = getActiveTSlices(t_source, Nt_forward, Nt_backward);
349 
350  const int decay_dir = Nd-1;
351  const int Lt = Layout::lattSize()[decay_dir];
352 
353  for(int spin_source=0; spin_source < Ns; ++spin_source)
354  {
355  for(int spin_sink=0; spin_sink < Ns; ++spin_sink)
356  {
357  for(int t=0; t < Lt; ++t)
358  {
359  if (! active_t_slices[t]) {continue;}
360 
362 
363  key.t_source = t_source;
364  key.t_slice = t;
365  key.colorvec_src = colorvec_src;
366  key.spin_src = spin_source;
367  key.spin_snk = spin_sink;
368  key.mass = mass;
369 
370  // QDPIO::cout << key << std::flush;
371 
372  keys.push_back(key);
373  } // for t
374  } // for spin_sink
375  } // for spin_source
376 
377  return keys;
378  }
379 
380  } // end anonymous
381 #endif
382 
383  } // end namespace
384 
385 
386  //----------------------------------------------------------------------------
387  namespace InlinePropDistillationEnv
388  {
389  namespace
390  {
391  AbsInlineMeasurement* createMeasurement(XMLReader& xml_in,
392  const std::string& path)
393  {
394  return new InlineMeas(Params(xml_in, path));
395  }
396 
397  //! Local registration flag
398  bool registered = false;
399  }
400 
401  const std::string name = "PROP_DISTILLATION";
402 
403  //! Register all the factories
404  bool registerAll()
405  {
406  bool success = true;
407  if (! registered)
408  {
410  success &= TheInlineMeasurementFactory::Instance().registerObject(name, createMeasurement);
411  registered = true;
412  }
413  return success;
414  }
415 
416 
417  //----------------------------------------------------------------------------
418  // Param stuff
420 
421  Params::Params(XMLReader& xml_in, const std::string& path)
422  {
423  try
424  {
425  XMLReader paramtop(xml_in, path);
426 
427  if (paramtop.count("Frequency") == 1)
428  read(paramtop, "Frequency", frequency);
429  else
430  frequency = 1;
431 
432  // Parameters for source construction
433  read(paramtop, "Param", param);
434 
435  // Read in the output propagator/source configuration info
436  read(paramtop, "NamedObject", named_obj);
437 
438  // Possible alternate XML file pattern
439  if (paramtop.count("xml_file") != 0)
440  {
441  read(paramtop, "xml_file", xml_file);
442  }
443  }
444  catch(const std::string& e)
445  {
446  QDPIO::cerr << __func__ << ": Caught Exception reading XML: " << e << std::endl;
447  QDP_abort(1);
448  }
449  }
450 
451 
452  //----------------------------------------------------------------------------
453  //----------------------------------------------------------------------------
454  // Function call
455  void
456  InlineMeas::operator()(unsigned long update_no,
457  XMLWriter& xml_out)
458  {
459  // If xml file not empty, then use alternate
460  if (params.xml_file != "")
461  {
462  std::string xml_file = makeXMLFileName(params.xml_file, update_no);
463 
464  push(xml_out, "PropDistillation");
465  write(xml_out, "update_no", update_no);
466  write(xml_out, "xml_file", xml_file);
467  pop(xml_out);
468 
469  XMLFileWriter xml(xml_file);
470  func(update_no, xml);
471  }
472  else
473  {
474  func(update_no, xml_out);
475  }
476  }
477 
478 
479  // Real work done here
480  void
481  InlineMeas::func(unsigned long update_no,
482  XMLWriter& xml_out)
483  {
484  START_CODE();
485 
486  StopWatch snoop;
487  snoop.reset();
488  snoop.start();
489 
490  // Test and grab a reference to the gauge field
491  multi1d<LatticeColorMatrix> u;
492  XMLBufferWriter gauge_xml;
493  try
494  {
495  u = TheNamedObjMap::Instance().getData< multi1d<LatticeColorMatrix> >(params.named_obj.gauge_id);
496  TheNamedObjMap::Instance().get(params.named_obj.gauge_id).getRecordXML(gauge_xml);
497  }
498  catch( std::bad_cast )
499  {
500  QDPIO::cerr << name << ": caught dynamic cast error" << std::endl;
501  QDP_abort(1);
502  }
503  catch (const std::string& e)
504  {
505  QDPIO::cerr << name << ": std::map call failed: " << e << std::endl;
506  QDP_abort(1);
507  }
508 
509  push(xml_out, "PropDistillation");
510  write(xml_out, "update_no", update_no);
511 
512  QDPIO::cout << name << ": propagator calculation" << std::endl;
513 
514  proginfo(xml_out); // Print out basic program info
515 
516  // Write out the input
517  write(xml_out, "Input", params);
518 
519  // Write out the config header
520  write(xml_out, "Config_info", gauge_xml);
521 
522  push(xml_out, "Output_version");
523  write(xml_out, "out_version", 1);
524  pop(xml_out);
525 
526  // Calculate some gauge invariant observables just for info.
527  MesPlq(xml_out, "Observables", u);
528 
529  // Will use TimeSliceSet-s a lot
530  const int decay_dir = params.param.contract.decay_dir;
531  const int Lt = Layout::lattSize()[decay_dir];
532 
533  // A sanity check
534  if (decay_dir != Nd-1)
535  {
536  QDPIO::cerr << __func__ << ": TimeSliceIO only supports decay_dir= " << Nd-1 << "\n";
537  QDP_abort(1);
538  }
539 
540  // The time-slice set
541  TimeSliceSet time_slice_set(decay_dir);
542 
543 
544  //
545  // Map-object-disk storage of the source file
546  //
547  //QDP::MapObjectDisk<KeyTimeSliceColorVec_t, TimeSliceIO<LatticeColorVectorF> > source_obj;
548  MODS_t source_obj;
549  source_obj.setDebug(0);
550 
551  try
552  {
553  QDPIO::cout << "Open source file" << std::endl;
554  source_obj.open(params.named_obj.colorvec_files);
555 
556  QDPIO::cout << "Finished opening solution file" << std::endl;
557 
558  std::string eigen_meta_data; // holds the eigenvalues
559  QDPIO::cout << "Get user data" << std::endl;
560  source_obj.getUserdata(eigen_meta_data);
561  }
562  catch (std::bad_cast) {
563  QDPIO::cerr << name << ": caught dynamic cast error" << std::endl;
564  QDP_abort(1);
565  }
566  catch (const std::string& e) {
567  QDPIO::cerr << name << ": error extracting source_header: " << e << std::endl;
568  QDP_abort(1);
569  }
570  catch( const char* e) {
571  QDPIO::cerr << name <<": Caught some char* exception:" << std::endl;
572  QDPIO::cerr << e << std::endl;
573  QDPIO::cerr << "Rethrowing" << std::endl;
574  throw;
575  }
576 
577  QDPIO::cout << "Source successfully read and parsed" << std::endl;
578 
579 
580  //
581  // Map-object-disk storage
582  //
583  QDP::MapObjectDisk<KeyPropDistillation_t, TimeSliceIO<LatticeColorVectorF> > prop_obj;
584  prop_obj.setDebug(0);
585 
586  QDPIO::cout << "Open solution file" << std::endl;
587 
588  if (! prop_obj.fileExists(params.named_obj.soln_file))
589  {
590  XMLBufferWriter file_xml;
591 
592  push(file_xml, "MODMetaData");
593  write(file_xml, "id", std::string("propDistillation"));
594  write(file_xml, "lattSize", QDP::Layout::lattSize());
595  write(file_xml, "decay_dir", decay_dir);
596  write(file_xml, "num_vecs", params.param.contract.num_vecs);
597  write(file_xml, "Nt_backward", params.param.contract.Nt_backward);
598  write(file_xml, "Nt_backward", params.param.contract.Nt_backward);
599  write(file_xml, "mass_label", params.param.contract.mass_label);
600  proginfo(file_xml); // Print out basic program info
601  write(file_xml, "Params", params.param);
602  write(file_xml, "Config_info", gauge_xml);
603  pop(file_xml);
604 
605  std::string file_str(file_xml.str());
606 
607  prop_obj.insertUserdata(file_xml.str());
608  prop_obj.open(params.named_obj.soln_file, std::ios_base::in | std::ios_base::out | std::ios_base::trunc);
609  }
610  else
611  {
612  prop_obj.open(params.named_obj.soln_file);
613  }
614 
615  QDPIO::cout << "Finished opening solution file" << std::endl;
616 
617 
618  // Total number of iterations
619  int ncg_had = 0;
620 
621 
622 #if 0
623  // NOTE: not doing this spin rotation, but leaving the code here for reference
624  // Rotation from DR to DP
625  SpinMatrix diracToDRMat(DiracToDRMat());
626  std::vector<MatrixSpinRep_t> diracToDrMatPlus = convertTwoQuarkSpin(diracToDRMat);
627  std::vector<MatrixSpinRep_t> diracToDrMatMinus = convertTwoQuarkSpin(adj(diracToDRMat));
628 #endif
629 
630 
631  //
632  // Try the factories
633  //
634  try
635  {
636  StopWatch swatch;
637  swatch.reset();
638  QDPIO::cout << "Try the various factories" << std::endl;
639 
640  // Typedefs to save typing
641  typedef LatticeFermion T;
642  typedef multi1d<LatticeColorMatrix> P;
643  typedef multi1d<LatticeColorMatrix> Q;
644 
645  //
646  // Initialize fermion action
647  //
648  std::istringstream xml_s(params.param.prop.fermact.xml);
649  XMLReader fermacttop(xml_s);
650  QDPIO::cout << "FermAct = " << params.param.prop.fermact.id << std::endl;
651 
652  // Generic Wilson-Type stuff
655  fermacttop,
657 
658  Handle< FermState<T,P,Q> > state(S_f->createState(u));
659 
662 
663  QDPIO::cout << "Suitable factory found: compute all the quark props" << std::endl;
664  swatch.start();
665 
666 
667  //
668  // Loop over the source color and spin, creating the source
669  // and calling the relevant propagator routines.
670  //
671  const int num_vecs = params.param.contract.num_vecs;
672  const multi1d<int>& t_sources = params.param.contract.t_sources;
673 
674  // Loop over each time source
675  for(int tt=0; tt < t_sources.size(); ++tt)
676  {
677  int t_source = t_sources[tt]; // This is the actual time-slice.
678  QDPIO::cout << "t_source = " << t_source << std::endl;
679 
680  // The space distillation loop
681  for(int colorvec_src=0; colorvec_src < num_vecs; ++colorvec_src)
682  {
683  StopWatch sniss1;
684  sniss1.reset();
685  sniss1.start();
686  QDPIO::cout << "colorvec_src = " << colorvec_src << std::endl;
687 
688  // Get the source std::vector
689  LatticeColorVector vec_srce = getSrc(source_obj, t_source, colorvec_src);
690 
691  //
692  // Loop over each spin source and invert.
693  // Use the same colorstd::vector source. No spin dilution will be used.
694  //
695  multi2d<LatticeColorVector> ferm_out(Ns,Ns);
696 
697  for(int spin_source=0; spin_source < Ns; ++spin_source)
698  {
699  QDPIO::cout << "spin_source = " << spin_source << std::endl;
700 
701  // Insert a ColorVector into spin index spin_source
702  // This only overwrites sections, so need to initialize first
703  LatticeFermion chi = zero;
704  CvToFerm(vec_srce, chi, spin_source);
705 
706  LatticeFermion quark_soln = zero;
707 
708  // Do the propagator inversion
709  SystemSolverResults_t res = (*PP)(quark_soln, chi);
710  ncg_had = res.n_count;
711 
712  // Extract into the temporary output array
713  for(int spin_sink=0; spin_sink < Ns; ++spin_sink)
714  {
715  ferm_out(spin_sink,spin_source) = peekSpin(quark_soln, spin_sink);
716  }
717  } // for spin_source
718 
719 
720 #if 0
721  // NOTE: not doing this spin rotation, but leaving the code here for reference
722  // Instead, leave the code in the DR basis - the convention for all the perambulators
723  // written
724 
725  // Rotate from DeGrand-Rossi (DR) to Dirac-Pauli (DP)
726  {
727  multi2d<LatticeColorVector> ferm_tmp;
728 
729  multiplyRep(ferm_tmp, diracToDrMatMinus, ferm_out);
730  multiplyRep(ferm_out, ferm_tmp, diracToDrMatPlus);
731  }
732 #endif
733 
734  sniss1.stop();
735  QDPIO::cout << "Time to assemble and transmogrify propagators for colorvec_src= " << colorvec_src << " time = "
736  << sniss1.getTimeInSeconds()
737  << " secs" << std::endl;
738 
739 
740  // Write out each time-slice chunk of a lattice colorvec soln to disk
741  QDPIO::cout << "Write propagator solutions to disk" << std::endl;
742  StopWatch sniss2;
743  sniss2.reset();
744  sniss2.start();
745 
746  // Write the solutions
747  QDPIO::cout << "Write propagator solution to disk" << std::endl;
748  std::list<KeyPropDistillation_t> snk_keys(getSnkKeys(t_source, colorvec_src,
752 
753  for(std::list<KeyPropDistillation_t>::const_iterator key= snk_keys.begin();
754  key != snk_keys.end();
755  ++key)
756  {
757  LatticeColorVectorF tmptmp = ferm_out(key->spin_snk,key->spin_src);
758 
759  prop_obj.insert(*key, TimeSliceIO<LatticeColorVectorF>(tmptmp, key->t_slice));
760  } // for key
761 
762  sniss2.stop();
763  QDPIO::cout << "Time to write propagators for colorvec_src= " << colorvec_src << " time = "
764  << sniss2.getTimeInSeconds()
765  << " secs" << std::endl;
766 
767  } // for colorvec_src
768  } // for tt
769 
770  swatch.stop();
771  QDPIO::cout << "Propagators computed: time= "
772  << swatch.getTimeInSeconds()
773  << " secs" << std::endl;
774  }
775  catch (const std::string& e)
776  {
777  QDPIO::cout << name << ": caught exception around qprop: " << e << std::endl;
778  QDP_abort(1);
779  }
780 
781  push(xml_out,"Relaxation_Iterations");
782  write(xml_out, "ncg_had", ncg_had);
783  pop(xml_out);
784 
785  pop(xml_out); // prop_dist
786 
787  snoop.stop();
788  QDPIO::cout << name << ": total time = "
789  << snoop.getTimeInSeconds()
790  << " secs" << std::endl;
791 
792  QDPIO::cout << name << ": ran successfully" << std::endl;
793 
794  END_CODE();
795  }
796 
797  }
798 
799 } // namespace Chroma
Inline measurement factory.
Class for counted reference semantics.
Definition: handle.h:33
void operator()(const unsigned long update_no, XMLWriter &xml_out)
Do the measurement.
void func(const unsigned long update_no, XMLWriter &xml_out)
Do the measurement.
static T & Instance()
Definition: singleton.h:432
Builds time slice subsets.
EXTERN Real dt
Basis rotation matrix from Dirac to Degrand-Rossi (and reverse)
Class structure for fermion actions.
Fermion action factories.
All Wilson-type fermion actions.
SpinMatrixD DiracToDRMat()
The Dirac to Degrand-Rossi spin transformation matrix (and reverse)
Definition: diractodr.cc:22
void CvToFerm(const LatticeColorVectorF &a, LatticeFermionF &b, int spin_index)
Convert (insert) a LatticeColorVector into a LatticeFermion.
Definition: transf.cc:18
void proginfo(XMLWriter &xml)
Print out basic information about this program.
Definition: proginfo.cc:24
std::string makeXMLFileName(std::string xml_file, unsigned long update_no)
Return a xml file name for inline measurements.
TimeSliceSet time_slice_set
Compute the propagator from distillation.
Key for vanilla distillation propagator sources and solutions.
Key for time-sliced color eigenvectors.
Make xml file writer.
int t
Definition: meslate.cc:37
Nd
Definition: meslate.cc:74
Named object function std::map.
static bool registered
Local registration flag.
multi1d< LatticeColorMatrix > P
QDP::MapObjectDiskMultiple< KeyTimeSliceColorVec_t, TimeSliceIO< LatticeColorVectorF > > MODS_t
bool registerAll()
Register all the factories.
void write(XMLWriter &xml, const std::string &path, const Params::NamedObject_t &input)
Propagator output.
void read(XMLReader &xml, const std::string &path, Params::NamedObject_t &input)
Propagator input.
bool registerAll()
Register all the factories.
Asqtad Staggered-Dirac operator.
Definition: klein_gord.cc:10
static multi1d< LatticeColorMatrix > u
LinOpSysSolverMGProtoClover::Q Q
LatticeFermion tmp
Definition: mespbg5p_w.cc:36
push(xml_out,"Condensates")
LinOpSysSolverMGProtoClover::T T
Double mass
Definition: pbg5p_w.cc:54
std::vector< MatrixSpinRep_t > convertTwoQuarkSpin(const SpinMatrix &in)
Convert a generic spin matrix into a sparse spin representation.
Definition: spin_rep.cc:67
multi1d< LatticeFermion > chi(Ncb)
QDP::MapObjectDisk< KeyPropDistillution_t, TimeSliceIO< LatticeColorVectorF > > MOD_t
pop(xml_out)
void MesPlq(const multi1d< LatticeColorMatrixF3 > &u, multi2d< Double > &plane_plaq, Double &link)
Definition: mesplq.cc:83
static QDP_ColorVector * out
Constructor.
START_CODE()
const WilsonTypeFermAct< multi1d< LatticeFermion > > Handle< const ConnectState > state
Definition: pbg5p_w.cc:28
const WilsonTypeFermAct< multi1d< LatticeFermion > > & S_f
Definition: pbg5p_w.cc:27
Double zero
Definition: invbicg.cc:106
void multiplyRep(multi2d< LatticeColorVector > &dist_rep, const std::vector< MatrixSpinRep_t > &spin, const multi2d< LatticeColorVector > &prop_rep)
Dist(t2) = SpinMatrix*Prop(t2)
static QDP_ColorVector * in
::std::string string
Definition: gtest.h:1979
Print out basic info about this program.
Sparse matrix representation of spin matrices.
GroupXML_t invParam
Definition: qprop_io.h:84
GroupXML_t fermact
Definition: qprop_io.h:80
Distillation propagators.
Holds return info from SystemSolver call.
Definition: syssolver.h:17
Convenience for building time-slice subsets.
Contraction operators for two quarks.