CHROMA
two_flavor_monomial5d_w.h
Go to the documentation of this file.
1 // -*- C++ -*-
2 
3 /*! @file
4  * @brief Two flavor Monomials - gauge action or fermion binlinear contributions for HMC
5  */
6 
7 #ifndef __two_flavor_monomial5d_w_h__
8 #define __two_flavor_monomial5d_w_h__
9 
15 
16 #include <typeinfo>
17 
18 namespace Chroma
19 {
20  //-------------------------------------------------------------------------------------------
21  //! Exact 2 degen flavor fermact monomial in extra dimensions
22  /*! @ingroup monomial
23  *
24  * Exact 2 degen flavor fermact monomial. Preconditioning is not
25  * specified yet.
26  * Can supply a default dsdq and pseudoferm refresh algorithm
27  *
28  * CAVEAT: I assume there is only 1 pseudofermion field in the following
29  * so called TwoFlavorExact actions.
30  */
31  template<typename P, typename Q, typename Phi>
33  {
34  public:
35  //! virtual destructor:
37 
38  //! Compute the total action
39  virtual Double S(const AbsFieldState<P,Q>& s) = 0;
40 
41  //! Compute dsdq for the system...
42  /*! Actions of the form chi^dag*(M^dag*M)*chi */
43  virtual void dsdq(P& F, const AbsFieldState<P,Q>& s)
44  {
45  START_CODE();
46 
47  // SelfIdentification/Encapsultaion Rule
48  XMLWriter& xml_out = TheXMLLogWriter::Instance();
49  push(xml_out, "TwoFlavorExactWilsonTypeFermMonomial5D");
50 
51  /**** Identical code for unprec and even-odd prec case *****/
52 
53  // S_f = chi^dag*(M^dag*M)^(-1)*chi
54  // Here, M is some 5D operator
55  //
56  // Need
57  // dS_f/dU = -chi^dag * (M^dag*M)^(-1) * [d(M^dag)*M + M^dag*dM] * (M^dag*M)^(-1) * chi
58  //
59  // where psi = (M^dag*M)^(-1) * chi
60  //
61  // In Balint's notation, the result is
62  // \dot{S} = chi^dag*X - X^dag*\dot{M}^\dag*Y - Y^dag\dot{M}*X + X*chi
63  // where
64  // X = (M^dag*M)^(-1)*chi Y = M*X = (M^dag)^(-1)*chi
65  // In Robert's notation, X -> psi .
66  //
68 
69  // Create a state for linop
71 
72  // Get linear operator
74 
75  // Get system solver
77 
78  // Chrono predictor and inversion
79  multi1d<Phi> X(FA.size());
80  int n_count;
81  {
82  // CG Chrono predictor needs MdagM
85 
86  // Do the inversion
87  SystemSolverResults_t res = (*invMdagM)(X, getPhi());
88  n_count = res.n_count;
89 
90  // Register the new std::vector
91  (getMDSolutionPredictor()).newVector(X);
92  }
93 
94  // Get/construct the pseudofermion solution
95  multi1d<Phi> Y(FA.size());
96 
97  (*M)(Y, X, PLUS);
98 
99  M->deriv(F, X, Y, MINUS);
100 
101  // fold M^dag into X^dag -> Y !!
102  P F_tmp;
103  M->deriv(F_tmp, Y, X, PLUS);
104  F += F_tmp;
105 
106  for(int mu=0; mu < F.size(); ++mu) {
107  F[mu] *= Real(-1);
108  }
109 
110  state->deriv(F);
111  write(xml_out, "n_count", n_count);
112  monitorForces(xml_out, "Forces", F);
113  pop(xml_out);
114 
115  END_CODE();
116  }
117 
118  //! Refresh pseudofermions
119  virtual void refreshInternalFields(const AbsFieldState<P,Q>& field_state)
120  {
121  START_CODE();
122 
123  // Heatbath all the fields
124 
125  // Get at the fermion action
127 
128  // Create a Connect State, apply fermionic boundaries
129  Handle< FermState<Phi,P,Q> > f_state(FA.createState(field_state.getQ()));
130 
131  // Create a linear operator
133 
134  const int N5 = FA.size();
135  multi1d<Phi> eta(N5);
136  eta = zero;
137 
138  // Fill the eta field with gaussian noise
139  for(int s=0; s < N5; ++s)
140  gaussian(eta[s], M->subset());
141 
142  // Account for fermion BC by modifying the proposed field
143  FA.getFermBC().modifyF(eta);
144 
145  // Temporary: Move to correct normalisation
146  for(int s=0; s < N5; ++s)
147  eta[s][M->subset()] *= sqrt(0.5);
148 
149  // Build phi = M^dag * eta
150  (*M)(getPhi(), eta, MINUS);
151 
152  // Reset the chronological predictor
153  QDPIO::cout << "TwoFlavWilson5DMonomial: resetting Predictor at end of field refresh" << std::endl;
154  getMDSolutionPredictor().reset();
155 
156  END_CODE();
157  }
158 
159  virtual void setInternalFields(const Monomial<P,Q>& m)
160  {
161  START_CODE();
162 
163  try
164  {
166 
167  // Do a resize here -- otherwise if the fields have not yet
168  // been refreshed there may be trouble
169  getPhi().resize(fm.getPhi().size());
170 
171  for(int i=0 ; i < fm.getPhi().size(); i++) {
172  (getPhi())[i] = (fm.getPhi())[i];
173  }
174  }
175  catch(std::bad_cast) {
176  QDPIO::cerr << "Failed to cast input Monomial to TwoFlavorExactWilsonTypeFermMonomial5D" << std::endl;
177  QDP_abort(1);
178  }
179 
180 
181  END_CODE();
182  }
183 
184  //! Reset predictors
185  virtual void resetPredictors() {
186  getMDSolutionPredictor().reset();
187 
188  }
189 
190  protected:
191  //! Accessor for pseudofermion with Pf index i (read only)
192  virtual const multi1d<Phi>& getPhi() const = 0;
193 
194  //! mutator for pseudofermion with Pf index i
195  virtual multi1d<Phi>& getPhi() = 0;
196 
197  //! Get at fermion action
198  virtual const WilsonTypeFermAct5D<Phi,P,Q>& getFermAct() const = 0;
199 
200  //! Get inverter params
201  virtual const GroupXML_t& getInvParams() const = 0;
202 
203  //! Get the initial guess predictor
205  };
206 
207 
208  //-------------------------------------------------------------------------------------------
209  //! Exact 2 degen flavor unpreconditioned fermact monomial living in extra dimensions
210  /*! @ingroup monomial
211  *
212  * Exact 2 degen flavor unpreconditioned fermact monomial.
213  *
214  * CAVEAT: I assume there is only 1 pseudofermion field in the following
215  * so called TwoFlavorExact actions.
216  */
217  template<typename P, typename Q, typename Phi>
219  {
220  public:
221  //! virtual destructor:
223 
224  //! Compute the total action
225  virtual Double S(const AbsFieldState<P,Q>& s)
226  {
227  START_CODE();
228 
229  // SelfEncapsulation/Identification Rule
230  XMLWriter& xml_out = TheXMLLogWriter::Instance();
231  push(xml_out, "TwoFlavorExactUnprecWilsonTypeFermMonomial5D");
232 
234 
236 
237  // Get system solver
239 
240  multi1d<Phi> X(FA.size());
241 
242  // Energy calc does not use chrono predictor
243  X = zero;
244 
245  // Reset the chrono predictor.
246  QDPIO::cout << "TwoFlavWilson5DMonomial: energy calc solve" << std::endl;
247  getMDSolutionPredictor().reset();
248 
249  // Do the inversion
250  SystemSolverResults_t res = (*invMdagM)(X, getPhi());
251  int n_count = res.n_count;
252 
253  // Action on the entire lattice
254  Double action = zero;
255  for(int s=0; s < FA.size(); ++s)
256  action += innerProductReal(getPhi()[s], X[s]);
257 
258  write(xml_out, "n_count", n_count);
259  write(xml_out, "S", action);
260  pop(xml_out);
261 
262  END_CODE();
263 
264  return action;
265  }
266 
267 
268  protected:
269  //! Accessor for pseudofermion with Pf index i (read only)
270  virtual const multi1d<Phi>& getPhi() const = 0;
271 
272  //! mutator for pseudofermion with Pf index i
273  virtual multi1d<Phi>& getPhi() = 0;
274 
275  //! Get at fermion action
276  virtual const UnprecWilsonTypeFermAct5D<Phi,P,Q>& getFermAct() const = 0;
277 
278  //! Get inverter params
279  virtual const GroupXML_t& getInvParams() const = 0;
280 
281  //! Get the initial guess predictor
283  };
284 
285 
286  //-------------------------------------------------------------------------------------------
287  //! Exact 2 degen flavor even-odd preconditioned fermact monomial living in extra dimensions
288  /*! @ingroup monomial
289  *
290  * Exact 2 degen flavor even-odd preconditioned fermact monomial.
291  * Can supply a default dsdq algorithm
292  */
293  template<typename P, typename Q, typename Phi>
295  {
296  public:
297  //! virtual destructor:
299 
300  //! Even even contribution (eg ln det Clover)
301  virtual Double S_even_even(const AbsFieldState<P,Q>& s) = 0;
302 
303  //! Compute the odd odd contribution (eg
305  {
306  START_CODE();
307 
308  XMLWriter& xml_out = TheXMLLogWriter::Instance();
309  push(xml_out, "S_odd_odd");
310 
312 
314 
315  // Get system solver
317 
318  // Need way to get gauge state from AbsFieldState<P,Q>
320 
321  // Get the X fields
322  multi1d<Phi> X(FA.size());
323 
324  // Chrono predictor not used in energy calculation
325  X = zero;
326 
327  // Reset the chrono predictor.
328  QDPIO::cout << "TwoFlavWilson5DMonomial: energy calc solve" << std::endl;
329  getMDSolutionPredictor().reset();
330 
331  // Do the inversion
332  SystemSolverResults_t res = (*invMdagM)(X, getPhi());
333  int n_count = res.n_count;
334 
335  // Total odd-subset action. NOTE: QDP has norm2(multi1d) but not innerProd
336  Double action = innerProductReal(getPhi(), X, M->subset());
337 
338  write(xml_out, "n_count", n_count);
339  write(xml_out, "S_oo", action);
340  pop(xml_out);
341 
342  END_CODE();
343 
344  return action;
345  }
346 
347  //! Compute the total action
349  {
350  START_CODE();
351 
352  XMLWriter& xml_out=TheXMLLogWriter::Instance();
353  push(xml_out, "TwoFlavorExactEvenOddPrecWilsonTypeFermMonomial5D");
354 
355  Double action = S_even_even(s) + S_odd_odd(s);
356 
357  write(xml_out, "S", action);
358  pop(xml_out);
359 
360  END_CODE();
361 
362  return action;
363  }
364 
365  protected:
366  //! Get at fermion action
368 
369  //! Get inverter params
370  virtual const GroupXML_t& getInvParams() const = 0;
371 
372  //! Get the initial guess predictor
374 
375  //! Accessor for pseudofermion with Pf index i (read only)
376  virtual const multi1d<Phi>& getPhi() const = 0;
377 
378  //! mutator for pseudofermion with Pf index i
379  virtual multi1d<Phi>& getPhi() = 0;
380  };
381 
382  //-------------------------------------------------------------------------------------------
383  //! Exact 2 degen flavor even-odd preconditioned fermact monomial living in extra dimensions
384  /*! @ingroup monomial
385  *
386  * Exact 2 degen flavor even-odd preconditioned fermact monomial.
387  * Can supply a default dsdq algorithm
388  */
389  template<typename P, typename Q, typename Phi>
391  {
392  public:
393  //! virtual destructor:
395 
396  //! Even even contribution (eg ln det Clover)
398  return Double(0);
399  }
400 
401  protected:
402  //! Get at fermion action
404  };
405 
406 
407 
408 }
409 
410 
411 #endif
Monomials - gauge action or fermion binlinear contributions for HMC.
Chronological predictor for HMC.
Abstract interface for a Chronological Solution predictor in 5D.
Abstract field state.
Definition: field_state.h:27
virtual const Q & getQ(void) const =0
virtual DiffLinearOperatorArray< T, P, Q > * linOp(Handle< FermState< T, P, Q > > state) const =0
Produce a linear operator for this action.
Even-odd preconditioned Wilson-like fermion actions including derivatives.
Even-odd preconditioned Wilson-like fermion actions including derivatives.
virtual EvenOddPrecLinearOperatorArray< T, P, Q > * linOp(Handle< FermState< T, P, Q > > state) const =0
Override to produce an even-odd prec. linear operator for this action.
Fermionic monomials (binlinears in fermion fields)
Definition: abs_monomial.h:280
virtual int size() const =0
Expected length of array index.
virtual const FermBC< T, P, Q > & getFermBC() const
Return the fermion BC object for this action.
Definition: fermact.h:79
virtual FermState< T, P, Q > * createState(const Q &q) const
Given links (coordinates Q) create the state needed for the linear operators.
Definition: fermact.h:59
Class for counted reference semantics.
Definition: handle.h:33
An abstract monomial class, for inexact algorithms.
Definition: abs_monomial.h:43
static T & Instance()
Definition: singleton.h:432
Exact 2 degen flavor even-odd preconditioned fermact monomial living in extra dimensions.
virtual const EvenOddPrecConstDetWilsonTypeFermAct5D< Phi, P, Q > & getFermAct() const =0
Get at fermion action.
virtual Double S_even_even(const AbsFieldState< P, Q > &s)
Even even contribution (eg ln det Clover)
Exact 2 degen flavor even-odd preconditioned fermact monomial living in extra dimensions.
virtual Double S_even_even(const AbsFieldState< P, Q > &s)=0
Even even contribution (eg ln det Clover)
virtual Double S_odd_odd(const AbsFieldState< P, Q > &s)
Compute the odd odd contribution (eg.
virtual multi1d< Phi > & getPhi()=0
mutator for pseudofermion with Pf index i
virtual const EvenOddPrecWilsonTypeFermAct5D< Phi, P, Q > & getFermAct() const =0
Get at fermion action.
virtual AbsChronologicalPredictor5D< Phi > & getMDSolutionPredictor()=0
Get the initial guess predictor.
virtual const GroupXML_t & getInvParams() const =0
Get inverter params.
Double S(const AbsFieldState< P, Q > &s)
Compute the total action.
virtual const multi1d< Phi > & getPhi() const =0
Accessor for pseudofermion with Pf index i (read only)
Exact 2 degen flavor unpreconditioned fermact monomial living in extra dimensions.
virtual Double S(const AbsFieldState< P, Q > &s)
Compute the total action.
virtual AbsChronologicalPredictor5D< Phi > & getMDSolutionPredictor()=0
Get the initial guess predictor.
virtual const UnprecWilsonTypeFermAct5D< Phi, P, Q > & getFermAct() const =0
Get at fermion action.
virtual const GroupXML_t & getInvParams() const =0
Get inverter params.
virtual multi1d< Phi > & getPhi()=0
mutator for pseudofermion with Pf index i
virtual const multi1d< Phi > & getPhi() const =0
Accessor for pseudofermion with Pf index i (read only)
Exact 2 degen flavor fermact monomial in extra dimensions.
virtual void setInternalFields(const Monomial< P, Q > &m)
Copy pseudofermions if any.
virtual void refreshInternalFields(const AbsFieldState< P, Q > &field_state)
Refresh pseudofermions.
virtual const WilsonTypeFermAct5D< Phi, P, Q > & getFermAct() const =0
Get at fermion action.
virtual const GroupXML_t & getInvParams() const =0
Get inverter params.
virtual void dsdq(P &F, const AbsFieldState< P, Q > &s)
Compute dsdq for the system...
virtual multi1d< Phi > & getPhi()=0
mutator for pseudofermion with Pf index i
virtual AbsChronologicalPredictor5D< Phi > & getMDSolutionPredictor()=0
Get the initial guess predictor.
virtual Double S(const AbsFieldState< P, Q > &s)=0
Compute the total action.
virtual const multi1d< Phi > & getPhi() const =0
Accessor for pseudofermion with Pf index i (read only)
Unpreconditioned Wilson-like fermion actions in extra dims with derivatives.
Definition: fermact.orig.h:571
Wilson-like fermion actions.
Definition: fermact.orig.h:403
virtual MdagMSystemSolverArray< T > * invMdagM(Handle< FermState< T, P, Q > > state, const GroupXML_t &invParam) const
Return a linear operator solver for this action to solve MdagM*psi=chi.
virtual DiffLinearOperatorArray< T, P, Q > * lMdagM(Handle< FermState< T, P, Q > > state) const
Produce a linear operator M^dag.M for this action.
Definition: fermact.orig.h:410
int mu
Definition: cool.cc:24
Even-odd const determinant Wilson-like fermact.
Helper function for calculating forces.
void write(XMLWriter &xml, const std::string &path, const AsqtadFermActParams &param)
Writer parameters.
void monitorForces(XMLWriter &xml_out, const std::string &path, const multi1d< LatticeColorMatrix > &F)
Calculate and write out forces.
static int m[4]
Definition: make_seeds.cc:16
Asqtad Staggered-Dirac operator.
Definition: klein_gord.cc:10
gaussian(aux)
int n_count
Definition: invbicg.cc:78
push(xml_out,"Condensates")
int i
Definition: pbg5p_w.cc:55
@ MINUS
Definition: chromabase.h:45
@ PLUS
Definition: chromabase.h:45
pop(xml_out)
LatticeFermion eta
Definition: mespbg5p_w.cc:37
START_CODE()
const WilsonTypeFermAct< multi1d< LatticeFermion > > Handle< const ConnectState > state
Definition: pbg5p_w.cc:28
Double zero
Definition: invbicg.cc:106
multi1d< LatticeFermion > s(Ncb)
FloatingPoint< double > Double
Definition: gtest.h:7351
Hold group xml and type id.
Holds return info from SystemSolver call.
Definition: syssolver.h:17
Handle< LinearOperator< T > > MdagM
multi1d< LatticeColorMatrix > P
Definition: t_clover.cc:13
Wilson-like fermion actions.
static INTERNAL_PRECISION F