LandauGinzburg
Loading...
Searching...
No Matches
field_potential.cpp
Go to the documentation of this file.
1
7#include "field_potential.hpp"
8
9
10MatrixXcd MatrixComp::ResizeMatrix (const MatrixXcd &m,
11 const int rw, const int cl) {
12 MatrixXcd res = MatrixXcd::Zero(rw, cl);
13 int r2 = m.rows(); int c2 = m.cols();
14 if (r2 >= rw) {
15 if (c2 >= cl)
16 res = m.block(0,0,rw,cl);
17 else
18 res.block(0,0,rw,c2) = m.block(0,0,rw,c2);
19 } else {
20 if (c2 >= cl)
21 res.block(0,0,r2,cl) = m.block(0,0,r2,cl);
22 else
23 res.block(0,0,r2,c2) = m;
24 }
25 return res;
26}
27
28MatrixXcd MatrixComp::DeleteZeroLine (const MatrixXcd &m) {
29 int r = m.rows();
30 int c = m.cols();
31 try {
32 if (r % 2 == 0 || c % 2 == 0)
33 throw "ERROR<DeleteZeroLine> Size of rows() or cols() is not ODD.";
34 int Qr = r / 2;
35 int Qc = c / 2;
36 MatrixXcd T(r - 1, c - 1);
37 T.block(0, 0, Qr, Qc) = m.block(0, 0, Qr, Qc);
38 T.block(0, Qc, Qr, Qc) = m.block(0, Qc+1, Qr, Qc);
39 T.block(Qr, 0, Qr, Qc) = m.block(Qr+1, 0, Qr, Qc);
40 T.block(Qr, Qc, Qr, Qc) = m.block(Qr+1, Qc+1, Qr, Qc);
41 return T;
42 } catch (const char *str) {cerr << str << endl; return MatrixXcd::Zero(r-1, c-1);}
43}
44double MatrixComp::DiffLU (const MatrixXcd &m,
45 const PartialPivLU<MatrixXcd> &lu) {
46 int N = m.rows();
47 try {
48 if (N != m.cols())
49 throw "ERROR<DiffLU> Not a square matrix.";
50 MatrixXcd P, L, U;
51 P = lu.permutationP();
52 L = lu.matrixLU().triangularView<StrictlyLower>();
53 L += MatrixXcd::Identity(N, N);
54 U = lu.matrixLU().triangularView<Upper>();
55 MatrixXcd Prod; Prod = P.inverse() * L * U;
56 MatrixXcd Diff; Diff = m - Prod;
57 return Diff.norm();
58 } catch (const char *str) {cerr << str << endl; return -1.0;}
59}
60int MatrixComp::SignDet (const MatrixXcd &m, const int debug) {
61 try {
62 if (m.rows() != m.cols())
63 throw "ERROR<SignDet> Not a square matrix.";
64 PartialPivLU<MatrixXcd> lu(m);
65 VectorXcd diag; diag = lu.matrixLU().diagonal();
66 std::complex<double> phase = 1;
67 for (int i = 0; i < diag.size(); ++i) {
68 phase *= diag(i) / std::abs(diag(i));
69 }
70 if ( lu.permutationP().determinant() < 0 ) phase *= -1;
71 int sign = 1;
72 if (phase.real() < 0) sign *= -1;
73 if (debug != 0) {
74 double diff = MatrixComp::DiffLU(m, lu);
75 VectorXd im = lu.matrixLU().diagonal().imag();
76 double im_norm = im.norm();
77 double im_max = 0;
78 for (int i = 0; i < im.size(); ++i) {
79 if (im_max < std::abs(im(i))) im_max = std::abs(im(i));
80 }
81 cerr << "PartialPiv: " << "Difference " << diff << endl
82 << " " << "Imag Norm " << im_norm << endl
83 << " " << " Max " << im_max << endl
84 << " " << "Phase of Det " << phase << endl;
85 }
86 return sign;
87 } catch (const char *str) {cerr << str << endl; return 0;}
88}
89
90
91
92void Potential::setpotential(const MatrixXcd &m) {
93 int N = m.rows(); int mcol = m.cols();
94 int num_field; int num_col;
95 for (int i = 1; i <= mcol; ++i) {
96 int tmp = int(i * (i + 1) / 2);
97 if (tmp >= mcol) {num_field = i; num_col = tmp; break;}
98 }
99 if (num_col == mcol) setconf(m, num_field);
100 else {
101 MatrixXcd mbig = MatrixXcd::Zero(N, num_col);
102 mbig.block(0,0, N, mcol) = m;
103 setconf(mbig, num_field);
104 }
105 return;
106}
107void Potential::setpotential(const MatrixXcd &m, const int num_field) {
108 if (num_field < 1) {setpotential(m); return;}
109 int num_col = int(num_field*(num_field+1)/2);
110 int mcol = m.cols();
111 try {
112 if (num_col > mcol) throw 1;
113 else if (num_col < mcol) throw 2;
114 setconf(m, num_field);
115 return;
116 } catch (const int i) {
117 int N = m.rows();
118 int num_field2; int num_col2;
119 switch (i) {
120 case 1:
121 num_field2 = num_field; num_col2 = num_col;
122 break;
123 case 2:
124 for (int j = num_field; j <= mcol; ++j) {
125 int tmp = int(j * (j+1) / 2);
126 if (tmp >= mcol) {num_field2 = j; num_col2 = tmp; break;}
127 }
128 break;
129 }
130 MatrixXcd mbig = MatrixXcd::Zero(N, num_col2);
131 mbig.block(0,0, N, mcol) = m;
132 setconf(mbig, num_field2);
133 return;
134 }
135}
136
137MatrixXcd Potential::unifield_tilde_type () const {
138 int N = (Li+1) * (Li+1);
139 try {
140 if (N % 2 == 0)
141 throw "ERROR<unifield_tilde_type> Size of Vector is not ODD.";
142 if (num_f != 1)
143 throw "ERROR<unifield_tilde_type> Number of superfields is not unity.";
144 MatrixXcd wm; wm = field2matrix(0) - dfield2matrix(0);
145 return MatrixComp::DeleteZeroLine(wm);
146 } catch (const char *str) {cerr << str << endl; return MatrixXcd::Zero(N-1, N-1);}
147}
148
149MatrixXcd Potential::unifield_jacobian () const {
150 int N = (Li+1) * (Li+1);
151 try {
152 if (N % 2 == 0)
153 throw "ERROR<unifield_jacobian> Size of Vector is not ODD.";
154 if (num_f != 1)
155 throw "ERROR<unifield_jacobian> Number of superfields is not unity.";
156 int L = Li+1;
157 MatrixXcd wt = unifield_tilde_type();
158 MatrixXcd p, pi;
159 p = MatrixXcd::Zero(N-1, N-1);
160 pi = MatrixXcd::Zero(N-1, N-1);
161 for (int i = 0; i < N/2; ++i) {
162 p(i, i) = 2.0 * M_PI * double(L-1)
163 * (double(L/2 - i/L) + IM * double(L/2 - i%L));
164 p(N-2 - i, N-2 - i) = - p(i, i);
165 pi(i, i) = p(i, i) / std::norm(p(i,i));
166 pi(N-2 - i, N-2 - i) = - pi(i,i);
167 }
168 return p + wt * pi * wt.adjoint();
169 } catch (const char *str) {cerr << str << endl; return MatrixXcd::Zero(N-1, N-1);}
170}
171MatrixXcd Potential::dualfield_jacobian () const {
172 int N = (Li+1) * (Li+1);
173 try {
174 if (N % 2 == 0)
175 throw "ERROR<dualfield_jacobian> Size of Vector is not ODD.";
176 if (num_f != 2)
177 throw "ERROR<dualfield_jacobian> Number of superfields is not 2.";
178 int L = Li+1;
179 MatrixXcd w12; w12 = field2matrix(1);
180 MatrixXcd m1 = MatrixXcd::Zero(2*N, 2*N);
181 m1.block(0,0, N,N) = w12; m1.block(N,N, N,N) = w12.adjoint();
182 MatrixXcd m2 = MatrixXcd::Zero(2*N, 2*N);
183 {
184 MatrixXcd p = MatrixXcd::Zero(N, N);
185 for (int i = 0; i < N; ++i) {
186 p(i, i) = 2.0 * M_PI * double(Li) * (IM * double(L/2 - i/L) + double(L/2 - i%L));
187 }
188 MatrixXcd m21 = MatrixXcd::Zero(2*N, 2*N);
189 {
190 MatrixXcd w22; w22 = field2matrix(2);
191 m21.block(0,0, N,N) = w22; m21.block(N,N, N,N) = w22.adjoint();
192 m21.block(N,0, N,N) = p; m21.block(0,N, N,N) = - p.adjoint();
193 }
194 MatrixXcd m22 = MatrixXcd::Zero(2*N, 2*N);
195 {
196 PartialPivLU<MatrixXcd> lu(w12);
197 MatrixXcd m_inv; m_inv = lu.inverse();
198 m22.block(0,0, N,N) = m_inv; m22.block(N,N, N,N) = m_inv.adjoint();
199 }
200 MatrixXcd m23 = MatrixXcd::Zero(2*N, 2*N);
201 {
202 MatrixXcd w11; w11 = field2matrix(0);
203 m23.block(0,0, N,N) = w11; m23.block(N,N, N,N) = w11.adjoint();
204 m23.block(N,0, N,N) = p; m23.block(0,N, N,N) = - p.adjoint();
205 }
206 m2 = m21 * m22 * m23;
207 }
208 return m1 - m2;
209 } catch (const char *str) {cerr << str << endl; return MatrixXcd::Zero(2*N,2*N);}
210}
211MatrixXcd Potential::triplefield_jacobian1 () const {
212 int N = (Li+1) * (Li+1);
213 try {
214 if (N % 2 == 0)
215 throw "ERROR<triplefield_jacobian> Size of Vector is not ODD.";
216 if (num_f != 3)
217 throw "ERROR<triplefield_jacobian> Number of superfields is not 3.";
218
219 MatrixXcd w11 = MatrixXcd::Zero(2*N, 2*N);
220 { MatrixXcd w; w = field2matrix(0);
221 w11.block(0,0, N,N) = w; w11.block(N,N, N,N) = w.adjoint();}
222 MatrixXcd w12 = MatrixXcd::Zero(2*N, 2*N);
223 { MatrixXcd w; w = field2matrix(1);
224 w12.block(0,0, N,N) = w; w12.block(N,N, N,N) = w.adjoint();}
225 MatrixXcd w13 = MatrixXcd::Zero(2*N, 2*N);
226 { MatrixXcd w; w = field2matrix(2);
227 w13.block(0,0, N,N) = w; w13.block(N,N, N,N) = w.adjoint();}
228 MatrixXcd w22 = MatrixXcd::Zero(2*N, 2*N);
229 { MatrixXcd w; w = field2matrix(3);
230 w22.block(0,0, N,N) = w; w22.block(N,N, N,N) = w.adjoint();}
231 MatrixXcd w23 = MatrixXcd::Zero(2*N, 2*N);
232 MatrixXcd w23_inv = MatrixXcd::Zero(2*N, 2*N);
233 { MatrixXcd w; w = field2matrix(4);
234 PartialPivLU<MatrixXcd> lu(w);
235 MatrixXcd w_inv; w_inv = lu.inverse();
236 w23.block(0,0, N,N) = w; w23.block(N,N, N,N) = w.adjoint();
237 w23_inv.block(0,0, N,N) = w_inv; w23_inv.block(N,N, N,N) = w_inv.adjoint();}
238 MatrixXcd w33 = MatrixXcd::Zero(2*N, 2*N);
239 { MatrixXcd w; w = field2matrix(5);
240 w33.block(0,0, N,N) = w; w33.block(N,N, N,N) = w.adjoint();}
241 {
242 int L = Li+1;
243 MatrixXcd p = MatrixXcd::Zero(N, N);
244 for (int i = 0; i < N; ++i) {
245 p(i, i) = 2.0 * M_PI * double(Li) * (IM * double(L/2 - i/L) + double(L/2 - i%L));
246 }
247 w11.block(N,0, N,N) = p; w11.block(0,N, N,N) = - p.adjoint();
248 w22.block(N,0, N,N) = p; w22.block(0,N, N,N) = - p.adjoint();
249 w33.block(N,0, N,N) = p; w33.block(0,N, N,N) = - p.adjoint();
250 }
251
252 MatrixXcd m1, m2;
253 m1 = w23 - w33 * w23_inv * w22;
254 {
255 MatrixXcd m2_1, m2_2, m2_3, m2_4;
256 m2_1 = w11 - w13 * w23_inv * w12;
257 m2_2 = w12 - w13 * w23_inv * w22;
258 {PartialPivLU<MatrixXcd> lu(m1); m2_3 = lu.inverse();}
259 m2_4 = w13 - w33 * w23_inv * w12;
260 m2 = m2_1 - m2_2 * m2_3 * m2_4;
261 }
262 return m1 * m2;
263 } catch (const char *str) {cerr << str << endl; return MatrixXcd::Zero(4*N,4*N);}
264}
265MatrixXcd Potential::triplefield_jacobian2 () const {
266 int N = (Li+1) * (Li+1);
267 try {
268 if (N % 2 == 0)
269 throw "ERROR<triplefield_jacobian> Size of Vector is not ODD.";
270 if (num_f != 3)
271 throw "ERROR<triplefield_jacobian> Number of superfields is not 3.";
272
273 MatrixXcd w11 = MatrixXcd::Zero(2*N, 2*N);
274 { MatrixXcd w; w = field2matrix(0);
275 w11.block(0,0, N,N) = w; w11.block(N,N, N,N) = w.adjoint();}
276 MatrixXcd w12 = MatrixXcd::Zero(2*N, 2*N);
277 MatrixXcd w12_inv = MatrixXcd::Zero(2*N, 2*N);
278 { MatrixXcd w; w = field2matrix(1);
279 PartialPivLU<MatrixXcd> lu(w);
280 MatrixXcd w_inv; w_inv = lu.inverse();
281 w12.block(0,0, N,N) = w; w12.block(N,N, N,N) = w.adjoint();
282 w12_inv.block(0,0, N,N) = w_inv; w12_inv.block(N,N, N,N) = w_inv.adjoint();}
283 MatrixXcd w13 = MatrixXcd::Zero(2*N, 2*N);
284 { MatrixXcd w; w = field2matrix(2);
285 w13.block(0,0, N,N) = w; w13.block(N,N, N,N) = w.adjoint();}
286 MatrixXcd w22 = MatrixXcd::Zero(2*N, 2*N);
287 { MatrixXcd w; w = field2matrix(3);
288 w22.block(0,0, N,N) = w; w22.block(N,N, N,N) = w.adjoint();}
289 MatrixXcd w23 = MatrixXcd::Zero(2*N, 2*N);
290 { MatrixXcd w; w = field2matrix(4);
291 w23.block(0,0, N,N) = w; w23.block(N,N, N,N) = w.adjoint();}
292 MatrixXcd w33 = MatrixXcd::Zero(2*N, 2*N);
293 { MatrixXcd w; w = field2matrix(5);
294 w33.block(0,0, N,N) = w; w33.block(N,N, N,N) = w.adjoint();}
295 {
296 int L = Li+1;
297 MatrixXcd p = MatrixXcd::Zero(N, N);
298 for (int i = 0; i < N; ++i) {
299 p(i, i) = 2.0 * M_PI * double(Li) * (IM * double(L/2 - i/L) + double(L/2 - i%L));
300 }
301 w11.block(N,0, N,N) = p; w11.block(0,N, N,N) = - p.adjoint();
302 w22.block(N,0, N,N) = p; w22.block(0,N, N,N) = - p.adjoint();
303 w33.block(N,0, N,N) = p; w33.block(0,N, N,N) = - p.adjoint();
304 }
305
306 MatrixXcd m1 = MatrixXcd::Zero(4*N, 4*N);
307 MatrixXcd m2 = MatrixXcd::Zero(4*N, 4*N);
308 {
309 m1.block(0, 0, 2*N,2*N) = w12; m1.block( 0,2*N, 2*N,2*N) = w13;
310 m1.block(2*N,0, 2*N,2*N) = w23; m1.block(2*N,2*N, 2*N,2*N) = w33;
311 }
312 {
313 MatrixXcd m2_1 = MatrixXcd::Zero(4*N, 2*N);
314 MatrixXcd m2_2 = MatrixXcd::Zero(2*N, 4*N);
315 m2_1.block(0,0, 2*N,2*N) = w11; m2_1.block(2*N,0, 2*N,2*N) = w13;
316 m2_2.block(0,0, 2*N,2*N) = w22; m2_2.block(0,2*N, 2*N,2*N) = w23;
317 m2 = m2_1 * w12_inv * m2_2;
318 }
319 return m1 - m2;
320 } catch (const char *str) {cerr << str << endl; return MatrixXcd::Zero(4*N,4*N);}
321}
322MatrixXcd Potential::multifield_jacobian () const {
323 int N = (Li+1) * (Li+1);
324 int Ns = 2 * N * num_f;
325 try {
326 if (N % 2 == 0)
327 throw "ERROR<multifield_jacobian> Size of Vector is not ODD.";
328 if (num_f < 1)
329 throw "ERROR<multifield_jacobian> Number of superfields is not positive.";
330 int L = Li+1;
331 MatrixXcd jac = MatrixXcd::Zero(Ns, Ns);
332 {
333 MatrixXcd iP = MatrixXcd::Zero(2*N, 2*N);
334 {
335 MatrixXcd p = MatrixXcd::Zero(N, N);
336 for (int i = 0; i < N/2; ++i) {
337 p(i, i) = 2.0 * M_PI * double(Li) * (IM * double(L/2 - i/L) + double(L/2 - i%L));
338 p(N-1 - i, N-1 - i) = - p(i,i);
339 }
340 iP.block(0,0, N,N) = p; iP.block(N,N, N,N) = - p.conjugate();
341 }
342 for (int i = 0; i < num_f; ++i) {
343 int tmp = 2*N * i;
344 jac.block(tmp, tmp, 2*N, 2*N) = iP;
345 }
346 }
347
348 int column_num = 0;
349 for (int i = 0; i < num_f; ++i) {
350 int tmp1 = 2*N * i;
351 for (int j = i; j < num_f; ++j) {
352 int tmp2 = 2*N * j;
353 MatrixXcd w; w = field2matrix(column_num);
354 jac.block(tmp1+N, tmp2, N,N) = w;
355 jac.block(tmp1, tmp2+N, N,N) = w.adjoint();
356 if (i != j) {
357 jac.block(tmp2+N, tmp1, N,N) = w;
358 jac.block(tmp2, tmp1+N, N,N) = w.adjoint();
359 }
360 ++column_num;
361 }
362 }
363 return jac;
364 } catch (const char *str) {cerr << str << endl; return MatrixXcd::Zero(Ns,Ns);}
365}
366
367MatrixXcd Potential::jacobian () const {
368 if (num_f == 1) {return unifield_jacobian();}
369 else if (num_f == 2) {return dualfield_jacobian();}
370 else if (num_f == 3) {return triplefield_jacobian2();}
371 return MatrixXcd::Zero( (Li+1)*(Li+1), (Li+1)*(Li+1) );
372}
373
374int Potential::sign_det (const int debug) const {
375 if (num_f < 4 && num_f > 0) {
376 MatrixXcd jac; jac = jacobian();
377 int sign = MatrixComp::SignDet(jac, debug);
378 if (num_f == 1) {return (Li%4==0)?(sign):(- sign);}
379 else if (num_f == 2) {return sign;}
380 else if (num_f == 3) {return sign;}
381 }
382 else if (num_f >=4) {
383 return sign_det4multifield(debug);
384 }
385 return 0;
386}
387int Potential::sign_det4multifield (const int debug) const {
388 MatrixXcd jac = multifield_jacobian();
389 int ph = (num_f % 2 == 0)?1:(-1);
390 return ph * MatrixComp::SignDet(jac, debug);
391}
392
393void Potential::show () const {
394 cout << "# class Potential" << endl
395 << "Li = " << Li << endl
396 << "num_f = " << num_f << endl
397 << "field(i, j) = " << endl
398 << field << endl;
399 return;
400}
401
402
403
404MatrixXd PotentialNR::nr_mat () const {
405 int L = Li + 1;
406 int N = L * L;
407
408 MatrixXd p0 = MatrixXd::Zero(N, N);
409 MatrixXd p1 = MatrixXd::Zero(N, N);
410 for (int i = 0; i < N; ++i) {
411 p0(i, i) = 2 * M_PI / double(Li) * double(L/2 - i/L);
412 p1(i, i) = 2 * M_PI / double(Li) * double(L/2 - i%L);
413 }
414
415 MatrixXd NR = MatrixXd::Zero(2*N * num_f, 2*N * num_f);
416 for (int i = 0; i < num_f; ++i) {
417 int tmp = 2*N * i;
418 //
419 // /p1 -p0\
420 // \p0 p1/
421 //
422 NR.block(tmp, tmp, N,N) = NR.block(tmp+N,tmp+N,N,N) = p1;
423 NR.block(tmp, tmp+N,N,N) = - (NR.block(tmp+N,tmp, N,N) = p0);
424 }
425 int column_num = 0;
426 for (int i = 0; i < num_f; ++i) {
427 int tmp1 = 2*N * i;
428 for (int j = i; j < num_f; ++j) {
429 int tmp2 = 2*N * j;
430 MatrixXcd w2f = field2matrix(column_num) / double( Li*Li );
431 MatrixXcd w2r = MatrixXcd::Zero(N, N);
432 for (int l = 0; l < N; ++l) {w2r.row(l) = w2f.row(N-1 - l);}
433 MatrixXd w2r_real = w2r.real();
434 MatrixXd w2r_imag = w2r.imag();
435 //
436 // /re -im\
437 // \-im -re/
438 //
439 NR.block(tmp1, tmp2, N,N) += w2r_real;
440 NR.block(tmp1, tmp2+N,N,N) -= w2r_imag;
441 NR.block(tmp1+N,tmp2, N,N) -= w2r_imag;
442 NR.block(tmp1+N,tmp2+N,N,N) -= w2r_real;
443 if (i != j) {
444 NR.block(tmp2, tmp1, N,N) += w2r_real;
445 NR.block(tmp2, tmp1+N,N,N) -= w2r_imag;
446 NR.block(tmp2+N,tmp1, N,N) -= w2r_imag;
447 NR.block(tmp2+N,tmp1+N,N,N) -= w2r_real;
448 }
449 ++column_num;
450 }
451 }
452 return NR;
453}
454VectorXd PotentialNR::nr_pvec () const {
455 int N = (Li+1) * (Li+1);
456 VectorXd NR = VectorXd::Zero(2*N * num_f);
457 for (int j = 0; j < num_f; ++j) {
458 int tmp = 2*N * j;
459 VectorXcd w1; w1 = field2.col(j).reverse();
460 VectorXcd w2; w2 = field3.col(j).reverse();
461 for (int i = 0; i < N; ++i) {
462 // NR(i) = double(k-2) * w1(i).real() + n(i).real();
463 // NR(i + N) = -double(k-2) * w1(i).imag() + n(i).imag();
464 NR(tmp + i) = w2(i).real() - w1(i).real();
465 NR(tmp + i + N) = - w2(i).imag() + w1(i).imag();
466 }
467 }
468 return NR;
469}
470
471VectorXcd PotentialNR::nrerr_pvec () const {
472 int N = (Li+1) * (Li+1);
473 int num_col = field2.cols();
474 VectorXcd NR = VectorXcd::Zero(N * num_col);
475 for (int j = 0; j < num_col; ++j) {
476 int tmp = N * j;
477 for (int i = 0; i < N; ++i) {NR(tmp + i) = field2(N-1 - i, j);}
478 }
479 return NR;
480}
481
482void PotentialNR::show () const {
483 cout << "# class PotentialNR" << endl
484 << "Li = " << Li << endl
485 << "num_f = " << num_f << endl
486 << "field(i, j)|_{i <= j} = " << endl
487 << field << endl
488 << "field2(i) = " << endl
489 << field2 << endl
490 << "field3(i) = " << endl
491 << field3 << endl;
492 return;
493}
494
Field::field
MatrixXcd field
Superfields.
Definition field_class.hpp:58
Field::field2matrix
MatrixXcd field2matrix(const int n) const
Definition field_class.cpp:154
Field::setconf
void setconf(const int num_col, const Distribution n, const double mean=0.0, const double dev=1.0)
Definition field_class.cpp:46
Field::Li
int Li
Physical box size, N_0=N_1.
Definition field_class.hpp:51
Field::dfield2matrix
MatrixXcd dfield2matrix(const int n) const
Definition field_class.cpp:175
Field::num_f
int num_f
Number of superfields.
Definition field_class.hpp:53
Potential::show
virtual void show() const
Output Li, field.
Definition field_potential.cpp:393
Potential::sign_det
int sign_det(const int debug=0) const
Sign determinant with a fast algorithm (num_f < 4)
Definition field_potential.cpp:374
Potential::sign_det4multifield
int sign_det4multifield(const int debug=0) const
Sign determinant for mutli-superfield.
Definition field_potential.cpp:387
PotentialNR::nr_pvec
VectorXd nr_pvec() const
Compute Vector for NR method (Real type)
Definition field_potential.cpp:454
PotentialNR::nr_mat
MatrixXd nr_mat() const
Compute Matrix for NR method (Real type)
Definition field_potential.cpp:404
PotentialNR::show
void show() const
Output Li, field, field2, field3.
Definition field_potential.cpp:482
PotentialNR::nrerr_pvec
VectorXcd nrerr_pvec() const
Compute Vector for NR error estimate (complex type)
Definition field_potential.cpp:471
IM
const std::complex< double > IM(0, 1.0)
Imaginary unit.
field_potential.hpp
Definition of the class Potential.
MatrixComp::DeleteZeroLine
MatrixXcd DeleteZeroLine(const MatrixXcd &m)
Delete a row and a column with p=0, q=0.
Definition field_potential.cpp:28
MatrixComp::DiffLU
double DiffLU(const MatrixXcd &m, const PartialPivLU< MatrixXcd > &lu)
Difference between a matrix and its LU decomposition.
Definition field_potential.cpp:44
MatrixComp::SignDet
int SignDet(const MatrixXcd &m, const int debug=0)
Sign of determinant with LU Decomposition.
Definition field_potential.cpp:60
MatrixComp::ResizeMatrix
MatrixXcd ResizeMatrix(const MatrixXcd &m, const int rw, const int cl)
Output m as MatrixXcd(rw, cl) with zeros.
Definition field_potential.cpp:10
Generated by doxygen 1.10.0