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"
39 #ifndef QDP_IS_QDPJIT_NO_NVPTX
41 #ifdef BUILD_JIT_CONTRACTION_KERNELS
42 #include "custom_kernels/custom_kernels.h"
52 multi1d<SubsetVectorWeight_t> readEigVals(
const std::string& meta)
54 std::istringstream xml_l(meta);
55 XMLReader eigtop(xml_l);
59 multi1d< multi1d<Real> > eigens;
62 read(eigtop, pat, eigens);
66 QDPIO::cerr << __func__ <<
": Caught Exception reading meta= XX" << meta <<
"XX with path= " << pat <<
" error= " << e << std::endl;
72 multi1d<SubsetVectorWeight_t> eigenvalues(eigens.size());
74 for(
int i=0;
i < eigens.size(); ++
i)
75 eigenvalues[
i].weights = eigens[
i];
84 namespace InlineMatElemDistillationEnv
89 XMLReader inputtop(xml, path);
112 XMLReader inputtop(xml, path);
124 if( inputtop.count(
"zero_colorvecs") == 1 ) {
128 QDPIO::cout <<
"zero_colorvecs found, *** timing mode activated ***\n";
133 if( inputtop.count(
"fuse_timeloop") == 1 ) {
137 QDPIO::cout <<
"fuse_timeloop found!\n";
162 XMLReader inputtop(xml, path);
200 namespace InlineMatElemDistillationEnv
204 typedef QDP::MapObjectDisk<KeyTimeSliceColorVec_t, TimeSliceIO<LatticeColorVectorF> >
MOD_t;
207 typedef QDP::MapObjectDiskMultiple<KeyTimeSliceColorVec_t, TimeSliceIO<LatticeColorVectorF> >
MODS_t;
210 typedef QDP::MapObjectMemory<KeyTimeSliceColorVec_t, SubLatticeColorVectorF>
SUB_MOD_t;
226 const SubLatticeColorVectorF& getVec(
int t_source,
int colorvec_src)
const;
248 const SubLatticeColorVectorF& SubEigenMap::getVec(
int t_source,
int colorvec_src)
const
256 QDPIO::cout << __func__ <<
": on t_source= " <<
t_source <<
" colorvec_src= " << colorvec_src << std::endl;
259 LatticeColorVectorF vec_srce;
263 TimeSliceIO<LatticeColorVectorF> time_slice_io(vec_srce,
t_source );
267 SubLatticeColorVectorF
tmp(getSet()[
t_source], vec_srce);
277 std::vector<bool> getActiveTSlices(
int t_source,
int Nt_forward,
int Nt_backward)
280 const int decay_dir =
Nd-1;
281 const int Lt = Layout::lattSize()[decay_dir];
283 std::vector<bool> active_t_slices(Lt);
284 for(
int t=0;
t < Lt; ++
t)
286 active_t_slices[
t] =
false;
290 for(
int dt=0;
dt < Nt_forward; ++
dt)
293 active_t_slices[
t % Lt] =
true;
297 for(
int dt=0;
dt < Nt_backward; ++
dt)
300 while (
t < 0) {
t += Lt;}
302 active_t_slices[
t % Lt] =
true;
305 return active_t_slices;
311 std::list<KeyPropElementalOperator_t> getSnkKeys(
int t_source,
int spin_source,
int Nt_forward,
int Nt_backward,
const std::string mass)
313 std::list<KeyPropElementalOperator_t> keys;
315 std::vector<bool> active_t_slices = getActiveTSlices(
t_source, Nt_forward, Nt_backward);
317 const int decay_dir =
Nd-1;
318 const int Lt = Layout::lattSize()[decay_dir];
320 for(
int spin_sink=0; spin_sink < Ns; ++spin_sink)
322 for(
int t=0;
t < Lt; ++
t)
324 if (! active_t_slices[
t]) {
continue;}
346 namespace InlineMatElemDistillationEnv
353 return new InlineMeas(Params(xml_in, path));
384 XMLReader paramtop(xml_in, path);
386 if (paramtop.count(
"Frequency") == 1)
398 if (paramtop.count(
"xml_file") != 0)
405 QDPIO::cerr << __func__ <<
": Caught Exception reading XML: " << e << std::endl;
423 push(xml_out,
"PropDistillation");
424 write(xml_out,
"update_no", update_no);
425 write(xml_out,
"xml_file", xml_file);
428 XMLFileWriter xml(xml_file);
429 func(update_no, xml);
433 func(update_no, xml_out);
449 push(xml_out,
"PropDistillation");
450 write(xml_out,
"update_no", update_no);
457 push(xml_out,
"Output_version");
458 write(xml_out,
"out_version", 1);
463 const int Lt = Layout::lattSize()[decay_dir];
468 if (decay_dir !=
Nd-1)
470 QDPIO::cerr <<
name <<
": TimeSliceIO only supports decay_dir= " <<
Nd-1 <<
"\n";
475 QDP::MapObjectDisk<KeyPropDistillation_t, TimeSliceIO<LatticeColorVectorF> > prop_obj;
476 prop_obj.setDebug(0);
478 QDPIO::cout <<
"Open source prop file" << std::endl;
490 QDPIO::cout <<
"Finished opening solution file" << std::endl;
496 QDPIO::cout <<
"Snarf the source from a std::map object disk file" << std::endl;
512 QDPIO::cout <<
"Get user data" << std::endl;
521 catch (std::bad_cast) {
522 QDPIO::cerr <<
name <<
": caught dynamic cast error" << std::endl;
526 QDPIO::cerr <<
name <<
": error extracting source_header: " << e << std::endl;
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;
536 QDPIO::cout <<
"Source successfully read and parsed" << std::endl;
543 QDPIO::cerr <<
name <<
": number of available eigenvectors is too small\n";
548 QDPIO::cout <<
"Number of vecs available is large enough" << std::endl;
551 QDPIO::cout <<
"Initialize sub-lattice std::map" << std::endl;
553 QDPIO::cout <<
"Finished initializing sub-lattice std::map" << std::endl;
566 XMLBufferWriter file_xml;
567 push(file_xml,
"DBMetaData");
570 auto tmp = QDP::Layout::lattSize();
578 write(file_xml,
"Weights", readEigVals(eigen_meta_data));
582 qdp_db.setMaxUserInfoLen(file_str.size());
586 qdp_db.insertUserdata(file_str);
593 QDPIO::cout <<
"Finished opening peram file" << std::endl;
606 QDPIO::cout <<
"Try the various factories" << std::endl;
610 #ifdef BUILD_JIT_CONTRACTION_KERNELS
612 QDPIO::cout <<
"Using JIT contraction kernels (fused timeloop)\n";
614 QDPIO::cout <<
"Using JIT contraction kernels (non-fused timeloop)\n";
626 for(
int tt=0; tt < t_sources.size(); ++tt)
629 QDPIO::cout <<
"t_source = " <<
t_source << std::endl;
632 QDPIO::cout <<
"new t_source in current lattice = " <<
t_source << std::endl;
635 for(
int spin_source=0; spin_source < Ns; ++spin_source)
637 QDPIO::cout <<
"spin_source = " << spin_source << std::endl;
640 std::list<KeyPropElementalOperator_t> snk_keys(getSnkKeys(
t_source,
648 QDP::MapObjectMemory<KeyPropElementalOperator_t, ValPropElementalOperator_t> peram;
651 for(std::list<KeyPropElementalOperator_t>::const_iterator key = snk_keys.begin();
652 key != snk_keys.end();
657 peram.insert(*key,
tmp);
659 peram[*key].mat.resize(num_vecs,num_vecs);
660 peram[*key].mat =
zero;
662 QDPIO::cout <<
"peram initialized! " << std::endl;
667 for(
int colorvec_src=0; colorvec_src < num_vecs; ++colorvec_src)
676 QDPIO::cout <<
"Do spin_source= " << spin_source <<
" colorvec_src= " << colorvec_src << std::endl;
683 multi1d<LatticeColorVector> ferm_out(Ns);
688 for(
int spin_sink=0; spin_sink < Ns; ++spin_sink)
690 for (
int t = 0 ;
t < Lt ; ++
t )
701 LatticeColorVectorF
tmp;
703 TimeSliceIO<LatticeColorVectorF> time_slice_io(
tmp,
t);
705 prop_obj.get( key , time_slice_io );
707 ferm_out(spin_sink)[ sub_eigen_map.getSet()[
t] ] =
tmp;
714 for (
int s = 0 ;
s < Ns ; ++
s )
716 zero_rep( ferm_out[
s] );
722 QDPIO::cout <<
"Time to read prop for spin_source= " << spin_source <<
" colorvec_src= " << colorvec_src <<
" time = "
723 << snarss1.getTimeInSeconds()
724 <<
" secs" << std::endl;
729 #ifndef BUILD_JIT_CONTRACTION_KERNELS
733 for(std::list<KeyPropElementalOperator_t>::const_iterator key= snk_keys.begin();
734 key != snk_keys.end();
737 if (key->t_slice !=
t_slice) {
continue;}
740 for(
int colorvec_sink=0; colorvec_sink < num_vecs; ++colorvec_sink)
742 peram[*key].mat(colorvec_sink,colorvec_src) =
innerProduct(sub_eigen_map.getVec(
t_slice, colorvec_sink),
743 ferm_out(key->spin_snk));
753 for(
int spin_snk = 0; spin_snk < Ns; ++spin_snk)
756 for(std::list<KeyPropElementalOperator_t>::const_iterator key= snk_keys.begin(); key != snk_keys.end(); ++key)
758 if (key->spin_snk != spin_snk) {
continue;}
764 multi1d<SubLatticeColorVectorF*> vec_ptr(
count );
765 multi1d<ComplexD*> contr_ptr(
count );
770 for(std::list<KeyPropElementalOperator_t>::const_iterator key= snk_keys.begin();
771 key != snk_keys.end();
774 if (key->spin_snk != spin_snk) {
continue;}
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 );
789 multi_innerProduct( contr_ptr , vec_ptr , ferm_out(spin_snk) );
801 for(std::list<KeyPropElementalOperator_t>::const_iterator key= snk_keys.begin();
802 key != snk_keys.end();
805 if (key->t_slice !=
t_slice) {
continue;}
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 );
819 multi_innerProduct( contr_ptr , vec_ptr , ferm_out(key->spin_snk) );
826 QDPIO::cout <<
"Time to compute and assemble peram for spin_source= " << spin_source <<
" colorvec_src= " << colorvec_src <<
" time = "
827 << sniss1.getTimeInSeconds()
829 <<
" (for contraction: " << sniss1.getTimeInSeconds() - snarss1.getTimeInSeconds() <<
")"
837 QDPIO::cout <<
"Write perambulator for spin_source= " << spin_source <<
" to disk" << std::endl;
843 for(std::list<KeyPropElementalOperator_t>::const_iterator key= snk_keys.begin();
844 key != snk_keys.end();
852 qdp_db.insert(
tmp , peram[*key]);
857 QDPIO::cout <<
"Time to write perambulators for spin_src= " << spin_source <<
" time = "
858 << sniss2.getTimeInSeconds()
859 <<
" secs" << std::endl;
867 QDPIO::cout <<
"Propagators computed: time= "
868 << swatch.getTimeInSeconds()
869 <<
" secs" << std::endl;
873 QDPIO::cout <<
name <<
": caught exception around qprop: " << e << std::endl;
877 push(xml_out,
"Relaxation_Iterations");
878 write(xml_out,
"ncg_had", ncg_had);
884 QDPIO::cout <<
name <<
": total time = "
885 << snoop.getTimeInSeconds()
886 <<
" secs" << std::endl;
888 QDPIO::cout <<
name <<
": ran successfully" << std::endl;
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.
Builds time slice subsets.
const Set & getSet() const
The set to be used in sumMulti.
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 ¶m)
Read parameters.
void proginfo(XMLWriter &xml)
Print out basic information about this program.
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.
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.
push(xml_out,"Condensates")
static QDP_ColorVector * in
multi1d< LatticeFermion > s(Ncb)
Print out basic info about this program.
Fourier transform phase factor support.
Sparse matrix representation of spin matrices.
std::vector< std::string > colorvec_files
Distillation propagators.
Holds of vectors and weights.
Convenience for building time-slice subsets.
Contraction operators for two quarks.