LORENE
vector_poisson.C
1/*
2 * Methods for solving vector Poisson equation.
3 *
4 * (see file vector.h for documentation).
5 *
6 */
7
8/*
9 * Copyright (c) 2003 Eric Gourgoulhon & Jerome Novak
10 *
11 * This file is part of LORENE.
12 *
13 * LORENE is free software; you can redistribute it and/or modify
14 * it under the terms of the GNU General Public License version 2
15 * as published by the Free Software Foundation.
16 *
17 * LORENE is distributed in the hope that it will be useful,
18 * but WITHOUT ANY WARRANTY; without even the implied warranty of
19 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
20 * GNU General Public License for more details.
21 *
22 * You should have received a copy of the GNU General Public License
23 * along with LORENE; if not, write to the Free Software
24 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
25 *
26 */
27
28char vector_poisson_C[] = "$Header: /cvsroot/Lorene/C++/Source/Tensor/vector_poisson.C,v 1.24 2014/10/13 08:53:45 j_novak Exp $" ;
29
30/*
31 * $Id: vector_poisson.C,v 1.24 2014/10/13 08:53:45 j_novak Exp $
32 * $Log: vector_poisson.C,v $
33 * Revision 1.24 2014/10/13 08:53:45 j_novak
34 * Lorene classes and functions now belong to the namespace Lorene.
35 *
36 * Revision 1.23 2014/10/06 15:13:21 j_novak
37 * Modified #include directives to use c++ syntax.
38 *
39 * Revision 1.22 2005/02/14 13:01:50 j_novak
40 * p_eta and p_mu are members of the class Vector. Most of associated functions
41 * have been moved from the class Vector_divfree to the class Vector.
42 *
43 * Revision 1.21 2005/01/10 15:36:09 j_novak
44 * In method 5: added transformation back from the Ylm base.
45 *
46 * Revision 1.20 2004/12/23 10:23:06 j_novak
47 * Improved method 5 in the case when some components are zero.
48 * Changed Vector_divfree::poisson to deduce eta from the equation.
49 *
50 * Revision 1.19 2004/08/24 09:14:50 p_grandclement
51 * Addition of some new operators, like Poisson in 2d... It now requieres the
52 * GSL library to work.
53 *
54 * Also, the way a variable change is stored by a Param_elliptic is changed and
55 * no longer uses Change_var but rather 2 Scalars. The codes using that feature
56 * will requiere some modification. (It should concern only the ones about monopoles)
57 *
58 * Revision 1.18 2004/07/27 09:40:05 j_novak
59 * Yet another method for solving vector Poisson eq. with spherical coordinates.
60 *
61 * Revision 1.17 2004/06/14 15:15:58 j_novak
62 * New method (No.4) to solve the vector Poisson eq. The output is continuous
63 * for all (spherical) components.
64 *
65 * Revision 1.16 2004/05/26 07:30:59 j_novak
66 * Version 1.15 was not good.
67 *
68 * Revision 1.15 2004/05/25 15:15:47 f_limousin
69 * Function Vector::poisson with parameters now returns a Vector (the
70 * result) instead of void.
71 *
72 * Revision 1.12 2004/03/26 17:05:24 j_novak
73 * Added new method n.3 using Tenseur::poisson_vect_oohara
74 *
75 * Revision 1.11 2004/03/11 08:48:45 f_limousin
76 * Implement method Vector::poisson with parameters, only with method
77 * 2 yet.
78 *
79 * Revision 1.10 2004/03/10 16:38:38 e_gourgoulhon
80 * Modified the prototype of poisson with param. to let it
81 * agree with declaration in vector.h.
82 *
83 * Revision 1.9 2004/03/03 09:07:03 j_novak
84 * In Vector::poisson(double, int), the flat metric is taken from the mapping.
85 *
86 * Revision 1.8 2004/02/24 17:00:25 j_novak
87 * Added a forgotten term.
88 *
89 * Revision 1.7 2004/02/24 09:46:20 j_novak
90 * Correction to cope with SGI compiler's warnings.
91 *
92 * Revision 1.6 2004/02/22 15:47:46 j_novak
93 * Added 2 more methods to solve the vector poisson equation. Method 1 is not
94 * tested yet.
95 *
96 * Revision 1.5 2004/02/20 10:53:41 j_novak
97 * Minor modifs.
98 *
99 * Revision 1.4 2004/02/16 17:40:14 j_novak
100 * Added a version of poisson with the flat metric as argument (avoids
101 * unnecessary calculations by decompose_div)
102 *
103 * Revision 1.3 2003/10/29 11:04:34 e_gourgoulhon
104 * dec2_dpzuis() replaced by dec_dzpuis(2).
105 * inc2_dpzuis() replaced by inc_dzpuis(2).
106 *
107 * Revision 1.2 2003/10/22 13:08:06 j_novak
108 * Better handling of dzpuis flags
109 *
110 * Revision 1.1 2003/10/20 15:15:42 j_novak
111 * New method Vector::poisson().
112 *
113 *
114 * $Headers: $
115 *
116 */
117
118//C headers
119#include <cstdlib>
120#include <cmath>
121
122// Lorene headers
123#include "metric.h"
124#include "tenseur.h"
125#include "param.h"
126#include "param_elliptic.h"
127#include "proto.h"
128
129namespace Lorene {
131 const {
132
133 bool nullite = true ;
134 for (int i=0; i<3; i++) {
135 assert(cmp[i]->check_dzpuis(4)) ;
136 if (cmp[i]->get_etat() != ETATZERO) nullite = false ;
137 }
138 assert ((method>=0) && (method<7)) ;
139
140 Vector resu(*mp, CON, triad) ;
141 if (nullite)
142 resu.set_etat_zero() ;
143 else {
144
145 switch (method) {
146
147 case 0 : {
148
149 Scalar poten(*mp) ;
150 if (fabs(lambda+1) < 1.e-8)
151 poten.set_etat_zero() ;
152 else {
153 poten = (potential(met_f) / (lambda + 1)).poisson() ;
154 }
155
156 Vector grad = poten.derive_con(met_f) ;
157 grad.dec_dzpuis(2) ;
158
159 return ( div_free(met_f).poisson() + grad) ;
160 break ;
161 }
162
163 case 1 : {
164
165 Scalar divf(*mp) ;
166 if (fabs(lambda+1) < 1.e-8)
167 divf.set_etat_zero() ;
168 else {
169 divf = (potential(met_f) / (lambda + 1)) ;
170 }
171
172 int nz = mp->get_mg()->get_nzone() ;
173
174 //-----------------------------------
175 // Removal of the l=0 part of div(F)
176 //-----------------------------------
178 div_lnot0.div_r_dzpuis(4) ;
179 Scalar source_r(*mp) ;
180 Valeur& va_div = div_lnot0.set_spectral_va() ;
181 if (div_lnot0.get_etat() != ETATZERO) {
182 va_div.coef() ;
183 va_div.ylm() ;
184 for (int lz=0; lz<nz; lz++) {
185 int np = mp->get_mg()->get_np(lz) ;
186 int nt = mp->get_mg()->get_nt(lz) ;
187 int nr = mp->get_mg()->get_nr(lz) ;
188 if (va_div.c_cf->operator()(lz).get_etat() != ETATZERO)
189 for (int k=0; k<np+1; k++)
190 for (int j=0; j<nt; j++) {
191 int l_q, m_q, base_r ;
192 if (nullite_plm(j, nt, k, np, va_div.base) == 1) {
193 donne_lm(nz, lz, j, k, va_div.base, m_q, l_q, base_r) ;
194 if (l_q == 0)
195 for (int i=0; i<nr; i++)
196 va_div.c_cf->set(lz, k, j, i) = 0. ;
197 }
198 }
199 }
200 source_r.set_etat_qcq() ;
201 source_r.set_spectral_va() = 2*(*va_div.c_cf) ; //2*div(F)
202 source_r.set_spectral_va().ylm_i() ;
203 source_r.set_dzpuis(4) ;
204 }
205 else
206 source_r.set_etat_zero() ;
207
208 //------------------------
209 // Other source terms ....
210 //------------------------
211 source_r += *(cmp[0]) - lambda*divf.dsdr() ;
212 Scalar f_r(*mp) ;
213 if (source_r.get_etat() != ETATZERO) {
214
215 //------------------------
216 // The elliptic operator
217 //------------------------
219 for (int lz=0; lz<nz; lz++)
220 param_fr.set_poisson_vect_r(lz) ;
221
222 f_r = source_r.sol_elliptic(param_fr) ;
223 }
224 else
225 f_r.set_etat_zero() ;
226
227 divf.dec_dzpuis(1) ;
228 Scalar source_eta = divf - f_r.dsdr() ;
229 source_eta.mult_r_dzpuis(0) ;
230 source_eta -= 2*f_r ;
231 resu.set_vr_eta_mu(f_r, source_eta.poisson_angu(), mu().poisson());
232
233 break ;
234
235 }
236
237 case 2 : {
238
239 Tenseur source_p(*mp, 1, CON, mp->get_bvect_spher() ) ;
240 source_p.set_etat_qcq() ;
241 for (int i=0; i<3; i++) {
242 source_p.set(i) = Cmp(*cmp[i]) ;
243 }
244 source_p.change_triad(mp->get_bvect_cart()) ;
245 Tenseur vect_auxi (*mp, 1, CON, mp->get_bvect_cart()) ;
246 vect_auxi.set_etat_qcq() ;
247 Tenseur scal_auxi (*mp) ;
248 scal_auxi.set_etat_qcq() ;
249
251 resu_p.change_triad(mp->get_bvect_spher() ) ;
252
253 for (int i=1; i<=3; i++)
254 resu.set(i) = resu_p(i-1) ;
255
256 break ;
257 }
258
259 case 3 : {
260
261 Tenseur source_p(*mp, 1, CON, mp->get_bvect_spher() ) ;
262 source_p.set_etat_qcq() ;
263 for (int i=0; i<3; i++) {
264 source_p.set(i) = Cmp(*cmp[i]) ;
265 }
266 source_p.change_triad(mp->get_bvect_cart()) ;
267 Tenseur scal_auxi (*mp) ;
268 scal_auxi.set_etat_qcq() ;
269
270 Tenseur resu_p(source_p.poisson_vect_oohara(lambda, scal_auxi)) ;
271 resu_p.change_triad(mp->get_bvect_spher() ) ;
272
273 for (int i=1; i<=3; i++)
274 resu.set(i) = resu_p(i-1) ;
275
276 break ;
277 }
278
279 case 4 : {
280 Scalar divf(*mp) ;
281 if (fabs(lambda+1) < 1.e-8)
282 divf.set_etat_zero() ;
283 else {
284 divf = (potential(met_f) / (lambda + 1)) ;
285 }
286
287 int nz = mp->get_mg()->get_nzone() ;
288
289 //-----------------------------------
290 // Removal of the l=0 part of div(F)
291 //-----------------------------------
293 div_tmp.div_r_dzpuis(4) ;
294 Scalar source_r(*mp) ;
295 Valeur& va_div = div_tmp.set_spectral_va() ;
296 if (div_tmp.get_etat() != ETATZERO) {
297 va_div.ylm() ;
298 for (int lz=0; lz<nz; lz++) {
299 int np = mp->get_mg()->get_np(lz) ;
300 int nt = mp->get_mg()->get_nt(lz) ;
301 int nr = mp->get_mg()->get_nr(lz) ;
302 if (va_div.c_cf->operator()(lz).get_etat() != ETATZERO)
303 for (int k=0; k<np+1; k++)
304 for (int j=0; j<nt; j++) {
305 int l_q, m_q, base_r ;
306 if (nullite_plm(j, nt, k, np, va_div.base) == 1) {
307 donne_lm(nz, lz, j, k, va_div.base, m_q, l_q, base_r) ;
308 if (l_q == 0)
309 for (int i=0; i<nr; i++)
310 va_div.c_cf->set(lz, k, j, i) = 0. ;
311 }
312 }
313 }
314 source_r.set_etat_qcq() ;
315 source_r.set_spectral_va() = 2*(*va_div.c_cf) ; //2*div(F)
316 source_r.set_spectral_va().ylm_i() ;
317 source_r.set_dzpuis(4) ;
318 }
319 else
320 source_r.set_etat_zero() ;
321
322 //------------------------
323 // Other source terms ....
324 //------------------------
325 source_r += *(cmp[0]) - lambda*divf.dsdr() ;
326 Scalar f_r(*mp) ;
327 if (source_r.get_etat() != ETATZERO) {
328
329 //------------------------
330 // The elliptic operator
331 //------------------------
333 for (int lz=0; lz<nz; lz++)
334 param_fr.set_poisson_vect_r(lz) ;
335
336 f_r = source_r.sol_elliptic(param_fr) ;
337 }
338 else
339 f_r.set_etat_zero() ;
340
341 //--------------------------
342 // Equation for eta
343 //--------------------------
344 Scalar source_eta = *cmp[1] ;
345 source_eta.mult_sint() ;
346 source_eta = source_eta.dsdt() ;
347 source_eta.div_sint() ;
348 source_eta = (source_eta + cmp[2]->stdsdp()).poisson_angu() ;
349
350 Scalar dfrsr = f_r.dsdr() ;
351 dfrsr.div_r_dzpuis(4) ;
352 Scalar frsr2 = f_r ;
353 frsr2.div_r_dzpuis(2) ;
354 frsr2.div_r_dzpuis(4) ;
355
356 Valeur& va_dfr = dfrsr.set_spectral_va() ;
357 Valeur& va_fsr = frsr2.set_spectral_va() ;
358 va_dfr.ylm() ;
359 va_fsr.ylm() ;
360
361 //Since the operator depends on the domain, the source
362 //must be transformed accordingly.
363 Valeur& va_eta = source_eta.set_spectral_va() ;
364 if (source_eta.get_etat() == ETATZERO) source_eta.annule_hard() ;
365 va_eta.ylm() ;
366 for (int lz=0; lz<nz; lz++) {
367 bool ced = (mp->get_mg()->get_type_r(lz) == UNSURR) ;
368 int np = mp->get_mg()->get_np(lz) ;
369 int nt = mp->get_mg()->get_nt(lz) ;
370 int nr = mp->get_mg()->get_nr(lz) ;
371 if (va_eta.c_cf->operator()(lz).get_etat() != ETATZERO)
372 for (int k=0; k<np+1; k++)
373 for (int j=0; j<nt; j++) {
374 int l_q, m_q, base_r ;
375 if (nullite_plm(j, nt, k, np, va_eta.base) == 1) {
376 donne_lm(nz, lz, j, k, va_eta.base, m_q, l_q, base_r) ;
377 if (l_q > 0)
378 for (int i=0; i<nr; i++) {
379 if (va_div.c_cf->operator()(lz).get_etat() != ETATZERO)
380 va_eta.c_cf->set(lz, k, j, i)
381 -= (lambda + 2. / double(ced ? -l_q : (l_q+1) ))
382 * va_div.c_cf->operator()(lz, k, j, i) ;
383 if (va_fsr.c_cf->operator()(lz).get_etat() != ETATZERO) {
384 va_eta.c_cf->set(lz, k, j, i)
385 += 2. / double(ced ? -l_q : (l_q+1) )
386 * va_dfr.c_cf->operator()(lz, k, j, i) ;
387 va_eta.c_cf->set(lz, k, j, i)
388 -= 2.*( ced ? double(l_q+2)/double(l_q)
389 : double(l_q-1)/double(l_q+1) )
390 * va_fsr.c_cf->set(lz, k, j, i) ;
391 }
392 } //Loop on r
393 } //nullite_plm
394 } //Loop on theta
395 } // Loop on domains
396 if (va_eta.c != 0x0) {
397 delete va_eta.c;
398 va_eta.c = 0x0 ;
399 }
400 va_eta.ylm_i() ;
401
402 //------------------------
403 // The elliptic operator
404 //------------------------
406 for (int lz=0; lz<nz; lz++)
407 param_eta.set_poisson_vect_eta(lz) ;
408
409 resu.set_vr_eta_mu(f_r, source_eta.sol_elliptic(param_eta), mu().poisson()) ;
410 break ;
411
412 }
413
414 case 5 : {
415
416 Scalar divf(*mp) ;
417 if (fabs(lambda+1) < 1.e-8)
418 divf.set_etat_zero() ;
419 else {
420 divf = (potential(met_f) / (lambda + 1)) ;
421 }
422
423 int nz = mp->get_mg()->get_nzone() ;
424
425 //-----------------------------------
426 // Removal of the l=0 part of div(F)
427 //-----------------------------------
429 div_lnot0.div_r_dzpuis(4) ;
430 Scalar source_r(*mp) ;
431 Valeur& va_div = div_lnot0.set_spectral_va() ;
432 if (div_lnot0.get_etat() != ETATZERO) {
433 va_div.coef() ;
434 va_div.ylm() ;
435 for (int lz=0; lz<nz; lz++) {
436 int np = mp->get_mg()->get_np(lz) ;
437 int nt = mp->get_mg()->get_nt(lz) ;
438 int nr = mp->get_mg()->get_nr(lz) ;
439 if (va_div.c_cf->operator()(lz).get_etat() != ETATZERO)
440 for (int k=0; k<np+1; k++)
441 for (int j=0; j<nt; j++) {
442 int l_q, m_q, base_r ;
443 if (nullite_plm(j, nt, k, np, va_div.base) == 1) {
444 donne_lm(nz, lz, j, k, va_div.base, m_q, l_q, base_r) ;
445 if (l_q == 0)
446 for (int i=0; i<nr; i++)
447 va_div.c_cf->set(lz, k, j, i) = 0. ;
448 }
449 }
450 }
451 source_r.set_etat_qcq() ;
452 source_r.set_spectral_va() = 2*(*va_div.c_cf) ; //2*div(F)
453 source_r.set_spectral_va().ylm_i() ;
454 source_r.set_dzpuis(4) ;
455 }
456 else
457 source_r.set_etat_zero() ;
458
459 //------------------------
460 // Other source terms ....
461 //------------------------
462 source_r += *(cmp[0]) - lambda*divf.dsdr() ;
463 Scalar f_r(*mp) ;
464 if (source_r.get_etat() != ETATZERO) {
465
466 //------------------------
467 // The elliptic operator
468 //------------------------
469
471 for (int lz=0; lz<nz; lz++)
472 param_fr.set_poisson_vect_r(lz) ;
473
474 f_r = source_r.sol_elliptic(param_fr) ;
475 }
476 else
477 f_r.set_etat_zero() ;
478
479 Scalar source_eta = - *(cmp[0]) ;
480 source_eta.mult_r_dzpuis(3) ;
481 source_eta -= (lambda+2.)*divf ;
482 source_eta.dec_dzpuis() ;
483 f_r.set_spectral_va().ylm() ;
484 Scalar tmp = 2*f_r + f_r.lapang() ;
485 tmp.div_r_dzpuis(2) ;
486 source_eta += tmp ;
487 tmp = (1.+lambda)*divf ;
488 tmp.mult_r_dzpuis(0) ;
489 tmp += f_r ;
490 source_eta = source_eta.primr() ;
491 f_r.set_spectral_va().ylm_i() ;
492 resu.set_vr_eta_mu(f_r, (tmp+source_eta).poisson_angu(), mu().poisson()) ;
493 break ;
494
495 }
496
497 case 6 : {
498
499 poisson_block(lambda, resu) ;
500 break ;
501
502 }
503 default : {
504 cout << "Vector::poisson : unexpected type of method !" << endl
505 << " method = " << method << endl ;
506 abort() ;
507 break ;
508 }
509
510 } // End of switch
511
512 } // End of non-null case
513
514 return resu ;
515
516}
517
519
520 const Base_vect_spher* tspher = dynamic_cast<const Base_vect_spher*>(triad) ;
521 const Base_vect_cart* tcart = dynamic_cast<const Base_vect_cart*>(triad) ;
522
523 assert ((tspher != 0x0) || (tcart != 0x0)) ;
524 const Metric_flat* met_f = 0x0 ;
525
526 if (tspher != 0x0) {
527 assert (tcart == 0x0) ;
528 met_f = &(mp->flat_met_spher()) ;
529 }
530
531 if (tcart != 0x0) {
532 assert (tspher == 0x0) ;
533 met_f = &(mp->flat_met_cart()) ;
534 }
535
536 return ( poisson(lambda, *met_f, method) );
537
538}
539
540// Version with parameters
541// -----------------------
542
543Vector Vector::poisson(const double lambda, Param& par, int method) const {
544
545
546 for (int i=0; i<3; i++)
547 assert(cmp[i]->check_dzpuis(4)) ;
548
549 assert ((method==0) || (method==2)) ;
550
551 Vector resu(*mp, CON, triad) ;
552
553 switch (method) {
554
555 case 0 : {
556
558 int nitermax = par.get_int(0) ;
559 int& niter = par.get_int_mod(0) ;
560 double relax = par.get_double(0) ;
561 double precis = par.get_double(1) ;
562 Cmp& ss_phi = par.get_cmp_mod(0) ;
563 Cmp& ss_khi = par.get_cmp_mod(1) ;
564 Cmp& ss_mu = par.get_cmp_mod(2) ;
565
566 Param par_phi ;
567 par_phi.add_int(nitermax, 0) ;
568 par_phi.add_int_mod(niter, 0) ;
569 par_phi.add_double(relax, 0) ;
570 par_phi.add_double(precis, 1) ;
571 par_phi.add_cmp_mod(ss_phi, 0) ;
572
573 Scalar poten(*mp) ;
574 if (fabs(lambda+1) < 1.e-8)
575 poten.set_etat_zero() ;
576 else {
577 Scalar tmp = potential(met_local) / (lambda + 1) ;
578 tmp.inc_dzpuis(2) ;
579 tmp.poisson(par_phi, poten) ;
580 }
581
582 Vector grad = poten.derive_con(met_local) ;
583 grad.dec_dzpuis(2) ;
584
586 par_free.add_int(nitermax, 0) ;
587 par_free.add_int_mod(niter, 0) ;
588 par_free.add_double(relax, 0) ;
589 par_free.add_double(precis, 1) ;
590 par_free.add_cmp_mod(ss_khi, 0) ;
591 par_free.add_cmp_mod(ss_mu, 1) ;
592
593 return div_free(met_local).poisson(par_free) + grad ;
594 break ;
595 }
596
597
598
599 case 2 : {
600
601 Tenseur source_p(*mp, 1, CON, mp->get_bvect_spher() ) ;
602 source_p.set_etat_qcq() ;
603 for (int i=0; i<3; i++) {
604 source_p.set(i) = Cmp(*cmp[i]) ;
605 }
606 source_p.change_triad(mp->get_bvect_cart()) ;
607
608 Tenseur vect_auxi (*mp, 1, CON, mp->get_bvect_cart()) ;
609 vect_auxi.set_etat_qcq() ;
610 for (int i=0; i<3 ;i++){
611 vect_auxi.set(i) = 0. ;
612 }
613 Tenseur scal_auxi (*mp) ;
614 scal_auxi.set_etat_qcq() ;
615 scal_auxi.set().annule_hard() ;
616 scal_auxi.set_std_base() ;
617
618 Tenseur resu_p(*mp, 1, CON, mp->get_bvect_cart() ) ;
619 resu_p.set_etat_qcq() ;
620 source_p.poisson_vect(lambda, par, resu_p, vect_auxi, scal_auxi) ;
621 resu_p.change_triad(mp->get_bvect_spher() ) ;
622
623 for (int i=1; i<=3; i++)
624 resu.set(i) = resu_p(i-1) ;
625
626 break ;
627 }
628
629 default : {
630 cout << "Vector::poisson : unexpected type of method !" << endl
631 << " method = " << method << endl ;
632
633
634 abort() ;
635 break ;
636 }
637
638 }// End of switch
639 return resu ;
640
641}
642
643
644
645
646}
Cartesian vectorial bases (triads).
Definition base_vect.h:201
Spherical orthonormal vectorial bases (triads).
Definition base_vect.h:308
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition cmp.h:446
Time evolution with partial storage (*** under development ***).
Definition evolution.h:371
Flat metric for tensor calculation.
Definition metric.h:261
This class contains the parameters needed to call the general elliptic solver.
Parameter storage.
Definition param.h:125
Tensor field of valence 0 (or component of a tensorial field).
Definition scalar.h:387
const Scalar & stdsdp() const
Returns of *this .
Tensor handling *** DEPRECATED : use class Tensor instead ***.
Definition tenseur.h:301
Values and coefficients of a (real-value) function.
Definition valeur.h:287
Tensor field of valence 1.
Definition vector.h:188
Vector poisson(double lambda, int method=6) const
Solves the vector Poisson equation with *this as a source.
virtual const Scalar & mu() const
Gives the field such that the angular components of the vector are written:
const Vector_divfree & div_free(const Metric &) const
Returns the div-free vector in the Helmholtz decomposition.
Definition vector.C:504
const Scalar & potential(const Metric &) const
Returns the potential in the Helmholtz decomposition.
Definition vector.C:492
const Map *const mp
Mapping on which the numerical values at the grid points are defined.
Definition tensor.h:295
Scalar ** cmp
Array of size n_comp of pointers onto the components.
Definition tensor.h:315
const Base_vect * triad
Vectorial basis (triad) with respect to which the tensor components are defined.
Definition tensor.h:303
Lorene prototypes.
Definition app_hor.h:64