Bridge++  Ver. 1.1.x
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
action_F_Standard.cpp
Go to the documentation of this file.
1 
14 #include "action_F_Standard.h"
15 
16 //====================================================================
18 {
20 
22 }
23 
24 
25 //====================================================================
27 {
28  int Nc = CommonParameters::Nc();
29  int Nvol = CommonParameters::Nvol();
30  int Ndim = CommonParameters::Ndim();
31  int NinG = 2 * Nc * Nc;
32 
33  vout.general(m_vl, "Action_F_Standard:\n");
34 
35  m_force.reset(NinG, Nvol, Ndim);
36 
37  int Niter = 2000;
38  double Stop_cond = 1.0e-24;
39 
40  string str_solver_type = "CG";
41  m_solver = Solver::New(str_solver_type, m_fopr);
42  m_solver->set_parameters(Niter, Stop_cond);
43 
44  // link variable update flag
45  m_status_linkv = 0;
46 }
47 
48 
49 //====================================================================
51 {
52  int Nvol = CommonParameters::Nvol();
53  int Ndim = CommonParameters::Ndim();
54 
55  int NinF = m_fopr->field_nin();
56  int NvolF = m_fopr->field_nvol();
57  int NexF = m_fopr->field_nex();
58  int size_psf = NinF * NvolF * NexF * CommonParameters::NPE();
59 
60  assert(NvolF == Nvol);
61  m_psf.reset(NinF, NvolF, NexF);
62 
63  vout.general(m_vl, " Action_F_Standard: %s\n", m_label.c_str());
64 
65  Field xi(NinF, NvolF, NexF);
66  rand->gauss_lex_global(xi);
67 
69  m_fopr->set_mode("Ddag");
70  m_psf = m_fopr->mult(xi);
71 
72  double xi2 = xi.norm();
73  double H_psf = xi2 * xi2;
74 
75  vout.general(m_vl, " H_Fstandard = %18.8f\n", H_psf);
76  vout.general(m_vl, " H_F/dof = %18.8f\n", H_psf / size_psf);
77 
78  return H_psf;
79 }
80 
81 
82 //====================================================================
84 {
85  int Nvol = CommonParameters::Nvol();
86  int Ndim = CommonParameters::Ndim();
87 
88  int NinF = m_fopr->field_nin();
89  int NvolF = m_fopr->field_nvol();
90  int NexF = m_fopr->field_nex();
91  int size_psf = NinF * NvolF * NexF * CommonParameters::NPE();
92 
93  Field v1(NinF, NvolF, NexF);
94 
95  vout.general(m_vl, " Action_F_Standard: %s\n", m_label.c_str());
96 
98  m_fopr->set_mode("DdagD");
99 
100  int Nconv;
101  double diff;
102  m_solver->solve(v1, m_psf, Nconv, diff);
103 
104  vout.general(m_vl, " Nconv = %d diff = %.8e\n", Nconv, diff);
105 
106  double H_psf = v1 * m_psf;
107 
108  vout.general(m_vl, " H_Fstandard = %18.8f\n", H_psf);
109  vout.general(m_vl, " H_F/dof = %18.8f\n", H_psf / size_psf);
110 
111  return H_psf;
112 }
113 
114 
115 //====================================================================
117 {
118  if (m_status_linkv == 0) {
119  int Nin = m_U->nin();
120  int Nvol = m_U->nvol();
121  int Nex = m_U->nex();
122  int Nc = CommonParameters::Nc();
123  int Ndim = CommonParameters::Ndim();
124 
125  int NinF = m_fopr->field_nin();
126  int NvolF = m_fopr->field_nvol();
127  int NexF = m_fopr->field_nex();
128  Field eta(NinF, NvolF, NexF);
129 
130  vout.general(m_vl, " Action_F_Standard: %s\n", m_label.c_str());
131 
132  //- fermion inversion for smeared gauge field
134  m_fopr->set_mode("DdagD");
135 
136  int Nconv;
137  double diff;
138  m_solver->solve(eta, m_psf, Nconv, diff);
139 
140  vout.general(m_vl, " Solver: Nconv = %6d diff = %12.6e\n", Nconv, diff);
141 
142  //- force of smeared fermion operator
144 
145  Field force(Nin, Nvol, Nex);
146  Field force1(Nin, Nvol, Nex);
147 
148  force = m_fopr_force->force_core(eta);
149 
150  m_force = force;
151  ++m_status_linkv;
152 
153  double Fave, Fmax, Fdev;
154  m_force.stat(Fave, Fmax, Fdev);
155  vout.general(m_vl, " Fstandard_ave = %12.6f Fstandard_max = %12.6f Fstandard_dev = %12.6f\n",
156  Fave, Fmax, Fdev);
157 
158  return m_force;
159  } else {
160  vout.general(m_vl, " Fstandard returns previous force.\n");
161  return m_force;
162  }
163 }
164 
165 
166 //====================================================================
167 //============================================================END=====