Bridge++  Ver. 2.0.2
fprop_alt_Standard_lex_Mixedprec-tmpl.h
Go to the documentation of this file.
1 
10 #include "Measurements/Fermion/fprop_alt_Standard_lex_Mixedprec.h"
11 
12 #include "lib/ResourceManager/threadManager.h"
13 
14 
15 template<typename AFIELD_d, typename AFIELD_f>
16 const std::string Fprop_alt_Standard_lex_Mixedprec<AFIELD_d, AFIELD_f>::class_name
17  = "Fprop_alt_Standard_lex_Mixedprec";
18 
19 //====================================================================
20 template<typename AFIELD_d, typename AFIELD_f>
21 void Fprop_alt_Standard_lex_Mixedprec<AFIELD_d, AFIELD_f>::init(const Parameters& params_fopr,
22  const Parameters& params_solver)
23 {
24  // this constructor assumes that the factories are available.
25  vout.general(m_vl, "\n");
26  vout.general(m_vl, "%s: being setup (without link smearing).\n",
27  class_name.c_str());
28 
29  m_dr_smear = 0;
30 
31  string fopr_type = params_fopr.get_string("fermion_type");
32 
33  m_kernel_d = 0;
34  m_fopr_d = AFopr<AFIELD_d>::New(fopr_type, params_fopr);
35 
36  m_kernel_f = 0;
37  m_fopr_f = AFopr<AFIELD_f>::New(fopr_type, params_fopr);
38 
39  // solver parameters for preconditioner solver
40  Parameters params_solver2 = params_solver;
41  params_solver2.set_double("convergence_criterion_squared", 1.0e-14);
42 
43  string solver_type = params_solver.get_string("solver_type");
44  m_solver_prec = ASolver<AFIELD_f>::New(solver_type, m_fopr_f);
45  m_solver_prec->set_parameters(params_solver2);
46 
47  m_precond
48  = new APrecond_Mixedprec<AFIELD_d, AFIELD_f>(m_solver_prec);
49 
50  m_solver
51  = new ASolver_FBiCGStab<AFIELD_d>(m_fopr_d, m_precond);
52  // ASolver_BiCGStab_Precond<AFIELD_d>(m_fopr_d, m_precond);
53  m_solver->set_parameters(params_solver);
54 
55  reset_performance();
56 
57  vout.general(m_vl, "%s: setup finished.\n", class_name.c_str());
58 }
59 
60 
61 //====================================================================
62 template<typename AFIELD_d, typename AFIELD_f>
63 void Fprop_alt_Standard_lex_Mixedprec<AFIELD_d, AFIELD_f>::init(const Parameters& params_fopr,
64  const Parameters& params_solver,
65  Director_Smear *dr_smear)
66 {
67  // this constructor assumes that the factories are available.
68  vout.general(m_vl, "\n");
69  vout.general(m_vl, "%s: being setup (with link smearing).\n",
70  class_name.c_str());
71 
72  m_dr_smear = dr_smear;
73 
74  string fopr_type = params_fopr.get_string("fermion_type");
75 
76  m_kernel_d = AFopr<AFIELD_d>::New(fopr_type, params_fopr);
77  m_fopr_d = new AFopr_Smeared<AFIELD_d>(m_kernel_d, m_dr_smear);
78 
79  m_kernel_f = AFopr<AFIELD_f>::New(fopr_type, params_fopr);
80  m_fopr_f = new AFopr_Smeared<AFIELD_f>(m_kernel_f, m_dr_smear);
81 
82  // solver parameters for preconditioner solver
83  Parameters params_solver2 = params_solver;
84  params_solver2.set_double("convergence_criterion_squared", 1.0e-14);
85 
86  string solver_type = params_solver.get_string("solver_type");
87  m_solver_prec = ASolver<AFIELD_f>::New(solver_type, m_fopr_f);
88  m_solver_prec->set_parameters(params_solver2);
89 
90  m_precond
91  = new APrecond_Mixedprec<AFIELD_d, AFIELD_f>(m_solver_prec);
92 
93  m_solver = new
94  ASolver_BiCGStab_Precond<AFIELD_d>(m_fopr_d, m_precond);
95  m_solver->set_parameters(params_solver);
96 
97  reset_performance();
98 
99  vout.general(m_vl, "%s: setup finished.\n", class_name.c_str());
100 }
101 
102 
103 //====================================================================
104 template<typename AFIELD_d, typename AFIELD_f>
105 void Fprop_alt_Standard_lex_Mixedprec<AFIELD_d, AFIELD_f>::tidyup()
106 {
107  delete m_solver;
108  delete m_precond;
109  delete m_solver_prec;
110  delete m_fopr_d;
111  delete m_fopr_f;
112  if (m_kernel_d != 0) delete m_kernel_d;
113  if (m_kernel_f != 0) delete m_kernel_f;
114 }
115 
116 
117 //====================================================================
118 template<typename AFIELD_d, typename AFIELD_f>
119 void Fprop_alt_Standard_lex_Mixedprec<AFIELD_d, AFIELD_f>::set_config(Field *U)
120 {
121  m_fopr_d->set_config(U);
122  m_fopr_f->set_config(U);
123 }
124 
125 
126 //====================================================================
127 template<typename AFIELD_d, typename AFIELD_f>
128 void Fprop_alt_Standard_lex_Mixedprec<AFIELD_d, AFIELD_f>::invert(Field& xq, const Field& b,
129  int& nconv, double& diff)
130 {
131  vout.paranoiac(m_vl, "%s: invert is called.\n", class_name.c_str());
132  vout.paranoiac(m_vl, "mode = %s.\n", m_mode.c_str());
133 
134  if (m_mode == "D") {
135  invert_D(xq, b, nconv, diff);
136  } else if (m_mode == "DdagD") {
137  invert_DdagD(xq, b, nconv, diff);
138  } else {
139  vout.crucial(m_vl, "%s: unsupported mode: %s\n",
140  class_name.c_str(), m_mode.c_str());
141  exit(EXIT_FAILURE);
142  }
143 }
144 
145 
146 //====================================================================
147 template<typename AFIELD_d, typename AFIELD_f>
148 void Fprop_alt_Standard_lex_Mixedprec<AFIELD_d, AFIELD_f>::invert_D(Field& xq, const Field& b,
149  int& nconv, double& diff)
150 {
151  m_timer.reset();
152  m_timer.start();
153 
154  int nin = m_fopr_d->field_nin();
155  int nvol = m_fopr_d->field_nvol();
156  int nex = m_fopr_d->field_nex();
157 
158  AFIELD_d axq(nin, nvol, nex);
159  AFIELD_d abq(nin, nvol, nex);
160 
161  AIndex_lex<real_t, AFIELD_d::IMPL> index_alt;
162 
163 #pragma omp parallel
164  {
165  if (m_fopr_d->needs_convert()) {
166  m_fopr_d->convert(abq, b);
167  } else {
168  convert(index_alt, abq, b);
169  }
170  }
171 
172  m_fopr_d->set_mode("D");
173  m_fopr_f->set_mode("D");
174 
175  double diff2;
176 
177 #pragma omp parallel
178  {
179  m_solver->solve(axq, abq, nconv, diff2);
180 #pragma omp barrier
181 
182  if (m_fopr_d->needs_convert()) {
183  m_fopr_d->reverse(xq, axq);
184  } else {
185  reverse(index_alt, xq, axq);
186  }
187  }
188  diff = double(diff2);
189 
190  m_timer.stop();
191  m_elapsed_time += m_timer.elapsed_sec();
192  m_flop_count += m_solver->flop_count();
193 }
194 
195 
196 //====================================================================
197 template<typename AFIELD_d, typename AFIELD_f>
198 void Fprop_alt_Standard_lex_Mixedprec<AFIELD_d, AFIELD_f>::invert_DdagD(
199  Field& xq, const Field& b,
200  int& nconv, double& diff)
201 {
202  m_timer.reset();
203  m_timer.start();
204 
205  int nin = m_fopr_d->field_nin();
206  int nvol = m_fopr_d->field_nvol();
207  int nex = m_fopr_d->field_nex();
208 
209  AFIELD_d axq(nin, nvol, nex);
210  AFIELD_d abq(nin, nvol, nex);
211 
212  AIndex_lex<real_t, AFIELD_d::IMPL> index_alt;
213 
214 #pragma omp parallel
215  {
216  if (m_fopr_d->needs_convert()) {
217  m_fopr_d->convert(abq, b);
218  } else {
219  convert(index_alt, abq, b);
220  }
221  }
222 
223  double diff2;
224 
225  m_fopr_d->set_mode("DdagD");
226  m_fopr_f->set_mode("DdagD");
227 
228 #pragma omp parallel
229  {
230  m_solver->solve(axq, abq, nconv, diff2);
231 #pragma omp barrier
232  if (m_fopr_d->needs_convert()) {
233  m_fopr_d->reverse(xq, axq);
234  } else {
235  reverse(index_alt, xq, axq);
236  }
237  }
238 
239  diff = double(diff2);
240 
241  m_timer.stop();
242  m_elapsed_time += m_timer.elapsed_sec();
243  m_flop_count += m_solver->flop_count();
244 }
245 
246 
247 //====================================================================
248 template<typename AFIELD_d, typename AFIELD_f>
249 void Fprop_alt_Standard_lex_Mixedprec<AFIELD_d, AFIELD_f>::invert(
250  AFIELD_d& xq, const AFIELD_d& b,
251  int& nconv, double& diff)
252 {
253  vout.paranoiac(m_vl, "%s: invert is called.\n", class_name.c_str());
254  vout.paranoiac(m_vl, "mode = %s.\n", m_mode.c_str());
255 
256  if (m_mode == "D") {
257  invert_D(xq, b, nconv, diff);
258  } else if (m_mode == "DdagD") {
259  invert_DdagD(xq, b, nconv, diff);
260  } else {
261  vout.crucial(m_vl, "%s: unsupported mode: %s\n",
262  class_name.c_str(), m_mode.c_str());
263  exit(EXIT_FAILURE);
264  }
265 }
266 
267 
268 //====================================================================
269 template<typename AFIELD_d, typename AFIELD_f>
270 void Fprop_alt_Standard_lex_Mixedprec<AFIELD_d, AFIELD_f>::invert_D(
271  AFIELD_d& xq, const AFIELD_d& b,
272  int& nconv, double& diff)
273 {
274  m_timer.reset();
275  m_timer.start();
276 
277  int nin = m_fopr_d->field_nin();
278  int nvol = m_fopr_d->field_nvol();
279  int nex = m_fopr_d->field_nex();
280 
281  m_fopr_d->set_mode("D");
282  m_fopr_f->set_mode("D");
283 
284  double diff2;
285 
286 #pragma omp parallel
287  {
288  m_solver->solve(xq, b, nconv, diff2);
289  }
290 
291  diff = double(diff2);
292 
293  m_timer.stop();
294  m_elapsed_time += m_timer.elapsed_sec();
295  m_flop_count += m_solver->flop_count();
296 }
297 
298 
299 //====================================================================
300 template<typename AFIELD_d, typename AFIELD_f>
301 void Fprop_alt_Standard_lex_Mixedprec<AFIELD_d, AFIELD_f>::invert_DdagD(
302  AFIELD_d& xq, const AFIELD_d& b,
303  int& nconv, double& diff)
304 {
305  m_timer.reset();
306  m_timer.start();
307 
308  int nin = m_fopr_d->field_nin();
309  int nvol = m_fopr_d->field_nvol();
310  int nex = m_fopr_d->field_nex();
311 
312  m_fopr_d->set_mode("DdagD");
313  m_fopr_f->set_mode("DdagD");
314 
315  double diff2;
316 
317 #pragma omp parallel
318  {
319  m_solver->solve(xq, b, nconv, diff2);
320  }
321 
322  diff = double(diff2);
323 
324  m_timer.stop();
325  m_elapsed_time += m_timer.elapsed_sec();
326  m_flop_count += m_solver->flop_count();
327 }
328 
329 
330 //====================================================================
331 template<typename AFIELD_d, typename AFIELD_f>
332 double Fprop_alt_Standard_lex_Mixedprec<AFIELD_d, AFIELD_f>::flop_count()
333 {
334  return m_solver->flop_count();
335 }
336 
337 
338 //====================================================================
339 template<typename AFIELD_d, typename AFIELD_f>
340 void Fprop_alt_Standard_lex_Mixedprec<AFIELD_d, AFIELD_f>::reset_performance()
341 {
342  m_flop_count = 0.0;
343  m_elapsed_time = 0.0;
344 }
345 
346 
347 //====================================================================
348 template<typename AFIELD_d, typename AFIELD_f>
349 void Fprop_alt_Standard_lex_Mixedprec<AFIELD_d, AFIELD_f>::report_performance()
350 {
351  double flops = m_flop_count / m_elapsed_time;
352  double gflops = flops * 1.0e-9;
353 
354  vout.general(m_vl, "\n");
355  vout.general(m_vl, "%s: solver performance:\n", class_name.c_str());
356  vout.general(m_vl, " Elapsed time = %14.6f sec\n", m_elapsed_time);
357  vout.general(m_vl, " Flop(total) = %18.0f\n", m_flop_count);
358  vout.general(m_vl, " Performance = %11.3f GFlops\n", gflops);
359 }
360 
361 
362 //====================================================================
363 template<typename AFIELD_d, typename AFIELD_f>
364 void Fprop_alt_Standard_lex_Mixedprec<AFIELD_d, AFIELD_f>::mult_performance(
365  const std::string mode,
366  const int Nrepeat)
367 {
368  int nin = m_fopr_d->field_nin();
369  int nvol = m_fopr_d->field_nvol();
370  int nex = m_fopr_d->field_nex();
371 
372  unique_ptr<Timer> timer(new Timer);
373 
374  std::string mode_prev_d = m_fopr_d->get_mode();
375  std::string mode_prev_f = m_fopr_f->get_mode();
376 
377  m_fopr_d->set_mode(mode);
378  m_fopr_f->set_mode(mode);
379 
380  {
381  AFIELD_d axq(nin, nvol, nex), abq(nin, nvol, nex);
382  abq.set(0.0);
383  abq.set(0, 1.0);
384 
385  timer->start();
386 
387 #pragma omp parallel
388  {
389  for (int i = 0; i < Nrepeat; ++i) {
390  m_fopr_d->mult(axq, abq);
391  m_fopr_d->mult(abq, axq);
392  }
393  }
394 
395  timer->stop();
396  }
397 
398  double flop_fopr = m_fopr_d->flop_count();
399  double flop_total = flop_fopr * double(2 * Nrepeat);
400 
401  double elapsed_time = timer->elapsed_sec();
402  double flops = flop_total / elapsed_time;
403  double gflops = flops * 1.0e-9;
404 
405  vout.general(m_vl, "\n");
406  vout.general(m_vl, "%s: mult performance:\n", class_name.c_str());
407  vout.general(m_vl, " Double precision:\n");
408  vout.general(m_vl, " Elapsed time = %14.6f sec\n", elapsed_time);
409  vout.general(m_vl, " Flop(Fopr) = %18.0f\n", flop_fopr);
410  vout.general(m_vl, " Flop(total) = %18.0f\n", flop_total);
411  vout.general(m_vl, " Performance = %11.3f GFlops\n", gflops);
412 
413  {
414  AFIELD_f axq(nin, nvol, nex), abq(nin, nvol, nex);
415  abq.set(0.0);
416  abq.set(0, 1.0);
417 
418  timer->reset();
419  timer->start();
420 
421 #pragma omp parallel
422  {
423  for (int i = 0; i < Nrepeat; ++i) {
424  m_fopr_f->mult(axq, abq);
425  m_fopr_f->mult(abq, axq);
426  }
427  }
428 
429  timer->stop();
430  }
431 
432  flop_fopr = m_fopr_f->flop_count();
433  flop_total = flop_fopr * double(2 * Nrepeat);
434 
435  elapsed_time = timer->elapsed_sec();
436  flops = flop_total / elapsed_time;
437  gflops = flops * 1.0e-9;
438 
439  vout.general(m_vl, " Single precision:\n");
440  vout.general(m_vl, " Elapsed time = %14.6f sec\n", elapsed_time);
441  vout.general(m_vl, " Flop(Fopr) = %18.0f\n", flop_fopr);
442  vout.general(m_vl, " Flop(total) = %18.0f\n", flop_total);
443  vout.general(m_vl, " Performance = %11.3f GFlops\n", gflops);
444 
445  m_fopr_d->set_mode(mode_prev_d);
446  m_fopr_f->set_mode(mode_prev_f);
447 }
448 
449 
450 //============================================================END=====
AFopr
Definition: afopr.h:48
ASolver
Definition: asolver.h:23
ASolver_FBiCGStab
Definition: asolver_FBiCGStab.h:32
Parameters
Class for parameters.
Definition: parameters.h:46
AIndex_lex
Definition: aindex_lex_base.h:17
Parameters::set_double
void set_double(const string &key, const double value)
Definition: parameters.cpp:33
convert
void convert(INDEX &index, AFIELD &v, const Field &w)
Definition: afield-inc.h:129
AField< double, QXS >
Fprop_alt_Standard_lex_Mixedprec::invert_D
void invert_D(Field &, const Field &, int &, double &)
Definition: fprop_alt_Standard_lex_Mixedprec-tmpl.h:148
Fprop_alt_Standard_lex_Mixedprec::tidyup
void tidyup()
Definition: fprop_alt_Standard_lex_Mixedprec-tmpl.h:105
Timer
Definition: timer.h:31
Bridge::BridgeIO::paranoiac
void paranoiac(const char *format,...)
Definition: bridgeIO.cpp:238
Fprop_alt_Standard_lex_Mixedprec::init
void init(const Parameters &params_fopr, const Parameters &params_solver)
Definition: fprop_alt_Standard_lex_Mixedprec-tmpl.h:21
AFopr_Smeared
smeared fermion operator: alternative version.
Definition: afopr_Smeared.h:41
threadManager.h
Fprop_alt_Standard_lex_Mixedprec::report_performance
void report_performance()
Definition: fprop_alt_Standard_lex_Mixedprec-tmpl.h:349
Director_Smear
Manager of smeared configurations.
Definition: director_Smear.h:39
Fprop_alt_Standard_lex_Mixedprec::invert_DdagD
void invert_DdagD(Field &, const Field &, int &, double &)
Definition: fprop_alt_Standard_lex_Mixedprec-tmpl.h:198
APrecond_Mixedprec
Definition: aprecond_Mixedprec.h:28
Fprop_alt_Standard_lex_Mixedprec::invert
void invert(Field &, const Field &, int &, double &)
invert accordingly to the mode. [22 Sep 2018 H.Matsufuru]
Definition: fprop_alt_Standard_lex_Mixedprec-tmpl.h:128
fprop_alt_Standard_lex_Mixedprec.h
Fprop_alt_Standard_lex_Mixedprec::flop_count
double flop_count()
Definition: fprop_alt_Standard_lex_Mixedprec-tmpl.h:332
Fprop_alt_Standard_lex_Mixedprec::reset_performance
void reset_performance()
Definition: fprop_alt_Standard_lex_Mixedprec-tmpl.h:340
reverse
void reverse(INDEX &index, Field &v, const AFIELD &w)
Definition: afield-inc.h:159
ASolver_BiCGStab_Precond
Definition: asolver_BiCGStab_Precond.h:24
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
Fprop_alt_Standard_lex_Mixedprec
Get quark propagator for Fopr with lexical site index: alternative version.
Definition: fprop_alt_Standard_lex_Mixedprec.h:31
Bridge::BridgeIO::general
void general(const char *format,...)
Definition: bridgeIO.cpp:200
Fprop_alt_Standard_lex_Mixedprec::mult_performance
void mult_performance(const std::string mode, const int Nrepeat)
Definition: fprop_alt_Standard_lex_Mixedprec-tmpl.h:364
Fprop_alt_Standard_lex_Mixedprec::set_config
void set_config(Field *)
Definition: fprop_alt_Standard_lex_Mixedprec-tmpl.h:119
Bridge::vout
BridgeIO vout
Definition: bridgeIO.cpp:512