Bridge++  Ver. 2.0.2
fprop_alt_Standard_lex-tmpl.h
Go to the documentation of this file.
1 
10 template<typename AFIELD>
11 const std::string Fprop_alt_Standard_lex<AFIELD>::class_name
12  = "Fprop_alt_Standard_lex";
13 //====================================================================
14 template<typename AFIELD>
15 void Fprop_alt_Standard_lex<AFIELD>::init(const Parameters& params_fopr,
16  const Parameters& params_solver)
17 { // this constructor assumes that the factories are available.
18  ThreadManager::assert_single_thread(class_name);
19 
20  // vout.general(m_vl, "\n");
21  vout.general(m_vl, "%s: being setup (without link smearing).\n",
22  class_name.c_str());
23 
24  // typedef AFopr<AField<real_t> > AltFopr;
25  // typedef ASolver<AField<real_t> > AltSolver;
26  typedef AFopr<AFIELD> AltFopr;
27  typedef ASolver<AFIELD> AltSolver;
28 
29  string fopr_type = params_fopr.get_string("fermion_type");
30  m_fopr = AltFopr::New(fopr_type, params_fopr);
31 
32  m_dr_smear = 0;
33 
34  m_kernel = 0;
35 
36  string solver_type = params_solver.get_string("solver_type");
37  m_solver = AltSolver::New(solver_type, m_fopr);
38  m_solver->set_parameters(params_solver);
39 
40  reset_performance();
41 
42  vout.general(m_vl, "%s: setup finished.\n", class_name.c_str());
43 }
44 
45 
46 //====================================================================
47 template<typename AFIELD>
48 void Fprop_alt_Standard_lex<AFIELD>::init(const Parameters& params_fopr,
49  const Parameters& params_solver,
50  Director_Smear *dr_smear)
51 { // this constructor assumes that the factories are available.
52  // vout.general(m_vl, "\n");
53  vout.general(m_vl, "%s: being setup (with link smearing).\n",
54  class_name.c_str());
55 
56  typedef AFopr<AFIELD> AltFopr;
57  typedef ASolver<AFIELD> AltSolver;
58 
59  m_dr_smear = dr_smear;
60 
61  string fopr_type = params_fopr.get_string("fermion_type");
62  m_kernel = AltFopr::New(fopr_type, params_fopr);
63 
64  // m_fopr = AltFopr::New("Smeared", m_kernel, m_dr_smear);
65  m_fopr = new AFopr_Smeared<AFIELD>(m_kernel, m_dr_smear);
66 
67  string solver_type = params_solver.get_string("solver_type");
68  m_solver = AltSolver::New(solver_type, m_fopr);
69  m_solver->set_parameters(params_solver);
70 
71  reset_performance();
72 
73  vout.general(m_vl, "%s: setup finished.\n", class_name.c_str());
74 }
75 
76 
77 //====================================================================
78 template<typename AFIELD>
79 void Fprop_alt_Standard_lex<AFIELD>::tidyup()
80 {
81  delete m_solver;
82  delete m_fopr;
83  if (m_kernel != 0) delete m_kernel;
84 }
85 
86 
87 //====================================================================
88 template<typename AFIELD>
89 void Fprop_alt_Standard_lex<AFIELD>::set_config(Field *U)
90 {
91  m_fopr->set_config(U);
92 }
93 
94 
95 //====================================================================
96 template<typename AFIELD>
97 void Fprop_alt_Standard_lex<AFIELD>::invert(Field& xq, const Field& b,
98  int& nconv, double& diff)
99 {
100  vout.paranoiac(m_vl, "%s: invert is called.\n", class_name.c_str());
101  vout.paranoiac(m_vl, "mode = %s.\n", m_mode.c_str());
102 
103  if (m_mode == "D") {
104  invert_D(xq, b, nconv, diff);
105  } else if (m_mode == "DdagD") {
106  invert_DdagD(xq, b, nconv, diff);
107  } else if (m_mode == "D_prec") {
108  invert_D_prec(xq, b, nconv, diff);
109  } else if (m_mode == "DdagD_prec") {
110  invert_DdagD_prec(xq, b, nconv, diff);
111  } else {
112  vout.crucial(m_vl, "%s: unsupported mode: %s\n",
113  class_name.c_str(), m_mode.c_str());
114  exit(EXIT_FAILURE);
115  }
116 }
117 
118 
119 //====================================================================
120 template<typename AFIELD>
121 void Fprop_alt_Standard_lex<AFIELD>::invert_D(Field& xq, const Field& b,
122  int& nconv, double& diff)
123 {
124  m_timer.reset();
125  m_timer.start();
126 
127  int nin = m_fopr->field_nin();
128  int nvol = m_fopr->field_nvol();
129  int nex = m_fopr->field_nex();
130 
131  vout.paranoiac(m_vl, "nin = %d nvol = %d nex = %d\n", nin, nvol, nex);
132 
133  AFIELD axq(nin, nvol, nex);
134  AFIELD abq(nin, nvol, nex);
135 
136  AIndex_lex<real_t, AFIELD::IMPL> index_alt;
137 
138 #pragma omp parallel
139  {
140  if (m_fopr->needs_convert()) {
141  vout.detailed(m_vl, "convert required.\n");
142  m_fopr->convert(abq, b);
143  } else {
144  vout.detailed(m_vl, "convert not required.\n");
145  convert(index_alt, abq, b);
146  }
147  }
148 
149  vout.detailed(m_vl, "%s: convert finished.\n", class_name.c_str());
150 
151  m_fopr->set_mode("D");
152 
153  real_t diff2;
154 
155 #pragma omp parallel
156  {
157  m_solver->solve(axq, abq, nconv, diff2);
158  }
159  diff = double(diff2);
160 
161 #pragma omp parallel
162  {
163  if (m_fopr->needs_convert()) {
164  m_fopr->reverse(xq, axq);
165  } else {
166  reverse(index_alt, xq, axq);
167  }
168  }
169 
170  m_timer.stop();
171  m_elapsed_time += m_timer.elapsed_sec();
172  m_flop_count += m_solver->flop_count();
173 }
174 
175 
176 //====================================================================
177 template<typename AFIELD>
178 void Fprop_alt_Standard_lex<AFIELD>::invert_DdagD(Field& xq, const Field& b,
179  int& nconv, double& diff)
180 {
181  m_timer.reset();
182  m_timer.start();
183 
184  int nin = m_fopr->field_nin();
185  int nvol = m_fopr->field_nvol();
186  int nex = m_fopr->field_nex();
187 
188  AFIELD axq(nin, nvol, nex);
189  AFIELD abq(nin, nvol, nex);
190 
191  AIndex_lex<real_t, AFIELD::IMPL> index_alt;
192 
193 #pragma omp parallel
194  {
195  if (m_fopr->needs_convert()) {
196  m_fopr->convert(abq, b);
197  } else {
198  convert(index_alt, abq, b);
199  }
200  }
201 
202  real_t diff2;
203 
204  m_fopr->set_mode("DdagD");
205 
206 #pragma omp parallel
207  {
208  m_solver->solve(axq, abq, nconv, diff2);
209  }
210  diff = double(diff2);
211 
212 #pragma omp parallel
213  {
214  if (m_fopr->needs_convert()) {
215  m_fopr->reverse(xq, axq);
216  } else {
217  reverse(index_alt, xq, axq);
218  }
219  }
220 
221  m_timer.stop();
222  m_elapsed_time += m_timer.elapsed_sec();
223  m_flop_count += m_solver->flop_count();
224 }
225 
226 
227 //====================================================================
228 template<typename AFIELD>
229 void Fprop_alt_Standard_lex<AFIELD>::invert_D_prec(Field& xq,
230  const Field& b,
231  int& nconv, double& diff)
232 {
233  m_timer.reset();
234  m_timer.start();
235 
236  int nin = m_fopr->field_nin();
237  int nvol = m_fopr->field_nvol();
238  int nex = m_fopr->field_nex();
239 
240  AFIELD axq(nin, nvol, nex);
241  AFIELD abq(nin, nvol, nex);
242 
243  AIndex_lex<real_t, AFIELD::IMPL> index_alt;
244 
245  if (m_fopr->needs_convert()) {
246  m_fopr->convert(abq, b);
247  } else {
248  convert(index_alt, abq, b);
249  }
250 
251  m_fopr->set_mode("D_prec");
252 
253  real_t diff2;
254 
255 #pragma omp parallel
256  {
257  m_solver->solve(axq, abq, nconv, diff2);
258  }
259  m_flop_count += m_solver->flop_count();
260 
261  m_fopr->mult(abq, axq, "Prec");
262 
263  if (m_fopr->needs_convert()) {
264  m_fopr->reverse(xq, abq);
265  } else {
266  reverse(index_alt, xq, abq);
267  }
268 
269  diff = double(diff2);
270 
271  m_timer.stop();
272  m_elapsed_time += m_timer.elapsed_sec();
273  // 1 solve + 2 mult; flop for the solve has been already accumulated
274  m_flop_count += 2 * m_fopr->flop_count();
275 }
276 
277 
278 //====================================================================
279 template<typename AFIELD>
280 void Fprop_alt_Standard_lex<AFIELD>::invert_DdagD_prec(Field& xq,
281  const Field& b,
282  int& nconv, double& diff)
283 {
284  vout.paranoiac(m_vl, "invert_DdagD_prec is called.\n");
285 
286  m_timer.reset();
287  m_timer.start();
288 
289  int nin = m_fopr->field_nin();
290  int nvol = m_fopr->field_nvol();
291  int nex = m_fopr->field_nex();
292 
293  AFIELD axq(nin, nvol, nex);
294  AFIELD abq(nin, nvol, nex);
295 
296  AIndex_lex<real_t, AFIELD::IMPL> index_alt;
297 
298  if (m_fopr->needs_convert()) {
299  m_fopr->convert(axq, b);
300  } else {
301  convert(index_alt, axq, b);
302  }
303 
304  m_fopr->mult(abq, axq, "Precdag");
305 
306  m_fopr->set_mode("DdagD_prec");
307 
308  real_t diff2;
309 
310 #pragma omp parallel
311  {
312  m_solver->solve(axq, abq, nconv, diff2);
313  }
314  m_flop_count += m_solver->flop_count();
315 
316  m_fopr->mult(abq, axq, "Prec");
317 
318  if (m_fopr->needs_convert()) {
319  m_fopr->reverse(xq, abq);
320  } else {
321  reverse(index_alt, xq, abq);
322  }
323 
324  diff = double(diff2);
325 
326  m_timer.stop();
327  m_elapsed_time += m_timer.elapsed_sec();
328  // 1 solve + 2 mult; flop for the solve has been already accumulated
329  m_flop_count += 2 * m_fopr->flop_count();
330 }
331 
332 
333 //====================================================================
334 template<typename AFIELD>
335 void Fprop_alt_Standard_lex<AFIELD>::invert(AFIELD& xq,
336  const AFIELD& b,
337  int& nconv, double& diff)
338 {
339  vout.paranoiac(m_vl, "%s: invert is called.\n", class_name.c_str());
340  vout.paranoiac(m_vl, "mode = %s.\n", m_mode.c_str());
341 
342  if (m_mode == "D") {
343  invert_D(xq, b, nconv, diff);
344  } else if (m_mode == "DdagD") {
345  invert_DdagD(xq, b, nconv, diff);
346  } else if (m_mode == "D_prec") {
347  invert_D_prec(xq, b, nconv, diff);
348  } else if (m_mode == "DdagD_prec") {
349  invert_DdagD_prec(xq, b, nconv, diff);
350  } else {
351  vout.crucial(m_vl, "%s: unsupported mode: %s\n",
352  class_name.c_str(), m_mode.c_str());
353  exit(EXIT_FAILURE);
354  }
355 }
356 
357 
358 //====================================================================
359 template<typename AFIELD>
360 void Fprop_alt_Standard_lex<AFIELD>::invert_D(AFIELD& axq,
361  const AFIELD& abq,
362  int& nconv, double& diff)
363 {
364  m_timer.reset();
365  m_timer.start();
366 
367  m_fopr->set_mode("D");
368 
369  real_t diff2;
370 
371 #pragma omp parallel
372  {
373  m_solver->solve(axq, abq, nconv, diff2);
374  }
375  diff = double(diff2);
376 
377  m_timer.stop();
378  m_elapsed_time += m_timer.elapsed_sec();
379  m_flop_count += m_solver->flop_count();
380 }
381 
382 
383 //====================================================================
384 template<typename AFIELD>
385 void Fprop_alt_Standard_lex<AFIELD>::invert_DdagD(AFIELD& axq,
386  const AFIELD& abq,
387  int& nconv, double& diff)
388 {
389  m_timer.reset();
390  m_timer.start();
391 
392  real_t diff2;
393 
394  m_fopr->set_mode("DdagD");
395 
396 #pragma omp parallel
397  {
398  m_solver->solve(axq, abq, nconv, diff2);
399  }
400  diff = double(diff2);
401 
402  m_timer.stop();
403  m_elapsed_time += m_timer.elapsed_sec();
404  m_flop_count += m_solver->flop_count();
405 }
406 
407 
408 //====================================================================
409 template<typename AFIELD>
410 void Fprop_alt_Standard_lex<AFIELD>::invert_D_prec(AFIELD& axq,
411  const AFIELD& abq,
412  int& nconv, double& diff)
413 {
414  m_timer.reset();
415  m_timer.start();
416 
417  int nin = m_fopr->field_nin();
418  int nvol = m_fopr->field_nvol();
419  int nex = m_fopr->field_nex();
420  AFIELD atmp(nin, nvol, nex);
421 
422  m_fopr->set_mode("D_prec");
423 
424  real_t diff2;
425 
426 #pragma omp parallel
427  {
428  m_solver->solve(atmp, abq, nconv, diff2);
429  }
430  m_flop_count += m_solver->flop_count();
431  m_fopr->mult(axq, atmp, "Prec");
432 
433  diff = double(diff2);
434 
435  m_timer.stop();
436  m_elapsed_time += m_timer.elapsed_sec();
437  // 1 solve + 2 mult; flop for the solve has been already accumulated
438  m_flop_count += 2 * m_fopr->flop_count();
439 }
440 
441 
442 //====================================================================
443 template<typename AFIELD>
444 void Fprop_alt_Standard_lex<AFIELD>::invert_DdagD_prec(AFIELD& axq,
445  const AFIELD& abq,
446  int& nconv, double& diff)
447 {
448  vout.paranoiac(m_vl, "invert_DdagD_prec is called.\n");
449 
450  m_timer.reset();
451  m_timer.start();
452 
453  int nin = m_fopr->field_nin();
454  int nvol = m_fopr->field_nvol();
455  int nex = m_fopr->field_nex();
456  AFIELD atmp(nin, nvol, nex);
457 
458  m_fopr->mult(axq, abq, "Precdag");
459 
460  m_fopr->set_mode("DdagD_prec");
461 
462  real_t diff2;
463 
464 #pragma omp parallel
465  {
466  m_solver->solve(atmp, axq, nconv, diff2);
467  }
468  m_flop_count += m_solver->flop_count();
469 
470  m_fopr->mult(axq, atmp, "Prec");
471 
472  diff = double(diff2);
473 
474  m_timer.stop();
475  m_elapsed_time += m_timer.elapsed_sec();
476  // 1 solve + 2 mult; flop for the solve has been already accumulated
477  m_flop_count += 2 * m_fopr->flop_count();
478 }
479 
480 
481 //====================================================================
482 template<typename AFIELD>
483 double Fprop_alt_Standard_lex<AFIELD>::flop_count()
484 {
485  return m_solver->flop_count();
486 }
487 
488 
489 //====================================================================
490 template<typename AFIELD>
491 void Fprop_alt_Standard_lex<AFIELD>::reset_performance()
492 {
493  m_flop_count = 0.0;
494  m_elapsed_time = 0.0;
495 }
496 
497 
498 //====================================================================
499 template<typename AFIELD>
500 void Fprop_alt_Standard_lex<AFIELD>::get_performance(double& flop_count,
501  double& elapsed_time)
502 {
503  flop_count = m_flop_count;
504  elapsed_time = m_elapsed_time;
505 }
506 
507 
508 //====================================================================
509 template<typename AFIELD>
510 void Fprop_alt_Standard_lex<AFIELD>::report_performance()
511 {
512  double flops = m_flop_count / m_elapsed_time;
513  double gflops = flops * 1.0e-9;
514 
515  // vout.general(m_vl, "\n");
516  vout.general(m_vl, "%s: solver performance:\n", class_name.c_str());
517  vout.general(m_vl, " Elapsed time = %14.6f sec\n", m_elapsed_time);
518  vout.general(m_vl, " Flop(total) = %18.0f\n", m_flop_count);
519  vout.general(m_vl, " Performance = %11.3f GFlops\n", gflops);
520 }
521 
522 
523 //====================================================================
524 template<typename AFIELD>
525 void Fprop_alt_Standard_lex<AFIELD>::mult_performance(
526  const std::string mode,
527  const int Nrepeat)
528 {
529  ThreadManager::assert_single_thread(class_name);
530 
531  int nin = m_fopr->field_nin();
532  int nvol = m_fopr->field_nvol();
533  int nex = m_fopr->field_nex();
534 
535  AFIELD axq(nin, nvol, nex), abq(nin, nvol, nex);
536  abq.set(0.0);
537  abq.set(0, 1.0);
538 
539  unique_ptr<Timer> timer(new Timer);
540 
541  std::string mode_prev = m_fopr->get_mode();
542 
543  m_fopr->set_mode(mode);
544 
545  timer->start();
546 
547 #pragma omp parallel
548  {
549  for (int i = 0; i < Nrepeat; ++i) {
550  m_fopr->mult(axq, abq);
551  m_fopr->mult(abq, axq);
552  }
553  }
554 
555  timer->stop();
556 
557  double flop_fopr = m_fopr->flop_count();
558  double flop_total = flop_fopr * double(2 * Nrepeat);
559 
560  double elapsed_time = timer->elapsed_sec();
561  double flops = flop_total / elapsed_time;
562  double gflops = flops * 1.0e-9;
563 
564  vout.general(m_vl, "\n");
565  vout.general(m_vl, "%s: mult performance:\n", class_name.c_str());
566  vout.general(m_vl, " mult mode = %s\n", mode.c_str());
567  vout.general(m_vl, " Elapsed time = %14.6f sec\n", elapsed_time);
568  vout.general(m_vl, " Flop(Fopr) = %18.0f\n", flop_fopr);
569  vout.general(m_vl, " Flop(total) = %18.0f\n", flop_total);
570  vout.general(m_vl, " Performance = %11.3f GFlops\n", gflops);
571 
572  m_fopr->set_mode(mode_prev);
573 }
574 
575 
576 //============================================================END=====
Fprop_alt_Standard_lex::invert_D
void invert_D(Field &, const Field &, int &, double &)
Definition: fprop_alt_Standard_lex-tmpl.h:121
Fprop_alt_Standard_lex::init
void init(const Parameters &params_fopr, const Parameters &params_solver)
Definition: fprop_alt_Standard_lex-tmpl.h:15
AFopr
Definition: afopr.h:48
ASolver
Definition: asolver.h:23
Fprop_alt_Standard_lex::tidyup
void tidyup()
Definition: fprop_alt_Standard_lex-tmpl.h:79
Field::set
void set(const int jin, const int site, const int jex, double v)
Definition: field.h:175
Parameters
Class for parameters.
Definition: parameters.h:46
AIndex_lex
Definition: aindex_lex_base.h:17
convert
void convert(INDEX &index, AFIELD &v, const Field &w)
Definition: afield-inc.h:129
Bridge::BridgeIO::detailed
void detailed(const char *format,...)
Definition: bridgeIO.cpp:219
real_t
double real_t
Definition: bridgeQXS_Clover_coarse_double.cpp:16
Fprop_alt_Standard_lex::mult_performance
void mult_performance(const std::string mode, const int Nrepeat)
Definition: fprop_alt_Standard_lex-tmpl.h:525
Timer
Definition: timer.h:31
Fprop_alt_Standard_lex::flop_count
double flop_count()
Definition: fprop_alt_Standard_lex-tmpl.h:483
Bridge::BridgeIO::paranoiac
void paranoiac(const char *format,...)
Definition: bridgeIO.cpp:238
Fprop_alt_Standard_lex
Get quark propagator for Fopr with lexical site index: alternative version.
Definition: fprop_alt_Standard_lex.h:30
Fprop_alt_Standard_lex::invert_DdagD_prec
void invert_DdagD_prec(Field &, const Field &, int &, double &)
Definition: fprop_alt_Standard_lex-tmpl.h:280
AFopr_Smeared
smeared fermion operator: alternative version.
Definition: afopr_Smeared.h:41
Fprop_alt_Standard_lex::report_performance
void report_performance()
Definition: fprop_alt_Standard_lex-tmpl.h:510
Fprop_alt::real_t
AFIELD::real_t real_t
Definition: fprop_alt.h:32
Fprop_alt_Standard_lex::get_performance
void get_performance(double &flop_count, double &elapsed_time)
Definition: fprop_alt_Standard_lex-tmpl.h:500
Director_Smear
Manager of smeared configurations.
Definition: director_Smear.h:39
Fprop_alt_Standard_lex::invert_DdagD
void invert_DdagD(Field &, const Field &, int &, double &)
Definition: fprop_alt_Standard_lex-tmpl.h:178
Fprop_alt_Standard_lex::invert_D_prec
void invert_D_prec(Field &, const Field &, int &, double &)
Definition: fprop_alt_Standard_lex-tmpl.h:229
Fprop_alt_Standard_lex::invert
void invert(Field &, const Field &, int &, double &)
invert accordingly to the mode. [22 Sep 2018 H.Matsufuru]
Definition: fprop_alt_Standard_lex-tmpl.h:97
reverse
void reverse(INDEX &index, Field &v, const AFIELD &w)
Definition: afield-inc.h:159
Fprop_alt_Standard_lex::reset_performance
void reset_performance()
Definition: fprop_alt_Standard_lex-tmpl.h:491
Fprop_alt_Standard_lex::set_config
void set_config(Field *)
Definition: fprop_alt_Standard_lex-tmpl.h:89
Parameters::get_string
string get_string(const string &key) const
Definition: parameters.cpp:221
Bridge::BridgeIO::crucial
void crucial(const char *format,...)
Definition: bridgeIO.cpp:180
Field
Container of Field-type object.
Definition: field.h:46
Bridge::BridgeIO::general
void general(const char *format,...)
Definition: bridgeIO.cpp:200
ThreadManager::assert_single_thread
static void assert_single_thread(const std::string &class_name)
assert currently running on single thread.
Definition: threadManager.cpp:372
Bridge::vout
BridgeIO vout
Definition: bridgeIO.cpp:512