Bridge++  Ver. 1.1.x
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
action_F_Ratio_lex.cpp
Go to the documentation of this file.
1 
14 #include "action_F_Ratio_lex.h"
15 
16 //====================================================================
17 void Action_F_Ratio_lex::set_parameters(const Parameters& params)
18 {
19  Bridge::VerboseLevel vl = params.get_VerboseLevel();
20 
21  set_parameter_verboselevel(vl);
22 }
23 
24 
25 //====================================================================
26 void Action_F_Ratio_lex::set_parameters()
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_Ratio_lex:\n");
34 
35  m_force.reset(NinG, Nvol, Ndim);
36 
37  // link variable update flag
38  m_status_linkv = 0;
39 }
40 
41 
42 //====================================================================
43 void Action_F_Ratio_lex::set_config(Field *U)
44 {
45  m_U = U;
46 
47  m_fopr_prec->set_config(U);
48  m_fopr_prec_force->set_config(U);
49 
50  m_fopr->set_config(U);
51  m_fopr_force->set_config(U);
52 }
53 
54 
55 //====================================================================
56 double Action_F_Ratio_lex::langevin(RandomNumbers *rand)
57 {
58  int Nvol = CommonParameters::Nvol();
59  int Ndim = CommonParameters::Ndim();
60 
61  int NinF = m_fopr_prec->field_nin();
62  int NvolF = m_fopr_prec->field_nvol();
63  int NexF = m_fopr_prec->field_nex();
64  int size_psf = NinF * NvolF * NexF * CommonParameters::NPE();
65 
66  assert(NvolF == Nvol);
67  m_psf.reset(NinF, NvolF, NexF);
68 
69  vout.general(m_vl, " Action_F_Ratio_lex: %s\n", m_label.c_str());
70 
71  Field xi(NinF, NvolF, NexF);
72  rand->gauss_lex_global(xi);
73 
74  m_fopr_prec->set_config(m_U);
75  m_fopr->set_config(m_U);
76 
77  Field v1(NinF, NvolF, NexF), v2(NinF, NvolF, NexF);
78 
79  m_fopr->set_mode("H");
80  v2 = m_fopr->mult_dag(xi);
81 
82  m_fopr_prec->set_mode("H");
83  v1 = m_fopr_prec->mult_dag(v2);
84 
85  int Nconv;
86  double diff;
87 
88  m_fprop_H_prec->set_config(m_U);
89  m_fprop_H_prec->invert_DdagD(m_psf, v1, Nconv, diff);
90  vout.general(m_vl, " Nconv = %d diff = %.8e\n", Nconv, diff);
91 
92  double xi2 = xi.norm();
93  double H_psf = xi2 * xi2;
94 
95  vout.general(m_vl, " H_Fratio = %18.8f\n", H_psf);
96  vout.general(m_vl, " H_F/dof = %18.8f\n", H_psf / size_psf);
97 
98  return H_psf;
99 }
100 
101 
102 //====================================================================
103 double Action_F_Ratio_lex::calcH()
104 {
105  int Nvol = CommonParameters::Nvol();
106  int Ndim = CommonParameters::Ndim();
107 
108  int NinF = m_fopr_prec->field_nin();
109  int NvolF = m_fopr_prec->field_nvol();
110  int NexF = m_fopr_prec->field_nex();
111  int size_psf = NinF * NvolF * NexF * CommonParameters::NPE();
112 
113  Field v1(NinF, NvolF, NexF), v2(NinF, NvolF, NexF);
114 
115  vout.general(m_vl, " Action_F_Ratio_lex: %s\n", m_label.c_str());
116 
117  m_fopr_prec->set_config(m_U);
118  m_fopr->set_config(m_U);
119 
120  m_fopr_prec->set_mode("H");
121  v1 = m_fopr_prec->mult(m_psf);
122 
123  int Nconv;
124  double diff;
125 
126  m_fprop_H->set_config(m_U);
127  m_fprop_H->invert_DdagD(v2, v1, Nconv, diff);
128  vout.general(m_vl, " Nconv = %d diff = %.8e\n", Nconv, diff);
129 
130  double H_psf = v1 * v2;
131 
132  vout.general(m_vl, " H_Fratio = %18.8f\n", H_psf);
133  vout.general(m_vl, " H_F/dof = %18.8f\n", H_psf / size_psf);
134 
135  return H_psf;
136 }
137 
138 
139 //====================================================================
140 const Field Action_F_Ratio_lex::force()
141 {
142  if (m_status_linkv == 0) {
143  int Nin = m_U->nin();
144  int Nvol = m_U->nvol();
145  int Nex = m_U->nex();
146  int Nc = CommonParameters::Nc();
147  int Ndim = CommonParameters::Ndim();
148 
149  int NinF = m_fopr_prec->field_nin();
150  int NvolF = m_fopr_prec->field_nvol();
151  int NexF = m_fopr_prec->field_nex();
152  Field eta(NinF, NvolF, NexF);
153 
154  vout.general(m_vl, " Action_F_Ratio_lex: %s\n", m_label.c_str());
155 
156  m_fopr_prec->set_config(m_U);
157  m_fopr->set_config(m_U);
158  m_fopr_prec_force->set_config(m_U);
159  m_fopr_force->set_config(m_U);
160 
161  Field v1(NinF, NvolF, NexF), v2(NinF, NvolF, NexF);
162 
163  m_fopr_prec->set_mode("H");
164  v1 = m_fopr_prec->mult(m_psf);
165 
166  int Nconv;
167  double diff;
168 
169  m_fprop_MD->set_config(m_U);
170  m_fprop_MD->invert_DdagD(v2, v1, Nconv, diff);
171  vout.general(m_vl, " Solver: Nconv = %6d diff = %12.6e\n", Nconv, diff);
172 
173  Field force(Nin, Nvol, Nex);
174 
175  force = m_fopr_force->force_core(v2);
176 
177  m_fopr_prec_force->set_mode("H");
178  force -= m_fopr_prec_force->force_core1(v2, m_psf);
179 
180  m_fopr_prec_force->set_mode("Hdag");
181  force -= m_fopr_prec_force->force_core1(m_psf, v2);
182 
183  m_force = force;
184  ++m_status_linkv;
185 
186  double Fave, Fmax, Fdev;
187  m_force.stat(Fave, Fmax, Fdev);
188  vout.general(m_vl, " Fratio_ave = %12.6f Fratio_max = %12.6f Fratio_dev = %12.6f\n",
189  Fave, Fmax, Fdev);
190 
191  return m_force;
192  } else {
193  vout.general(m_vl, " Action_F_Ratio_lex returns previous force.\n");
194  return m_force;
195  }
196 }
197 
198 
199 //====================================================================
200 //============================================================END=====