CHROMA
inline_matelem_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 "qdp.h"
8 #include "fermact.h"
11 #include "meas/glue/mesplq.h"
12 #include "qdp_map_obj.h"
13 #include "qdp_map_obj_disk.h"
14 #include "qdp_map_obj_disk_multiple.h"
15 #include "qdp_map_obj_memory.h"
16 #include "qdp_disk_map_slice.h"
22 #include "util/ferm/key_val_db.h"
23 #include "util/ferm/transf.h"
24 #include "util/ferm/spin_rep.h"
25 #include "util/ferm/diractodr.h"
27 #include "util/ft/sftmom.h"
28 #include "util/ft/time_slice_set.h"
29 #include "util/info/proginfo.h"
30 #include "util/info/proginfo.h"
34 
36 
37 #include "chroma_config.h"
38 
39 #ifndef QDP_IS_QDPJIT_NO_NVPTX
40 
41 #ifdef BUILD_JIT_CONTRACTION_KERNELS
42 #include "custom_kernels/custom_kernels.h"
43 #endif
44 
45 namespace Chroma
46 {
47  //----------------------------------------------------------------------------
48  // Utilities
49  namespace
50  {
51 
52  multi1d<SubsetVectorWeight_t> readEigVals(const std::string& meta)
53  {
54  std::istringstream xml_l(meta);
55  XMLReader eigtop(xml_l);
56 
57  std::string pat("/MODMetaData/Weights");
58  // multi1d<SubsetVectorWeight_t> eigenvalues;
59  multi1d< multi1d<Real> > eigens;
60  try
61  {
62  read(eigtop, pat, eigens);
63  }
64  catch(const std::string& e)
65  {
66  QDPIO::cerr << __func__ << ": Caught Exception reading meta= XX" << meta << "XX with path= " << pat << " error= " << e << std::endl;
67  QDP_abort(1);
68  }
69 
70  eigtop.close();
71 
72  multi1d<SubsetVectorWeight_t> eigenvalues(eigens.size());
73 
74  for(int i=0; i < eigens.size(); ++i)
75  eigenvalues[i].weights = eigens[i];
76 
77  return eigenvalues;
78  }
79  } // end anonymous namespace
80 
81 
82 
83  //----------------------------------------------------------------------------
84  namespace InlineMatElemDistillationEnv
85  {
86  //! Propagator input
87  void read(XMLReader& xml, const std::string& path, Params::NamedObject_t& input)
88  {
89  XMLReader inputtop(xml, path);
90 
91  read(inputtop, "colorvec_files", input.colorvec_files);
92  read(inputtop, "prop_file", input.prop_file);
93  read(inputtop, "peram_file", input.peram_file);
94  }
95 
96  //! Propagator output
97  void write(XMLWriter& xml, const std::string& path, const Params::NamedObject_t& input)
98  {
99  push(xml, path);
100 
101  write(xml, "colorvec_files", input.colorvec_files);
102  write(xml, "prop_file", input.prop_file);
103  write(xml, "peram_file", input.peram_file);
104 
105  pop(xml);
106  }
107 
108 
109  //! Propagator input
110  void read(XMLReader& xml, const std::string& path, Params::Param_t::Contract_t& input)
111  {
112  XMLReader inputtop(xml, path);
113 
114  read(inputtop, "num_vecs", input.num_vecs);
115  read(inputtop, "lt_orig", input.lt_orig);
116  read(inputtop, "t_offset", input.t_offset);
117  read(inputtop, "t_sources", input.t_sources);
118  read(inputtop, "decay_dir", input.decay_dir);
119  read(inputtop, "Nt_forward", input.Nt_forward);
120  read(inputtop, "Nt_backward", input.Nt_backward);
121  read(inputtop, "mass_label", input.mass_label);
122 
123  input.zero_colorvecs = false;
124  if( inputtop.count("zero_colorvecs") == 1 ) {
125  read(inputtop, "zero_colorvecs", input.zero_colorvecs );
126  if (input.zero_colorvecs)
127  {
128  QDPIO::cout << "zero_colorvecs found, *** timing mode activated ***\n";
129  }
130  }
131 
132  input.fuse_timeloop = false;
133  if( inputtop.count("fuse_timeloop") == 1 ) {
134  read(inputtop, "fuse_timeloop", input.fuse_timeloop );
135  if (input.fuse_timeloop)
136  {
137  QDPIO::cout << "fuse_timeloop found!\n";
138  }
139  }
140 
141  }
142 
143  //! Propagator output
144  void write(XMLWriter& xml, const std::string& path, const Params::Param_t::Contract_t& input)
145  {
146  push(xml, path);
147 
148  write(xml, "num_vecs", input.num_vecs);
149  write(xml, "t_sources", input.t_sources);
150  write(xml, "decay_dir", input.decay_dir);
151  write(xml, "Nt_forward", input.Nt_forward);
152  write(xml, "Nt_backward", input.Nt_backward);
153  write(xml, "mass_label", input.mass_label);
154 
155  pop(xml);
156  }
157 
158 
159  //! Propagator input
160  void read(XMLReader& xml, const std::string& path, Params::Param_t& input)
161  {
162  XMLReader inputtop(xml, path);
163 
164  read(inputtop, "Contractions", input.contract);
165  }
166 
167  //! Propagator output
168  void write(XMLWriter& xml, const std::string& path, const Params::Param_t& input)
169  {
170  push(xml, path);
171 
172  write(xml, "Contractions", input.contract);
173 
174  pop(xml);
175  }
176 
177 
178  //! Propagator input
179  void read(XMLReader& xml, const std::string& path, Params& input)
180  {
181  Params tmp(xml, path);
182  input = tmp;
183  }
184 
185  //! Propagator output
186  void write(XMLWriter& xml, const std::string& path, const Params& input)
187  {
188  push(xml, path);
189 
190  write(xml, "Param", input.param);
191  write(xml, "NamedObject", input.named_obj);
192 
193  pop(xml);
194  }
195  } // namespace InlinePropDistillationEnv
196 
197 
198  //----------------------------------------------------------------------------
199  //----------------------------------------------------------------------------
200  namespace InlineMatElemDistillationEnv
201  {
202  //----------------------------------------------------------------------------
203  // Convenience type
204  typedef QDP::MapObjectDisk<KeyTimeSliceColorVec_t, TimeSliceIO<LatticeColorVectorF> > MOD_t;
205 
206  // Convenience type
207  typedef QDP::MapObjectDiskMultiple<KeyTimeSliceColorVec_t, TimeSliceIO<LatticeColorVectorF> > MODS_t;
208 
209  // Convenience type
210  typedef QDP::MapObjectMemory<KeyTimeSliceColorVec_t, SubLatticeColorVectorF> SUB_MOD_t;
211 
212 
213 
214  // Anonymous namespace
215  namespace
216  {
217  //----------------------------------------------------------------------------
218  //! Class to hold std::map of eigenvectors
219  class SubEigenMap
220  {
221  public:
222  //! Constructor
223  SubEigenMap(MODS_t& eigen_source_, int decay_dir, bool zero_colorvecs, int t_offset, int lt_orig) : eigen_source(eigen_source_), time_slice_set(decay_dir), zero_colorvecs(zero_colorvecs), t_offset(t_offset), lt_orig(lt_orig) {}
224 
225  //! Getter
226  const SubLatticeColorVectorF& getVec(int t_source, int colorvec_src) const;
227 
228  //! The set to be used in sumMulti
229  const Set& getSet() const {return time_slice_set.getSet();}
230 
231  private:
232  //! Eigenvectors
234 
235  // The time-slice set
237 
238  private:
239  //! Where we store the sublattice versions
242  int t_offset;
243  int lt_orig;
244  };
245 
246  //----------------------------------------------------------------------------
247  //! Getter
248  const SubLatticeColorVectorF& SubEigenMap::getVec(int t_source, int colorvec_src) const
249  {
250  // The key
251  KeyTimeSliceColorVec_t src_key( (t_source + t_offset ) % lt_orig , colorvec_src);
252 
253  // If item does not exist, read from original std::map and put in memory std::map
254  if (! sub_eigen.exist(src_key))
255  {
256  QDPIO::cout << __func__ << ": on t_source= " << t_source << " colorvec_src= " << colorvec_src << std::endl;
257 
258  // No need to initialize with 'zero' - we are returning a subtype.
259  LatticeColorVectorF vec_srce;
260 
261  if (!zero_colorvecs)
262  {
263  TimeSliceIO<LatticeColorVectorF> time_slice_io(vec_srce, t_source );
264  eigen_source.get(src_key, time_slice_io);
265  }
266 
267  SubLatticeColorVectorF tmp(getSet()[t_source], vec_srce);
268 
269  sub_eigen.insert(src_key, tmp);
270  }
271  return sub_eigen[src_key];
272  }
273 
274 
275  //----------------------------------------------------------------------------
276  //! Get active time-slices
277  std::vector<bool> getActiveTSlices(int t_source, int Nt_forward, int Nt_backward)
278  {
279  // Initialize the active time slices
280  const int decay_dir = Nd-1;
281  const int Lt = Layout::lattSize()[decay_dir];
282 
283  std::vector<bool> active_t_slices(Lt);
284  for(int t=0; t < Lt; ++t)
285  {
286  active_t_slices[t] = false;
287  }
288 
289  // Forward
290  for(int dt=0; dt < Nt_forward; ++dt)
291  {
292  int t = t_source + dt;
293  active_t_slices[t % Lt] = true;
294  }
295 
296  // Backward
297  for(int dt=0; dt < Nt_backward; ++dt)
298  {
299  int t = t_source - dt;
300  while (t < 0) {t += Lt;}
301 
302  active_t_slices[t % Lt] = true;
303  }
304 
305  return active_t_slices;
306  }
307 
308 
309  //----------------------------------------------------------------------------
310  //! Get sink keys
311  std::list<KeyPropElementalOperator_t> getSnkKeys(int t_source, int spin_source, int Nt_forward, int Nt_backward, const std::string mass)
312  {
313  std::list<KeyPropElementalOperator_t> keys;
314 
315  std::vector<bool> active_t_slices = getActiveTSlices(t_source, Nt_forward, Nt_backward);
316 
317  const int decay_dir = Nd-1;
318  const int Lt = Layout::lattSize()[decay_dir];
319 
320  for(int spin_sink=0; spin_sink < Ns; ++spin_sink)
321  {
322  for(int t=0; t < Lt; ++t)
323  {
324  if (! active_t_slices[t]) {continue;}
325 
327 
328  key.t_source = t_source;
329  key.t_slice = t;
330  key.spin_src = spin_source;
331  key.spin_snk = spin_sink;
332  key.mass_label = mass;
333 
334  keys.push_back(key);
335  } // for t
336  } // for spin_sink
337 
338  return keys;
339  }
340 
341  } // end anonymous
342  } // end namespace
343 
344 
345  //----------------------------------------------------------------------------
346  namespace InlineMatElemDistillationEnv
347  {
348  namespace
349  {
350  AbsInlineMeasurement* createMeasurement(XMLReader& xml_in,
351  const std::string& path)
352  {
353  return new InlineMeas(Params(xml_in, path));
354  }
355 
356  //! Local registration flag
357  bool registered = false;
358  }
359 
360  const std::string name = "MATELEM_DISTILLATION";
361 
362  //! Register all the factories
363  bool registerAll()
364  {
365  bool success = true;
366  if (! registered)
367  {
369  success &= TheInlineMeasurementFactory::Instance().registerObject(name, createMeasurement);
370  registered = true;
371  }
372  return success;
373  }
374 
375 
376  //----------------------------------------------------------------------------
377  // Param stuff
379 
380  Params::Params(XMLReader& xml_in, const std::string& path)
381  {
382  try
383  {
384  XMLReader paramtop(xml_in, path);
385 
386  if (paramtop.count("Frequency") == 1)
387  read(paramtop, "Frequency", frequency);
388  else
389  frequency = 1;
390 
391  // Parameters for source construction
392  read(paramtop, "Param", param);
393 
394  // Read in the output propagator/source configuration info
395  read(paramtop, "NamedObject", named_obj);
396 
397  // Possible alternate XML file pattern
398  if (paramtop.count("xml_file") != 0)
399  {
400  read(paramtop, "xml_file", xml_file);
401  }
402  }
403  catch(const std::string& e)
404  {
405  QDPIO::cerr << __func__ << ": Caught Exception reading XML: " << e << std::endl;
406  QDP_abort(1);
407  }
408  }
409 
410 
411  //----------------------------------------------------------------------------
412  //----------------------------------------------------------------------------
413  // Function call
414  void
415  InlineMeas::operator()(unsigned long update_no,
416  XMLWriter& xml_out)
417  {
418  // If xml file not empty, then use alternate
419  if (params.xml_file != "")
420  {
421  std::string xml_file = makeXMLFileName(params.xml_file, update_no);
422 
423  push(xml_out, "PropDistillation");
424  write(xml_out, "update_no", update_no);
425  write(xml_out, "xml_file", xml_file);
426  pop(xml_out);
427 
428  XMLFileWriter xml(xml_file);
429  func(update_no, xml);
430  }
431  else
432  {
433  func(update_no, xml_out);
434  }
435  }
436 
437 
438  // Real work done here
439  void
440  InlineMeas::func(unsigned long update_no,
441  XMLWriter& xml_out)
442  {
443  START_CODE();
444 
445  StopWatch snoop;
446  snoop.reset();
447  snoop.start();
448 
449  push(xml_out, "PropDistillation");
450  write(xml_out, "update_no", update_no);
451 
452  proginfo(xml_out); // Print out basic program info
453 
454  // Write out the input
455  write(xml_out, "Input", params);
456 
457  push(xml_out, "Output_version");
458  write(xml_out, "out_version", 1);
459  pop(xml_out);
460 
461  // Will use TimeSliceSet-s a lot
462  const int decay_dir = params.param.contract.decay_dir;
463  const int Lt = Layout::lattSize()[decay_dir];
464  const int lt_orig = params.param.contract.lt_orig;
465  const int t_offset = params.param.contract.t_offset;
466 
467  // A sanity check
468  if (decay_dir != Nd-1)
469  {
470  QDPIO::cerr << name << ": TimeSliceIO only supports decay_dir= " << Nd-1 << "\n";
471  QDP_abort(1);
472  }
473 
474 
475  QDP::MapObjectDisk<KeyPropDistillation_t, TimeSliceIO<LatticeColorVectorF> > prop_obj;
476  prop_obj.setDebug(0);
477 
478  QDPIO::cout << "Open source prop file" << std::endl;
479 
480  if (! prop_obj.fileExists(params.named_obj.prop_file))
481  {
482  QDPIO::cerr << name << ": source file does not exist: src_file= " << params.named_obj.prop_file << std::endl;
483  QDP_abort(1);
484  }
485  else
486  {
488  }
489 
490  QDPIO::cout << "Finished opening solution file" << std::endl;
491 
492 
493  //
494  // Read in the source along with relevant information.
495  //
496  QDPIO::cout << "Snarf the source from a std::map object disk file" << std::endl;
497 
499  eigen_source.setDebug(0);
500 
501  std::string eigen_meta_data; // holds the eigenvalues
502 
504  {
505  try
506  {
507  // Open
508  QDPIO::cout << "Open file= " << params.named_obj.colorvec_files[0] << std::endl;
510 
511  // Snarf the source info.
512  QDPIO::cout << "Get user data" << std::endl;
513  eigen_source.getUserdata(eigen_meta_data);
514  // QDPIO::cout << "User data= " << eigen_meta_data << std::endl;
515 
516  // Write it
517  // QDPIO::cout << "Write to an xml file" << std::endl;
518  // XMLBufferWriter xml_buf(eigen_meta_data);
519  // write(xml_out, "Source_info", xml_buf);
520  }
521  catch (std::bad_cast) {
522  QDPIO::cerr << name << ": caught dynamic cast error" << std::endl;
523  QDP_abort(1);
524  }
525  catch (const std::string& e) {
526  QDPIO::cerr << name << ": error extracting source_header: " << e << std::endl;
527  QDP_abort(1);
528  }
529  catch( const char* e) {
530  QDPIO::cerr << name <<": Caught some char* exception:" << std::endl;
531  QDPIO::cerr << e << std::endl;
532  QDPIO::cerr << "Rethrowing" << std::endl;
533  throw;
534  }
535 
536  QDPIO::cout << "Source successfully read and parsed" << std::endl;
537  }
538 
539 #if 0
540  // Sanity check
542  {
543  QDPIO::cerr << name << ": number of available eigenvectors is too small\n";
544  QDP_abort(1);
545  }
546 #endif
547 
548  QDPIO::cout << "Number of vecs available is large enough" << std::endl;
549 
550  // The sub-lattice eigenstd::vector std::map
551  QDPIO::cout << "Initialize sub-lattice std::map" << std::endl;
552  SubEigenMap sub_eigen_map(eigen_source, decay_dir, params.param.contract.zero_colorvecs, t_offset , lt_orig );
553  QDPIO::cout << "Finished initializing sub-lattice std::map" << std::endl;
554 
555 
556  //
557  // DB storage
558  //
559  BinaryStoreDB< SerialDBKey<KeyPropElementalOperator_t>, SerialDBData<ValPropElementalOperator_t> > qdp_db;
560 
562  {
563  // Open the file, and write the meta-data and the binary for this operator
564  if (! qdp_db.fileExists(params.named_obj.peram_file))
565  {
566  XMLBufferWriter file_xml;
567  push(file_xml, "DBMetaData");
568  write(file_xml, "id", std::string("propElemOp"));
569 
570  auto tmp = QDP::Layout::lattSize();
571  tmp[Nd-1] = lt_orig;
572 
573  write(file_xml, "lattSize", tmp );
574  write(file_xml, "decay_dir", params.param.contract.decay_dir);
575  proginfo(file_xml); // Print out basic program info
576  write(file_xml, "Params", params.param);
578  write(file_xml, "Weights", readEigVals(eigen_meta_data));
579  pop(file_xml);
580 
581  std::string file_str(file_xml.str());
582  qdp_db.setMaxUserInfoLen(file_str.size());
583 
584  qdp_db.open(params.named_obj.peram_file, O_RDWR | O_CREAT, 0664);
585 
586  qdp_db.insertUserdata(file_str);
587  }
588  else
589  {
590  qdp_db.open(params.named_obj.peram_file, O_RDWR, 0664);
591  }
592 
593  QDPIO::cout << "Finished opening peram file" << std::endl;
594  }
595 
596  // Total number of iterations
597  int ncg_had = 0;
598 
599  //
600  // Try the factories
601  //
602  try
603  {
604  StopWatch swatch;
605  swatch.reset();
606  QDPIO::cout << "Try the various factories" << std::endl;
607 
608  swatch.start();
609 
610 #ifdef BUILD_JIT_CONTRACTION_KERNELS
612  QDPIO::cout << "Using JIT contraction kernels (fused timeloop)\n";
613  else
614  QDPIO::cout << "Using JIT contraction kernels (non-fused timeloop)\n";
615 #endif
616 
617  //
618  // Loop over the source color and spin, creating the source
619  // and calling the relevant propagator routines.
620  //
621  const int num_vecs = params.param.contract.num_vecs;
622  const multi1d<int>& t_sources = params.param.contract.t_sources;
623 
624 
625  // Loop over each time source
626  for(int tt=0; tt < t_sources.size(); ++tt)
627  {
628  int t_source = t_sources[tt]; // This is the actual time-slice.
629  QDPIO::cout << "t_source = " << t_source << std::endl;
630 
631  t_source -= t_offset;
632  QDPIO::cout << "new t_source in current lattice = " << t_source << std::endl;
633 
634  // Loop over each spin source
635  for(int spin_source=0; spin_source < Ns; ++spin_source)
636  {
637  QDPIO::cout << "spin_source = " << spin_source << std::endl;
638 
639  // These are the common parts of a perambulator that are needed for this time source
640  std::list<KeyPropElementalOperator_t> snk_keys(getSnkKeys(t_source,
641  spin_source,
645 
646  if (1) {
647  // The final perambulator
648  QDP::MapObjectMemory<KeyPropElementalOperator_t, ValPropElementalOperator_t> peram;
649 
650  // Initialize
651  for(std::list<KeyPropElementalOperator_t>::const_iterator key = snk_keys.begin();
652  key != snk_keys.end();
653  ++key)
654  {
655  // The perambulator value
657  peram.insert(*key, tmp);
658 
659  peram[*key].mat.resize(num_vecs,num_vecs);
660  peram[*key].mat = zero;
661  } // key
662  QDPIO::cout << "peram initialized! " << std::endl;
663 
664  //
665  // The space distillation loop
666  //
667  for(int colorvec_src=0; colorvec_src < num_vecs; ++colorvec_src)
668  {
669  StopWatch sniss1;
670  sniss1.reset();
671  sniss1.start();
672 
673  StopWatch snarss1;
674  snarss1.reset();
675  snarss1.start();
676  QDPIO::cout << "Do spin_source= " << spin_source << " colorvec_src= " << colorvec_src << std::endl;
677 
678 
679  //
680  // Loop over each spin source and invert.
681  // Use the same colorstd::vector source. No spin dilution will be used.
682  //
683  multi1d<LatticeColorVector> ferm_out(Ns);
684 
686  {
687  // Extract into the temporary output array
688  for(int spin_sink=0; spin_sink < Ns; ++spin_sink)
689  {
690  for ( int t = 0 ; t < Lt ; ++t )
691  {
693 
694  key.t_source = ( t_source + t_offset ) % lt_orig;
695  key.t_slice = ( t + t_offset ) % lt_orig;
696  key.colorvec_src = colorvec_src;
697  key.spin_src = spin_source;
698  key.spin_snk = spin_sink;
700 
701  LatticeColorVectorF tmp;
702 
703  TimeSliceIO<LatticeColorVectorF> time_slice_io(tmp, t);
704 
705  prop_obj.get( key , time_slice_io );
706 
707  ferm_out(spin_sink)[ sub_eigen_map.getSet()[t] ] = tmp;
708  }
709  }
710 
711  } // zero_colorvecs ??
712  else
713  {
714  for ( int s = 0 ; s < Ns ; ++s )
715  {
716  zero_rep( ferm_out[s] );
717  }
718  }
719 
720 
721  snarss1.stop();
722  QDPIO::cout << "Time to read prop for spin_source= " << spin_source << " colorvec_src= " << colorvec_src << " time = "
723  << snarss1.getTimeInSeconds()
724  << " secs" << std::endl;
725 
726  // The perambulator part
727  // Loop over time
728 
729 #ifndef BUILD_JIT_CONTRACTION_KERNELS
730  for(int t_slice = 0; t_slice < Lt; ++t_slice)
731  {
732  // Loop over all the keys
733  for(std::list<KeyPropElementalOperator_t>::const_iterator key= snk_keys.begin();
734  key != snk_keys.end();
735  ++key)
736  {
737  if (key->t_slice != t_slice) {continue;}
738 
739  // Loop over the sink colorvec, form the innerproduct and the resulting perambulator
740  for(int colorvec_sink=0; colorvec_sink < num_vecs; ++colorvec_sink)
741  {
742  peram[*key].mat(colorvec_sink,colorvec_src) = innerProduct(sub_eigen_map.getVec(t_slice, colorvec_sink),
743  ferm_out(key->spin_snk));
744  } // for colorvec_sink
745  } // for key
746  } // for t_slice
747 #else
749  {
750  //
751  // fused timeloop
752  //
753  for(int spin_snk = 0; spin_snk < Ns; ++spin_snk)
754  {
755  int count = 0;
756  for(std::list<KeyPropElementalOperator_t>::const_iterator key= snk_keys.begin(); key != snk_keys.end(); ++key)
757  {
758  if (key->spin_snk != spin_snk) { continue;}
759  count += num_vecs;
760  }
761 
762  //QDPIO::cout << "spin_snk = " << spin_snk << ", count = " << count << "\n";
763 
764  multi1d<SubLatticeColorVectorF*> vec_ptr( count );
765  multi1d<ComplexD*> contr_ptr( count );
766 
767  int run_count = 0;
768 
769  // Loop over all the keys
770  for(std::list<KeyPropElementalOperator_t>::const_iterator key= snk_keys.begin();
771  key != snk_keys.end();
772  ++key)
773  {
774  if (key->spin_snk != spin_snk) {continue;}
775 
776  //
777  // Pack pointers to the vectors and matrix elements
778  //
779  for (int i=0 ; i < num_vecs ; ++i ) {
780  vec_ptr[run_count] = const_cast<SubLatticeColorVectorF*>( &sub_eigen_map.getVec( key->t_slice , i ) );
781  contr_ptr[run_count] = &peram[*key].mat( i , colorvec_src );
782  ++run_count;
783  }
784  }
785 
786  //
787  // Big call-out
788  //
789  multi_innerProduct( contr_ptr , vec_ptr , ferm_out(spin_snk) );
790 
791  } // for spin_snk
792  }
793  else
794  {
795  //
796  // non-fused timeloop
797  //
798  for(int t_slice = 0; t_slice < Lt; ++t_slice)
799  {
800  // Loop over all the keys
801  for(std::list<KeyPropElementalOperator_t>::const_iterator key= snk_keys.begin();
802  key != snk_keys.end();
803  ++key)
804  {
805  if (key->t_slice != t_slice) {continue;}
806  //
807  // Pack pointers to the vectors and matrix elements
808  //
809  multi1d<SubLatticeColorVectorF*> vec_ptr( num_vecs );
810  multi1d<ComplexD*> contr_ptr( num_vecs );
811  for (int i=0 ; i < num_vecs ; ++i ) {
812  vec_ptr[i] = const_cast<SubLatticeColorVectorF*>( &sub_eigen_map.getVec( t_slice , i ) );
813  contr_ptr[i] = &peram[*key].mat( i , colorvec_src );
814  }
815 
816  //
817  // Big call-out
818  //
819  multi_innerProduct( contr_ptr , vec_ptr , ferm_out(key->spin_snk) );
820  } // for key
821  } // for t_slice
822  }
823 #endif
824 
825  sniss1.stop();
826  QDPIO::cout << "Time to compute and assemble peram for spin_source= " << spin_source << " colorvec_src= " << colorvec_src << " time = "
827  << sniss1.getTimeInSeconds()
828  << " secs"
829  << " (for contraction: " << sniss1.getTimeInSeconds() - snarss1.getTimeInSeconds() << ")"
830  << std::endl;
831 
832  } // for colorvec_src
833 
835  {
836  // Write out each time-slice chunk of a lattice colorvec soln to disk
837  QDPIO::cout << "Write perambulator for spin_source= " << spin_source << " to disk" << std::endl;
838  StopWatch sniss2;
839  sniss2.reset();
840  sniss2.start();
841 
842  // The perambulator is complete. Write it.
843  for(std::list<KeyPropElementalOperator_t>::const_iterator key= snk_keys.begin();
844  key != snk_keys.end();
845  ++key)
846  {
848  tmp.t_source = ( tmp.t_source + t_offset ) % lt_orig;
849  tmp.t_slice = ( tmp.t_slice + t_offset ) % lt_orig;
850 
851  // Insert/write to disk
852  qdp_db.insert( tmp , peram[*key]);
853 
854  } // for key
855 
856  sniss2.stop();
857  QDPIO::cout << "Time to write perambulators for spin_src= " << spin_source << " time = "
858  << sniss2.getTimeInSeconds()
859  << " secs" << std::endl;
860  }
861 
862  } // this ends the lifetime of qdp_db.
863  } // for spin_src
864  } // for tt
865 
866  swatch.stop();
867  QDPIO::cout << "Propagators computed: time= "
868  << swatch.getTimeInSeconds()
869  << " secs" << std::endl;
870  }
871  catch (const std::string& e)
872  {
873  QDPIO::cout << name << ": caught exception around qprop: " << e << std::endl;
874  QDP_abort(1);
875  }
876 
877  push(xml_out,"Relaxation_Iterations");
878  write(xml_out, "ncg_had", ncg_had);
879  pop(xml_out);
880 
881  pop(xml_out); // prop_dist
882 
883  snoop.stop();
884  QDPIO::cout << name << ": total time = "
885  << snoop.getTimeInSeconds()
886  << " secs" << std::endl;
887 
888  QDPIO::cout << name << ": ran successfully" << std::endl;
889 
890  END_CODE();
891  }
892 
893  }
894 
895 } // namespace Chroma
896 
897 #endif
Inline measurement factory.
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.
Serializable value harness.
Definition: key_val_db.h:69
static T & Instance()
Definition: singleton.h:432
Builds time slice subsets.
const Set & getSet() const
The set to be used in sumMulti.
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.
void read(XMLReader &xml, const std::string &path, AsqtadFermActParams &param)
Read parameters.
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
MODS_t & eigen_source
Eigenvectors.
SUB_MOD_t sub_eigen
Where we store the sublattice versions.
Compute the propagator from distillation.
Key for propagator colorstd::vector sources.
Key for vanilla distillation propagator sources and solutions.
Key for propagator colorstd::vector matrix elements.
Key for time-sliced color eigenvectors.
Key and values for DB.
Make xml file writer.
int t
Definition: meslate.cc:37
int t_slice
Definition: meslate.cc:23
Nd
Definition: meslate.cc:74
Named object function std::map.
static bool registered
Local registration flag.
QDP::MapObjectDiskMultiple< KeyTimeSliceColorVec_t, TimeSliceIO< LatticeColorVectorF > > MODS_t
QDP::MapObjectDisk< KeyTimeSliceColorVec_t, TimeSliceIO< LatticeColorVectorF > > MOD_t
void read(XMLReader &xml, const std::string &path, Params::NamedObject_t &input)
Propagator input.
QDP::MapObjectMemory< KeyTimeSliceColorVec_t, SubLatticeColorVectorF > SUB_MOD_t
bool registerAll()
Register all the factories.
void write(XMLWriter &xml, const std::string &path, const Params::NamedObject_t &input)
Propagator output.
BinaryReturn< C1, C2, FnInnerProduct >::Type_t innerProduct(const QDPSubType< T1, C1 > &s1, const QDPType< T2, C2 > &s2)
bool registerAll()
Register all the factories.
Asqtad Staggered-Dirac operator.
Definition: klein_gord.cc:10
LatticeFermion tmp
Definition: mespbg5p_w.cc:36
push(xml_out,"Condensates")
int i
Definition: pbg5p_w.cc:55
Double mass
Definition: pbg5p_w.cc:54
pop(xml_out)
START_CODE()
Double zero
Definition: invbicg.cc:106
static QDP_ColorVector * in
multi1d< LatticeFermion > s(Ncb)
int count
Definition: octave.h:14
::std::string string
Definition: gtest.h:1979
Print out basic info about this program.
Fourier transform phase factor support.
Sparse matrix representation of spin matrices.
Distillation propagators.
Holds of vectors and weights.
Convenience for building time-slice subsets.
Contraction operators for two quarks.