Bridge++  Ver.2.1.3
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.detailed(m_vl, "%s::invert_D() is called\n", class_name.c_str());
132 
133  vout.paranoiac(m_vl, "nin = %d nvol = %d nex = %d\n", nin, nvol, nex);
134 
135  AFIELD axq(nin, nvol, nex);
136  AFIELD abq(nin, nvol, nex);
137 
138  AIndex_lex<real_t, AFIELD::IMPL> index_alt;
139 
140 #pragma omp parallel
141  {
142  if (m_fopr->needs_convert()) {
143  vout.detailed(m_vl, "convert required.\n");
144  m_fopr->convert(abq, b);
145  } else {
146  vout.detailed(m_vl, "convert not required.\n");
147  convert(index_alt, abq, b);
148  }
149  }
150 
151  vout.detailed(m_vl, "%s: convert finished.\n", class_name.c_str());
152 
153  m_fopr->set_mode("D");
154 
155  real_t diff2;
156 
157 #pragma omp parallel
158  {
159  m_solver->solve(axq, abq, nconv, diff2);
160  }
161  diff = double(diff2);
162 
163 #pragma omp parallel
164  {
165  if (m_fopr->needs_convert()) {
166  m_fopr->reverse(xq, axq);
167  } else {
168  reverse(index_alt, xq, axq);
169  }
170  }
171 
172  m_timer.stop();
173  double elapsed_sec = m_timer.elapsed_sec();
174  double flop_count = m_solver->flop_count();
175 
176  m_elapsed_time += elapsed_sec;
177  m_flop_count += flop_count;
178 
179  vout.general(m_vl, "%s::invert_D():\n", class_name.c_str());
180  vout.general(m_vl, " elapase time = %e sec\n", elapsed_sec);
181  vout.general(m_vl, " flop count = %e\n", flop_count);
182  double gflops = 1.e-9 * flop_count/elapsed_sec;
183  vout.general(m_vl, " performance = %e GFlops\n", gflops);
184 
185 }
186 
187 
188 //====================================================================
189 template<typename AFIELD>
190 void Fprop_alt_Standard_lex<AFIELD>::invert_DdagD(Field& xq, const Field& b,
191  int& nconv, double& diff)
192 {
193  m_timer.reset();
194  m_timer.start();
195 
196  vout.detailed(m_vl, "%s::invert_DdagD() is called\n", class_name.c_str());
197 
198  int nin = m_fopr->field_nin();
199  int nvol = m_fopr->field_nvol();
200  int nex = m_fopr->field_nex();
201 
202  AFIELD axq(nin, nvol, nex);
203  AFIELD abq(nin, nvol, nex);
204 
205  AIndex_lex<real_t, AFIELD::IMPL> index_alt;
206 
207 #pragma omp parallel
208  {
209  if (m_fopr->needs_convert()) {
210  m_fopr->convert(abq, b);
211  } else {
212  convert(index_alt, abq, b);
213  }
214  }
215 
216  real_t diff2;
217 
218  m_fopr->set_mode("DdagD");
219 
220 #pragma omp parallel
221  {
222  m_solver->solve(axq, abq, nconv, diff2);
223  }
224  diff = double(diff2);
225 
226 #pragma omp parallel
227  {
228  if (m_fopr->needs_convert()) {
229  m_fopr->reverse(xq, axq);
230  } else {
231  reverse(index_alt, xq, axq);
232  }
233  }
234 
235  m_timer.stop();
236  double elapsed_sec = m_timer.elapsed_sec();
237  double flop_count = m_solver->flop_count();
238 
239  m_elapsed_time += elapsed_sec;
240  m_flop_count += flop_count;
241 
242  vout.general(m_vl, "%s::invert_DdagD():\n", class_name.c_str());
243  vout.general(m_vl, " elapase time = %e sec\n", elapsed_sec);
244  vout.general(m_vl, " flop count = %e\n", flop_count);
245  double gflops = 1.e-9 * flop_count/elapsed_sec;
246  vout.general(m_vl, " performance = %e GFlops\n", gflops);
247 
248 }
249 
250 
251 //====================================================================
252 template<typename AFIELD>
253 void Fprop_alt_Standard_lex<AFIELD>::invert_D_prec(Field& xq,
254  const Field& b,
255  int& nconv, double& diff)
256 {
257  m_timer.reset();
258  m_timer.start();
259 
260  int nin = m_fopr->field_nin();
261  int nvol = m_fopr->field_nvol();
262  int nex = m_fopr->field_nex();
263 
264  AFIELD axq(nin, nvol, nex);
265  AFIELD abq(nin, nvol, nex);
266 
267  AIndex_lex<real_t, AFIELD::IMPL> index_alt;
268 
269  if (m_fopr->needs_convert()) {
270  m_fopr->convert(abq, b);
271  } else {
272  convert(index_alt, abq, b);
273  }
274 
275  m_fopr->set_mode("D_prec");
276 
277  real_t diff2;
278 
279 #pragma omp parallel
280  {
281  m_solver->solve(axq, abq, nconv, diff2);
282  }
283  m_flop_count += m_solver->flop_count();
284 
285  m_fopr->mult(abq, axq, "Prec");
286 
287  if (m_fopr->needs_convert()) {
288  m_fopr->reverse(xq, abq);
289  } else {
290  reverse(index_alt, xq, abq);
291  }
292 
293  diff = double(diff2);
294 
295  m_timer.stop();
296  m_elapsed_time += m_timer.elapsed_sec();
297  // 1 solve + 2 mult; flop for the solve has been already accumulated
298  m_flop_count += 2 * m_fopr->flop_count();
299 }
300 
301 
302 //====================================================================
303 template<typename AFIELD>
304 void Fprop_alt_Standard_lex<AFIELD>::invert_DdagD_prec(Field& xq,
305  const Field& b,
306  int& nconv, double& diff)
307 {
308  vout.paranoiac(m_vl, "invert_DdagD_prec is called.\n");
309 
310  m_timer.reset();
311  m_timer.start();
312 
313  int nin = m_fopr->field_nin();
314  int nvol = m_fopr->field_nvol();
315  int nex = m_fopr->field_nex();
316 
317  AFIELD axq(nin, nvol, nex);
318  AFIELD abq(nin, nvol, nex);
319 
320  AIndex_lex<real_t, AFIELD::IMPL> index_alt;
321 
322  if (m_fopr->needs_convert()) {
323  m_fopr->convert(axq, b);
324  } else {
325  convert(index_alt, axq, b);
326  }
327 
328  m_fopr->mult(abq, axq, "Precdag");
329 
330  m_fopr->set_mode("DdagD_prec");
331 
332  real_t diff2;
333 
334 #pragma omp parallel
335  {
336  m_solver->solve(axq, abq, nconv, diff2);
337  }
338  m_flop_count += m_solver->flop_count();
339 
340  m_fopr->mult(abq, axq, "Prec");
341 
342  if (m_fopr->needs_convert()) {
343  m_fopr->reverse(xq, abq);
344  } else {
345  reverse(index_alt, xq, abq);
346  }
347 
348  diff = double(diff2);
349 
350  m_timer.stop();
351  m_elapsed_time += m_timer.elapsed_sec();
352  // 1 solve + 2 mult; flop for the solve has been already accumulated
353  m_flop_count += 2 * m_fopr->flop_count();
354 }
355 
356 
357 //====================================================================
358 template<typename AFIELD>
359 void Fprop_alt_Standard_lex<AFIELD>::invert(AFIELD& xq,
360  const AFIELD& b,
361  int& nconv, double& diff)
362 {
363  vout.paranoiac(m_vl, "%s: invert is called.\n", class_name.c_str());
364  vout.paranoiac(m_vl, "mode = %s.\n", m_mode.c_str());
365 
366  if (m_mode == "D") {
367  invert_D(xq, b, nconv, diff);
368  } else if (m_mode == "DdagD") {
369  invert_DdagD(xq, b, nconv, diff);
370  } else if (m_mode == "D_prec") {
371  invert_D_prec(xq, b, nconv, diff);
372  } else if (m_mode == "DdagD_prec") {
373  invert_DdagD_prec(xq, b, nconv, diff);
374  } else {
375  vout.crucial(m_vl, "%s: unsupported mode: %s\n",
376  class_name.c_str(), m_mode.c_str());
377  exit(EXIT_FAILURE);
378  }
379 }
380 
381 
382 //====================================================================
383 template<typename AFIELD>
384 void Fprop_alt_Standard_lex<AFIELD>::invert_D(AFIELD& axq,
385  const AFIELD& abq,
386  int& nconv, double& diff)
387 {
388  m_timer.reset();
389  m_timer.start();
390 
391  m_fopr->set_mode("D");
392 
393  real_t diff2;
394 
395 #pragma omp parallel
396  {
397  m_solver->solve(axq, abq, nconv, diff2);
398  }
399  diff = double(diff2);
400 
401  m_timer.stop();
402  m_elapsed_time += m_timer.elapsed_sec();
403  m_flop_count += m_solver->flop_count();
404 }
405 
406 
407 //====================================================================
408 template<typename AFIELD>
409 void Fprop_alt_Standard_lex<AFIELD>::invert_DdagD(AFIELD& axq,
410  const AFIELD& abq,
411  int& nconv, double& diff)
412 {
413  m_timer.reset();
414  m_timer.start();
415 
416  real_t diff2;
417 
418  m_fopr->set_mode("DdagD");
419 
420 #pragma omp parallel
421  {
422  m_solver->solve(axq, abq, nconv, diff2);
423  }
424  diff = double(diff2);
425 
426  m_timer.stop();
427  m_elapsed_time += m_timer.elapsed_sec();
428  m_flop_count += m_solver->flop_count();
429 }
430 
431 
432 //====================================================================
433 template<typename AFIELD>
434 void Fprop_alt_Standard_lex<AFIELD>::invert_D_prec(AFIELD& axq,
435  const AFIELD& abq,
436  int& nconv, double& diff)
437 {
438  m_timer.reset();
439  m_timer.start();
440 
441  int nin = m_fopr->field_nin();
442  int nvol = m_fopr->field_nvol();
443  int nex = m_fopr->field_nex();
444  AFIELD atmp(nin, nvol, nex);
445 
446  m_fopr->set_mode("D_prec");
447 
448  real_t diff2;
449 
450 #pragma omp parallel
451  {
452  m_solver->solve(atmp, abq, nconv, diff2);
453  }
454  m_flop_count += m_solver->flop_count();
455  m_fopr->mult(axq, atmp, "Prec");
456 
457  diff = double(diff2);
458 
459  m_timer.stop();
460  m_elapsed_time += m_timer.elapsed_sec();
461  // 1 solve + 2 mult; flop for the solve has been already accumulated
462  m_flop_count += 2 * m_fopr->flop_count();
463 }
464 
465 
466 //====================================================================
467 template<typename AFIELD>
468 void Fprop_alt_Standard_lex<AFIELD>::invert_DdagD_prec(AFIELD& axq,
469  const AFIELD& abq,
470  int& nconv, double& diff)
471 {
472  vout.paranoiac(m_vl, "invert_DdagD_prec is called.\n");
473 
474  m_timer.reset();
475  m_timer.start();
476 
477  int nin = m_fopr->field_nin();
478  int nvol = m_fopr->field_nvol();
479  int nex = m_fopr->field_nex();
480  AFIELD atmp(nin, nvol, nex);
481 
482  m_fopr->mult(axq, abq, "Precdag");
483 
484  m_fopr->set_mode("DdagD_prec");
485 
486  real_t diff2;
487 
488 #pragma omp parallel
489  {
490  m_solver->solve(atmp, axq, nconv, diff2);
491  }
492  m_flop_count += m_solver->flop_count();
493 
494  m_fopr->mult(axq, atmp, "Prec");
495 
496  diff = double(diff2);
497 
498  m_timer.stop();
499  m_elapsed_time += m_timer.elapsed_sec();
500  // 1 solve + 2 mult; flop for the solve has been already accumulated
501  m_flop_count += 2 * m_fopr->flop_count();
502 }
503 
504 
505 //====================================================================
506 template<typename AFIELD>
507 double Fprop_alt_Standard_lex<AFIELD>::flop_count()
508 {
509  return m_solver->flop_count();
510 }
511 
512 
513 //====================================================================
514 template<typename AFIELD>
515 void Fprop_alt_Standard_lex<AFIELD>::reset_performance()
516 {
517  m_flop_count = 0.0;
518  m_elapsed_time = 0.0;
519 }
520 
521 
522 //====================================================================
523 template<typename AFIELD>
524 void Fprop_alt_Standard_lex<AFIELD>::get_performance(double& flop_count,
525  double& elapsed_time)
526 {
527  flop_count = m_flop_count;
528  elapsed_time = m_elapsed_time;
529 }
530 
531 
532 //====================================================================
533 template<typename AFIELD>
534 void Fprop_alt_Standard_lex<AFIELD>::report_performance()
535 {
536  double flops = m_flop_count / m_elapsed_time;
537  double gflops = flops * 1.0e-9;
538 
539  // vout.general(m_vl, "\n");
540  vout.general(m_vl, "%s: solver performance:\n", class_name.c_str());
541  vout.general(m_vl, " Elapsed time = %14.6f sec\n", m_elapsed_time);
542  vout.general(m_vl, " Flop(total) = %18.0f\n", m_flop_count);
543  vout.general(m_vl, " Performance = %11.3f GFlops\n", gflops);
544 }
545 
546 
547 //====================================================================
548 template<typename AFIELD>
549 void Fprop_alt_Standard_lex<AFIELD>::mult_performance(
550  const std::string mode,
551  const int Nrepeat)
552 {
553  ThreadManager::assert_single_thread(class_name);
554 
555  int nin = m_fopr->field_nin();
556  int nvol = m_fopr->field_nvol();
557  int nex = m_fopr->field_nex();
558 
559  AFIELD axq(nin, nvol, nex), abq(nin, nvol, nex);
560  abq.set(0.0);
561  abq.set(0, 1.0);
562 
563  unique_ptr<Timer> timer(new Timer);
564 
565  std::string mode_prev = m_fopr->get_mode();
566 
567  m_fopr->set_mode(mode);
568 
569  timer->start();
570 
571 #pragma omp parallel
572  {
573  for (int i = 0; i < Nrepeat; ++i) {
574  m_fopr->mult(axq, abq);
575  m_fopr->mult(abq, axq);
576  }
577  }
578 
579  timer->stop();
580 
581  double flop_fopr = m_fopr->flop_count();
582  double flop_total = flop_fopr * double(2 * Nrepeat);
583 
584  double elapsed_time = timer->elapsed_sec();
585  double flops = flop_total / elapsed_time;
586  double gflops = flops * 1.0e-9;
587 
588  vout.general(m_vl, "\n");
589  vout.general(m_vl, "%s: mult performance:\n", class_name.c_str());
590  vout.general(m_vl, " mult mode = %s\n", mode.c_str());
591  vout.general(m_vl, " Elapsed time = %14.6f sec\n", elapsed_time);
592  vout.general(m_vl, " Flop(Fopr) = %18.0f\n", flop_fopr);
593  vout.general(m_vl, " Flop(total) = %18.0f\n", flop_total);
594  vout.general(m_vl, " Performance = %11.3f GFlops\n", gflops);
595 
596  m_fopr->set_mode(mode_prev);
597 }
598 
599 
600 //============================================================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
Bridge::BridgeIO::detailed
void detailed(const char *format,...)
Definition: bridgeIO.cpp:281
BridgeACC::convert
void convert(double *v, double *w, int nin, int nvol, int nvol_pad)
Fprop_alt_Standard_lex::mult_performance
void mult_performance(const std::string mode, const int Nrepeat)
Definition: fprop_alt_Standard_lex-tmpl.h:549
Timer
Definition: timer.h:31
Fprop_alt_Standard_lex::flop_count
double flop_count()
Definition: fprop_alt_Standard_lex-tmpl.h:507
Bridge::BridgeIO::paranoiac
void paranoiac(const char *format,...)
Definition: bridgeIO.cpp:300
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:304
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:534
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:524
Director_Smear
Manager of smeared configurations.
Definition: director_Smear.h:32
Fprop_alt_Standard_lex::invert_DdagD
void invert_DdagD(Field &, const Field &, int &, double &)
Definition: fprop_alt_Standard_lex-tmpl.h:190
Fprop_alt_Standard_lex::invert_D_prec
void invert_D_prec(Field &, const Field &, int &, double &)
Definition: fprop_alt_Standard_lex-tmpl.h:253
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
Fprop_alt_Standard_lex::reset_performance
void reset_performance()
Definition: fprop_alt_Standard_lex-tmpl.h:515
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:242
Field
Container of Field-type object.
Definition: field.h:46
Bridge::BridgeIO::general
void general(const char *format,...)
Definition: bridgeIO.cpp:262
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:572
BridgeACC::reverse
void reverse(double *v, double *w, int nin, int nvol, int nvol_pad)