CHROMA
inline_seqsource_w.cc
Go to the documentation of this file.
1 /*! \file
2  * \brief Inline construction of sequential sources
3  *
4  * Sequential source construction
5  */
6 
7 #include "handle.h"
10 #include "meas/glue/mesplq.h"
15 #include "util/ft/sftmom.h"
16 #include "util/info/proginfo.h"
17 #include "util/info/unique_id.h"
19 
20 namespace Chroma
21 {
22  //! Propagator input
23  void read(XMLReader& xml, const std::string& path, InlineSeqSourceEnv::Params::NamedObject_t& input)
24  {
25  XMLReader inputtop(xml, path);
26 
27  read(inputtop, "gauge_id", input.gauge_id);
28  read(inputtop, "prop_ids", input.prop_ids);
29  read(inputtop, "seqsource_id", input.seqsource_id);
30  }
31 
32  //! Propagator output
33  void write(XMLWriter& xml, const std::string& path, const InlineSeqSourceEnv::Params::NamedObject_t& input)
34  {
35  push(xml, path);
36 
37  write(xml, "gauge_id", input.gauge_id);
38  write(xml, "prop_ids", input.prop_ids);
39  write(xml, "seqsource_id", input.seqsource_id);
40 
41  pop(xml);
42  }
43 
44 
45  namespace InlineSeqSourceEnv
46  {
47  namespace
48  {
49  AbsInlineMeasurement* createMeasurement(XMLReader& xml_in,
50  const std::string& path)
51  {
52  return new InlineMeas(Params(xml_in, path));
53  }
54 
55  //! Local registration flag
56  bool registered = false;
57  }
58 
59  const std::string name = "SEQSOURCE";
60 
61  //! Register all the factories
62  bool registerAll()
63  {
64  bool success = true;
65  if (! registered)
66  {
69  success &= TheInlineMeasurementFactory::Instance().registerObject(name, createMeasurement);
70  registered = true;
71  }
72  return success;
73  }
74 
75 
76  // Param stuff
78 
79  Params::Params(XMLReader& xml_in, const std::string& path)
80  {
81  try
82  {
83  XMLReader paramtop(xml_in, path);
84 
85  if (paramtop.count("Frequency") == 1)
86  read(paramtop, "Frequency", frequency);
87  else
88  frequency = 1;
89 
90  // The parameters holds the version number
91  read(paramtop, "Param", param);
92 
93  // The parameters for smearing the sink
94  read(paramtop, "PropSink", sink_header);
95 
96  // Read in the forward_prop/seqsource info
97  read(paramtop, "NamedObject", named_obj);
98  }
99  catch(const std::string& e)
100  {
101  QDPIO::cerr << __func__ << ": Caught Exception reading XML: " << e << std::endl;
102  QDP_abort(1);
103  }
104  }
105 
106 
107  void
108  Params::writeXML(XMLWriter& xml_out, const std::string& path)
109  {
110  push(xml_out, path);
111 
112  write(xml_out, "Param", param);
113  write(xml_out, "PropSink", sink_header);
114  write(xml_out, "NamedObject", named_obj);
115 
116  pop(xml_out);
117  }
118 
119 
120  // Function call
121  void
122  InlineMeas::operator()(unsigned long update_no,
123  XMLWriter& xml_out)
124  {
125  START_CODE();
126 
127  StopWatch snoop;
128  snoop.reset();
129  snoop.start();
130 
131  // Test and grab a reference to the gauge field
132  XMLBufferWriter gauge_xml;
133  try
134  {
135  TheNamedObjMap::Instance().getData< multi1d<LatticeColorMatrix> >(params.named_obj.gauge_id);
136  TheNamedObjMap::Instance().get(params.named_obj.gauge_id).getRecordXML(gauge_xml);
137  }
138  catch( std::bad_cast )
139  {
140  QDPIO::cerr << name << ": caught dynamic cast error"
141  << std::endl;
142  QDP_abort(1);
143  }
144  catch (const std::string& e)
145  {
146  QDPIO::cerr << name << ": std::map call failed: " << e
147  << std::endl;
148  QDP_abort(1);
149  }
150  const multi1d<LatticeColorMatrix>& u =
151  TheNamedObjMap::Instance().getData< multi1d<LatticeColorMatrix> >(params.named_obj.gauge_id);
152 
153  push(xml_out, "seqsource");
154  write(xml_out, "update_no", update_no);
155 
156  QDPIO::cout << name << ": propagator sequential source constructor" << std::endl;
157  StopWatch swatch;
158 
159  proginfo(xml_out); // Print out basic program info
160 
161  // Write out the input
162  params.writeXML(xml_out, "Input");
163 
164  // Write out the config header
165  write(xml_out, "Config_info", gauge_xml);
166 
167  push(xml_out, "Output_version");
168  write(xml_out, "out_version", 1);
169  pop(xml_out);
170 
171  // Calculate some gauge invariant observables just for info.
172  MesPlq(xml_out, "Observables", u);
173 
174  // Sanity check
175  if (params.named_obj.prop_ids.size() == 0)
176  {
177  QDPIO::cerr << name << ": sanity error: " << std::endl;
178  QDP_abort(1);
179  }
180 
181  //
182  // Read the quark propagator and extract headers
183  //
184  multi1d<LatticePropagator> forward_props(params.named_obj.prop_ids.size());
185  multi1d<ForwardProp_t> forward_headers(params.named_obj.prop_ids.size());
186  push(xml_out, "Forward_prop_infos");
187  for(int loop=0; loop < params.named_obj.prop_ids.size(); ++loop)
188  {
189  push(xml_out, "elem");
190  try
191  {
192  // Snarf the data into a copy
193  forward_props[loop] =
194  TheNamedObjMap::Instance().getData<LatticePropagator>(params.named_obj.prop_ids[loop]);
195 
196  // Snarf the source info. This is will throw if the source_id is not there
197  XMLReader prop_file_xml, prop_record_xml;
198  TheNamedObjMap::Instance().get(params.named_obj.prop_ids[loop]).getFileXML(prop_file_xml);
199  TheNamedObjMap::Instance().get(params.named_obj.prop_ids[loop]).getRecordXML(prop_record_xml);
200 
201  // Try to invert this record XML into a ChromaProp struct
202  {
203  Propagator_t header;
204  read(prop_record_xml, "/Propagator", header);
205 
206  forward_headers[loop].prop_header = header.prop_header;
207  forward_headers[loop].source_header = header.source_header;
208  forward_headers[loop].gauge_header = header.gauge_header;
209  }
210 
211  // Save prop input
212  write(xml_out, "Propagator_info", prop_record_xml);
213  }
214  catch( std::bad_cast )
215  {
216  QDPIO::cerr << name << ": caught dynamic cast error"
217  << std::endl;
218  QDP_abort(1);
219  }
220  catch (const std::string& e)
221  {
222  QDPIO::cerr << name << ": std::map call failed: " << e
223  << std::endl;
224  QDP_abort(1);
225  }
226  pop(xml_out);
227  }
228  pop(xml_out);
229 
230  QDPIO::cout << "Forward propagator successfully read and parsed" << std::endl;
231 
232  // Derived from input prop
233  int j_decay = forward_headers[0].source_header.j_decay;
234 
235  // Initialize the slow Fourier transform phases
236  SftMom phases(0, true, j_decay);
237 
238  // Sanity check - write out the norm2 of the forward prop in the j_decay direction
239  // Use this for any possible verification
240  push(xml_out, "Forward_prop_correlators");
241  for(int loop=0; loop < params.named_obj.prop_ids.size(); ++loop)
242  {
243  multi1d<Double> forward_prop_corr = sumMulti(localNorm2(forward_props[loop]),
244  phases.getSet());
245 
246  push(xml_out, "elem");
247  write(xml_out, "forward_prop_corr", forward_prop_corr);
248  pop(xml_out);
249  }
250  pop(xml_out);
251 
252  // A sanity check
253  if (params.param.t_sink < 0 || params.param.t_sink >= QDP::Layout::lattSize()[j_decay])
254  {
255  QDPIO::cerr << "Sink time coordinate incorrect." << std::endl;
256  QDPIO::cerr << "t_sink = " << params.param.t_sink << std::endl;
257  QDP_abort(1);
258  }
259 
260 
261  //------------------ Start main body of calculations -----------------------------
262 
263  LatticePropagator quark_prop_src;
264 
265  try
266  {
267  // Sink smear the forward propagators
268  // NOTE: The smearing construction is pulled outside the loop
269  // for efficiency. However, I'm anticipating that we will
270  // have different smearings at the sink of the forward props.
271  // In that case, the loop needs to be in inverted.
272  std::istringstream xml_s(params.sink_header.sink.xml);
273  XMLReader sinktop(xml_s);
274  QDPIO::cout << "Sink = " << params.sink_header.sink.id << std::endl;
275 
277  sinkSmearing(ThePropSinkSmearingFactory::Instance().createObject(params.sink_header.sink.id,
278  sinktop,
280  u));
281 
282  // Do the sink smearing BEFORE the interpolating operator
283  for(int loop=0; loop < params.named_obj.prop_ids.size(); ++loop)
284  {
285  forward_headers[loop].sink_header = params.sink_header;
286  (*sinkSmearing)(forward_props[loop]);
287  }
288 
289 
290  //
291  // Construct the sequential source
292  //
293  QDPIO::cout << "Sequential source = " << params.param.seqsrc.xml << std::endl;
294 
295  std::istringstream xml_seq(params.param.seqsrc.xml);
296  XMLReader seqsrctop(xml_seq);
297  QDPIO::cout << "SeqSource = " << params.param.seqsrc.id << std::endl;
298 
300  hadSeqSource(TheWilsonHadronSeqSourceFactory::Instance().createObject(params.param.seqsrc.id,
301  seqsrctop,
303 
304  swatch.reset();
305  swatch.start();
306  quark_prop_src = (*hadSeqSource)(u, forward_headers, forward_props);
307 
308  swatch.stop();
309 
310  QDPIO::cout << "Hadron sequential source computed: time= "
311  << swatch.getTimeInSeconds()
312  << " secs" << std::endl;
313 
314 
315  // Do the sink smearing AFTER the interpolating operator
316  (*sinkSmearing)(quark_prop_src);
317 
318  }
319  catch(const std::string& e)
320  {
321  QDPIO::cerr << name << ": Caught Exception in sink: " << e << std::endl;
322  QDP_abort(1);
323  }
324 
325 
326  // Sanity check - write out the norm2 of the propagator source in the j_decay direction
327  // Use this for any possible verification
328  {
329  multi1d<Double> seqsource_corr = sumMulti(localNorm2(quark_prop_src),
330  phases.getSet());
331 
332  push(xml_out, "SeqSource_correlator");
333  write(xml_out, "seqsource_corr", seqsource_corr);
334  pop(xml_out);
335  }
336 
337 
338  /*
339  * Write the sequential source out to a named buffer
340  */
341  try
342  {
343  QDPIO::cout << "Attempt to store sequential source" << std::endl;
344 
345  XMLBufferWriter file_xml;
346  push(file_xml, "seqsource");
347  write(file_xml, "id", uniqueId()); // NOTE: new ID form
348  pop(file_xml);
349 
350  // Sequential source header
351  // Header composed of all forward prop headers
352  SequentialSource_t new_header;
353  new_header.sink_header = params.sink_header;
354  new_header.seqsource_header = params.param;
355  new_header.forward_props = forward_headers;
356  new_header.gauge_header = gauge_xml.printCurrentContext();
357 
358  XMLBufferWriter record_xml;
359  write(record_xml, "SequentialSource", new_header);
360 
361  // Store the seqsource
362  TheNamedObjMap::Instance().create<LatticePropagator>(params.named_obj.seqsource_id);
363  TheNamedObjMap::Instance().getData<LatticePropagator>(params.named_obj.seqsource_id) = quark_prop_src;
364  TheNamedObjMap::Instance().get(params.named_obj.seqsource_id).setFileXML(file_xml);
365  TheNamedObjMap::Instance().get(params.named_obj.seqsource_id).setRecordXML(record_xml);
366 
367  QDPIO::cout << "Sequential source successfully stored" << std::endl;
368  }
369  catch (std::bad_cast)
370  {
371  QDPIO::cerr << name << ": dynamic cast error"
372  << std::endl;
373  QDP_abort(1);
374  }
375  catch (const std::string& e)
376  {
377  QDPIO::cerr << name << ": error storing seqsource: " << e << std::endl;
378  QDP_abort(1);
379  }
380 
381  pop(xml_out); // seqsource
382 
383  snoop.stop();
384  QDPIO::cout << name << ": total time = "
385  << snoop.getTimeInSeconds()
386  << " secs" << std::endl;
387 
388  QDPIO::cout << name << ": ran successfully" << std::endl;
389 
390  END_CODE();
391  }
392 
393  }
394 
395 }
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.
Fourier transform phase factor support.
Definition: sftmom.h:35
const Set & getSet() const
The set to be used in sumMulti.
Definition: sftmom.h:57
static T & Instance()
Definition: singleton.h:432
void read(XMLReader &xml, const std::string &path, AsqtadFermActParams &param)
Read parameters.
void write(XMLWriter &xml, const std::string &path, const AsqtadFermActParams &param)
Writer parameters.
std::string uniqueId()
Return a unique id.
Definition: unique_id.cc:18
void proginfo(XMLWriter &xml)
Print out basic information about this program.
Definition: proginfo.cc:24
Class for counted reference semantics.
Inline construction of sequential sources.
int j_decay
Definition: meslate.cc:22
Named object function std::map.
static bool registered
Local registration flag.
bool registerAll()
Register all the factories.
bool registerAll()
Register all the factories.
bool registerAll()
Register all the factories.
multi1d< ForwardProp_t > & forward_headers
Asqtad Staggered-Dirac operator.
Definition: klein_gord.cc:10
static multi1d< LatticeColorMatrix > u
push(xml_out,"Condensates")
pop(xml_out)
void MesPlq(const multi1d< LatticeColorMatrixF3 > &u, multi2d< Double > &plane_plaq, Double &link)
Definition: mesplq.cc:83
START_CODE()
::std::string string
Definition: gtest.h:1979
Print out basic info about this program.
All sequential source constructors.
Factory for producing quark prop sinks.
Fourier transform phase factor support.
All make sink constructors.
Factory for producing quark prop sinks.
void writeXML(XMLWriter &xml_out, const std::string &path)
struct Chroma::InlineSeqSourceEnv::Params::NamedObject_t named_obj
Mega structure holding a full forward prop (not sink smeared)
Definition: qprop_io.h:111
ChromaProp_t prop_header
Definition: qprop_io.h:112
std::string gauge_header
Definition: qprop_io.h:114
PropSourceConst_t source_header
Definition: qprop_io.h:113
GroupXML_t seqsrc
Definition: qprop_io.h:93
Mega structure holding a full sequential source.
Definition: qprop_io.h:130
PropSinkSmear_t sink_header
Definition: qprop_io.h:131
multi1d< ForwardProp_t > forward_props
Definition: qprop_io.h:133
SeqSource_t seqsource_header
Definition: qprop_io.h:132
std::string gauge_header
Definition: qprop_io.h:134
Generate a unique id.