flexiblesusy is hosted by Hepforge, IPPP Durham
FlexibleSUSY
lattice_constraint.hpp
Go to the documentation of this file.
1// ====================================================================
2// This file is part of FlexibleSUSY.
3//
4// FlexibleSUSY is free software: you can redistribute it and/or modify
5// it under the terms of the GNU General Public License as published
6// by the Free Software Foundation, either version 3 of the License,
7// or (at your option) any later version.
8//
9// FlexibleSUSY is distributed in the hope that it will be useful, but
10// WITHOUT ANY WARRANTY; without even the implied warranty of
11// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12// General Public License for more details.
13//
14// You should have received a copy of the GNU General Public License
15// along with FlexibleSUSY. If not, see
16// <http://www.gnu.org/licenses/>.
17// ====================================================================
18
19#ifndef lattice_constraint_hpp
20#define lattice_constraint_hpp
21
22
23#include <Eigen/Dense>
24#include "constraint.hpp"
25#include "lattice_solver.hpp"
26
27namespace flexiblesusy {
28
29template <class T>
31
32class Lattice;
33
35public:
36 // Lattice_constraint() {}
37 virtual ~Lattice_constraint() = default;
38 virtual void init(RGFlow<Lattice> *flow) { f = flow; activate(); }
39 virtual void deactivate();
40 virtual void alloc_rows() = 0;
41 virtual void free_rows() { rfree(); }
42 virtual void operator()() = 0;
43 virtual void relocate(const std::vector<std::vector<size_t>>& site_maps)=0;
44 // size_t relocated_site(size_t m, size_t old_height, size_t new_height)
45 // { return size_t(0.5 + m*Real(new_height-1)/Real(old_height-1)); }
47protected:
48 virtual void activate();
49 Real& A(size_t r, size_t T, size_t m, size_t j)
50 { return f->A(rows[r]->rowSpec.realRow, T, m, j); }
51 Real y(size_t T, size_t m, size_t i) const { return f->y(T, m, i); }
52 Real& z(size_t r) { return f->z[rows[r]->rowSpec.realRow]; }
53 Real u(size_t T, size_t i) const { return f->u(T, i); }
54 Real x(size_t T, size_t m, size_t i) const { return f->x(T, m, i); }
55 void ralloc(size_t nrows, size_t T, size_t m, size_t span);
56 void rfree();
57private:
58 std::vector<RGFlow<Lattice>::EqRow *> rows;
59};
60
62
63template<>
64class Matching<Lattice> : public Lattice_constraint {
65public:
66 // InterTheoryConstraint() {}
67 virtual void init
68 (RGFlow<Lattice> *flow, size_t lower_theory)
69 { Lattice_constraint::init(flow); TL = lower_theory; }
70 virtual void relocate(const std::vector<std::vector<size_t>>&) {}
72protected:
73 Real& A(size_t r, size_t To, size_t j)
74 { return Lattice_constraint::A(r, TL+To, m(To), j); }
75 Real y(size_t To, size_t i) const
76 { return Lattice_constraint::y(TL+To, m(To), i); }
77 Real& z(size_t r) { return Lattice_constraint::z(r); }
78 Real u(size_t To, size_t i) const
79 { return Lattice_constraint::u(TL+To, i); }
80 Real x(size_t To, size_t i) const
81 { return Lattice_constraint::x(TL+To, m(To), i); }
82 size_t m(size_t To) const { return To ? 0 : f->efts[TL].height - 1; }
83 void ralloc(size_t nrows)
84 { Lattice_constraint::ralloc(nrows, TL, m(0), 2); }
85 size_t TL;
86private:
87};
88
90public:
91 // IntraTheoryConstraint() {}
92 virtual void init
93 (RGFlow<Lattice> *flow, size_t theory, size_t site, size_t span_) {
95 T = theory; mbegin = site; span = span_;
96 }
97 virtual void relocate(const std::vector<size_t>& site_map) {
98 size_t new_mbegin = site_map[mbegin];
99 size_t new_span = site_map[mbegin + span - 1] - new_mbegin + 1;
100 mbegin = new_mbegin;
101 span = new_span;
102 }
103 virtual void relocate(const std::vector<std::vector<size_t>>& site_maps)
104 { relocate(site_maps[T]); }
105 size_t mbegin;
107protected:
108 Real& A(size_t r, size_t m, size_t j)
109 { return Lattice_constraint::A(r, T, m, j); }
110 Real y(size_t m, size_t i) const { return Lattice_constraint::y(T, m, i); }
111 Real& z(size_t r) { return Lattice_constraint::z(r); }
112 Real u(size_t i) const { return Lattice_constraint::u(T, i); }
113 Real x(size_t m, size_t i) const { return Lattice_constraint::x(T, m, i); }
114 void ralloc(size_t nrows, size_t m, size_t span)
115 { Lattice_constraint::ralloc(nrows, T, m, span); }
116 size_t T;
117 size_t span;
118private:
119};
120
122
123template<>
124class Constraint<Lattice> : public IntraTheoryConstraint {
125public:
126 // SingleSiteConstraint() {}
127 virtual void init(RGFlow<Lattice> *flow, size_t theory, size_t site)
128 { IntraTheoryConstraint::init(flow, theory, site, 1); }
129 virtual double get_scale() const { return std::exp(x(0))*f->scl0; }
131protected:
132 Real& A(size_t r, size_t j)
133 { return IntraTheoryConstraint::A(r, mbegin, j); }
134 Real y(size_t i) const { return IntraTheoryConstraint::y(mbegin, i); }
135 Real& z(size_t r) { return IntraTheoryConstraint::z(r); }
136 Real x(size_t i) const { return IntraTheoryConstraint::x(mbegin, i); }
137 void ralloc(size_t nrows)
138 { IntraTheoryConstraint::ralloc(nrows, mbegin, 1); }
139private:
140};
141
143public:
145 (SingleSiteConstraint *c, size_t T_offset) :
146 ssc(c), To(T_offset) {}
147 virtual void init(RGFlow<Lattice> *flow, size_t lower_theory) {
148 InterTheoryConstraint::init(flow, lower_theory);
149 ssc->init(f, TL+To, m(To));
150 }
151 void alloc_rows() { ssc->alloc_rows(); }
152 void free_rows() { ssc->free_rows(); }
153 void operator()() { (*ssc)(); }
155private:
157 size_t To;
158};
159
161public:
163 (size_t nrows, std::function<void(AnySingleSiteConstraint *)> fxn) :
164 nr(nrows), fxn_(fxn)
165 {}
166 void alloc_rows() { ralloc(nr); }
167 void operator()() { fxn_(this); }
168private:
169 size_t nr;
170 std::function<void(AnySingleSiteConstraint *)> fxn_;
171};
172
174public:
175 // Lattice_RGE() {}
176 void init(RGFlow<Lattice> *flow, size_t theory, size_t site, size_t span) {
177 IntraTheoryConstraint::init(flow, theory, site, span);
178 x.resize(f->efts[T].w->width);
179 ddxm0i.resize(x.size());
180 ddxm1i.resize(x.size());
181 }
182 void alloc_rows() {
183 for (size_t n = 0; n < span - 1; n++)
184 ralloc(f->efts[T].w->width - 1, mbegin + n, 2);
185 }
186 void operator()();
188
189private:
193
194 void calc_dxmi_ddxmi(size_t m, size_t i, Real& dxmi, RVec& ddxmi);
195 void set_diff(size_t r, size_t m, size_t i);
196};
197
199public:
200 template<class A, class V, class M>
201 struct Adapter_ {
202 Adapter_() : x(nullptr,0), D(nullptr,0,0) {}
203 void set(A& xD, size_t width) {
204 v = &xD;
205 n = width;
206 new (&x) Eigen::Map<V>(v->data() , n);
207 new (&D) Eigen::Map<M>(v->data()+n, n, n);
208 }
209 size_t n;
210 Eigen::Map<V> x;
211 Eigen::Map<M> D;
212 A *v;
213 };
214
216 using const_Adapter = Adapter_<const Eigen::ArrayXd, const Eigen::VectorXd,
217 const Eigen::MatrixXd>;
218
219 // Lattice_RKRGE() {}
220 void init(RGFlow<Lattice> *flow, size_t theory, size_t site, size_t span) {
221 assert(span == 2);
222 IntraTheoryConstraint::init(flow, theory, site, span);
223 xD0.resize(f->efts[T].w->width * (1 + f->efts[T].w->width));
224 xD1.resize(xD0.size());
225 a0.set(xD0, f->efts[T].w->width);
226 a1.set(xD1, f->efts[T].w->width);
227 dx0.resize(f->efts[T].w->width);
228 dx1.resize(f->efts[T].w->width);
229 }
230 void alloc_rows() {
231 ralloc(f->efts[T].w->width - 1, mbegin, 2);
232 }
233 void operator()();
235
236private:
237 int evolve_to(Real t_new, Adapter& a, Real eps = -1);
238
239 Eigen::ArrayXd xD0, xD1;
242};
243
245public:
246 // Uniform_dt() {}
247 void alloc_rows() {
248 for (size_t n = 1; n < span - 1; n++)
249 ralloc(1, mbegin + n - 1, 3);
250 }
251 void operator()() {
252 for (size_t n = 1; n < span - 1; n++) {
253 size_t m = mbegin + n;
254 size_t r = n - 1;
255 A(r,m-1,0) = 1;
256 A(r,m ,0) = -2;
257 A(r,m+1,0) = 1;
258 z(r) = 0;
259 }
260 }
261private:
262};
263
265public:
266 // Match_t() {}
267 void alloc_rows() {
268 ralloc(1);
269 }
270 void operator()() {
271 A(0,0,0) = u(0,0);
272 A(0,1,0) = -u(1,0);
273 z(0) = 0;
274 }
275};
276
278public:
279 Unified_xi_xj(size_t i, size_t j) :
282 A(0,i) = u(i);
283 A(0,j) = -u(j);
284 z(0) = 0;
285 })
286 {}
287};
288
290public:
291 Fixed_x(size_t i, Real x) :
294 A(0,i) = u(i);
295 z(0) = x;
296 })
297 {}
298};
299
301public:
305 A(0,0) = u(0);
306 z(0) = std::log(mu/f->scl0);
307 })
308 {}
309};
310
311} // namespace flexiblesusy
312
313#endif // lattice_constraint_hpp
std::function< void(AnySingleSiteConstraint *)> fxn_
AnySingleSiteConstraint(size_t nrows, std::function< void(AnySingleSiteConstraint *)> fxn)
virtual void init(RGFlow< Lattice > *flow, size_t theory, size_t site)
virtual void init(RGFlow< Lattice > *flow, size_t theory, size_t site, size_t span_)
virtual void relocate(const std::vector< std::vector< size_t > > &site_maps)
virtual void relocate(const std::vector< size_t > &site_map)
void ralloc(size_t nrows, size_t m, size_t span)
Real & A(size_t r, size_t m, size_t j)
Real x(size_t m, size_t i) const
Real y(size_t m, size_t i) const
void init(RGFlow< Lattice > *flow, size_t theory, size_t site, size_t span)
void set_diff(size_t r, size_t m, size_t i)
void calc_dxmi_ddxmi(size_t m, size_t i, Real &dxmi, RVec &ddxmi)
Adapter_< Eigen::ArrayXd, Eigen::VectorXd, Eigen::MatrixXd > Adapter
void init(RGFlow< Lattice > *flow, size_t theory, size_t site, size_t span)
int evolve_to(Real t_new, Adapter &a, Real eps=-1)
Real x(size_t T, size_t m, size_t i) const
std::vector< RGFlow< Lattice >::EqRow * > rows
Real y(size_t T, size_t m, size_t i) const
Real u(size_t T, size_t i) const
virtual void relocate(const std::vector< std::vector< size_t > > &site_maps)=0
Real & A(size_t r, size_t T, size_t m, size_t j)
virtual void init(RGFlow< Lattice > *flow)
void ralloc(size_t nrows, size_t T, size_t m, size_t span)
virtual ~Lattice_constraint()=default
Real u(size_t To, size_t i) const
virtual void relocate(const std::vector< std::vector< size_t > > &)
virtual void init(RGFlow< Lattice > *flow, size_t lower_theory)
Real y(size_t To, size_t i) const
Real & A(size_t r, size_t To, size_t j)
Real x(size_t To, size_t i) const
Real & y(size_t T, size_t m, size_t i)
Real x(size_t T, size_t m, size_t i) const
Real & A(size_t r, size_t T, size_t m, size_t i)
Real u(size_t T, size_t i) const
SingleSiteInterTheoryConstraint(SingleSiteConstraint *c, size_t T_offset)
virtual void init(RGFlow< Lattice > *flow, size_t lower_theory)
std::complex< double > f(double tau) noexcept
double Real
Definition: mathdefs.hpp:12
std::vector< Real > RVec
Complex< T > log(const Complex< T > &z) noexcept
Definition: complex.hpp:54