LORENE
vector.h
1 /*
2 * Definition of Lorene class Vector
3 *
4 */
5
6/*
7 * Copyright (c) 2003 Eric Gourgoulhon & Jerome Novak
8 *
9 * This file is part of LORENE.
10 *
11 * LORENE is free software; you can redistribute it and/or modify
12 * it under the terms of the GNU General Public License version 2
13 * as published by the Free Software Foundation.
14 *
15 * LORENE is distributed in the hope that it will be useful,
16 * but WITHOUT ANY WARRANTY; without even the implied warranty of
17 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18 * GNU General Public License for more details.
19 *
20 * You should have received a copy of the GNU General Public License
21 * along with LORENE; if not, write to the Free Software
22 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
23 *
24 */
25
26#ifndef __VECTOR_H_
27#define __VECTOR_H_
28
29/*
30 * $Id: vector.h,v 1.42 2014/10/13 08:52:37 j_novak Exp $
31 * $Log: vector.h,v $
32 * Revision 1.42 2014/10/13 08:52:37 j_novak
33 * Lorene classes and functions now belong to the namespace Lorene.
34 *
35 * Revision 1.41 2008/12/03 10:18:56 j_novak
36 * Method 6 is now the default for calls to vector Poisson solver.
37 *
38 * Revision 1.40 2008/10/29 08:19:08 jl_cornou
39 * Typo in the doxygen documentation + spectral bases for pseudo vectors added
40 * and curl
41 *
42 * Revision 1.39 2008/08/27 08:45:21 jl_cornou
43 * Implemented routines to solve dirac systems for divergence-free vectors
44 *
45 * Revision 1.38 2007/12/21 16:06:16 j_novak
46 * Methods to filter Tensor, Vector and Sym_tensor objects.
47 *
48 * Revision 1.37 2007/05/11 11:38:16 n_vasset
49 * Added poisson_boundary2.C routine
50 *
51 * Revision 1.36 2007/01/16 15:05:59 n_vasset
52 * New constructor (taking a Scalar in mono-domain angular grid for
53 * boundary) for function sol_elliptic_boundary
54 *
55 * Revision 1.35 2005/06/09 07:56:25 f_limousin
56 * Implement a new function sol_elliptic_boundary() and
57 * Vector::poisson_boundary(...) which solve the vectorial poisson
58 * equation (method 6) with an inner boundary condition.
59 *
60 * Revision 1.34 2005/02/16 15:00:05 m_forot
61 * Add visu_streamlime function
62 *
63 * Revision 1.33 2005/02/14 13:01:48 j_novak
64 * p_eta and p_mu are members of the class Vector. Most of associated functions
65 * have been moved from the class Vector_divfree to the class Vector.
66 *
67 * Revision 1.32 2004/05/25 14:54:01 f_limousin
68 * Change arguments of method poisson with parameters.
69 *
70 * Revision 1.31 2004/05/09 20:54:22 e_gourgoulhon
71 * Added method flux (to compute the flux accross a sphere).
72 *
73 * Revision 1.30 2004/03/29 11:57:13 e_gourgoulhon
74 * Added methods ope_killing and ope_killing_conf.
75 *
76 * Revision 1.29 2004/03/28 21:25:14 e_gourgoulhon
77 * Minor modif comments (formula for V^\theta in Vector_divfree).
78 *
79 * Revision 1.28 2004/03/27 20:59:55 e_gourgoulhon
80 * Slight modif comment (doxygen \ingroup).
81 *
82 * Revision 1.27 2004/03/22 13:12:44 j_novak
83 * Modification of comments to use doxygen instead of doc++
84 *
85 * Revision 1.26 2004/03/10 12:52:19 f_limousin
86 * Add a new argument "method" in poisson method.
87 *
88 * Revision 1.25 2004/03/03 09:07:02 j_novak
89 * In Vector::poisson(double, int), the flat metric is taken from the mapping.
90 *
91 * Revision 1.24 2004/02/26 22:45:44 e_gourgoulhon
92 * Added method derive_lie.
93 *
94 * Revision 1.23 2004/02/22 15:47:45 j_novak
95 * Added 2 more methods to solve the vector Poisson equation. Method 1 is not
96 * tested yet.
97 *
98 * Revision 1.22 2004/02/21 16:27:53 j_novak
99 * Modif. comments
100 *
101 * Revision 1.21 2004/02/16 17:40:14 j_novak
102 * Added a version of poisson with the flat metric as argument (avoids
103 * unnecessary calculations by decompose_div)
104 *
105 * Revision 1.20 2003/12/14 21:47:24 e_gourgoulhon
106 * Added method visu_arrows for visualization through OpenDX.
107 *
108 * Revision 1.19 2003/11/06 14:43:37 e_gourgoulhon
109 * Gave a name to const arguments in certain method prototypes (e.g.
110 * constructors) to correct a bug of DOC++.
111 *
112 * Revision 1.18 2003/11/03 22:31:02 e_gourgoulhon
113 * Class Vector_divfree: parameters of methods set_vr_eta_mu and set_vr_mu
114 * are now all references.
115 *
116 * Revision 1.17 2003/11/03 14:02:17 e_gourgoulhon
117 * Class Vector_divfree: the members p_eta and p_mu are no longer declared
118 * "const".
119 *
120 * Revision 1.16 2003/10/20 15:12:17 j_novak
121 * New method Vector::poisson
122 *
123 * Revision 1.15 2003/10/20 14:44:49 e_gourgoulhon
124 * Class Vector_divfree: added method poisson().
125 *
126 * Revision 1.14 2003/10/20 09:32:10 j_novak
127 * Members p_potential and p_div_free of the Helmholtz decomposition
128 * + the method decompose_div(Metric).
129 *
130 * Revision 1.13 2003/10/17 16:36:53 e_gourgoulhon
131 * Class Vector_divfree: Added new methods set_vr_eta_mu and set_vr_mu.
132 *
133 * Revision 1.12 2003/10/16 21:36:01 e_gourgoulhon
134 * Added method Vector_divfree::update_vtvp().
135 *
136 * Revision 1.11 2003/10/16 15:25:00 e_gourgoulhon
137 * Changes in documentation.
138 *
139 * Revision 1.10 2003/10/16 14:21:33 j_novak
140 * The calculation of the divergence of a Tensor is now possible.
141 *
142 * Revision 1.9 2003/10/13 20:49:26 e_gourgoulhon
143 * Corrected typo in the comments.
144 *
145 * Revision 1.8 2003/10/13 13:52:39 j_novak
146 * Better managment of derived quantities.
147 *
148 * Revision 1.7 2003/10/12 20:33:15 e_gourgoulhon
149 * Added new derived class Vector_divfree (preliminary).
150 *
151 * Revision 1.6 2003/10/06 13:58:45 j_novak
152 * The memory management has been improved.
153 * Implementation of the covariant derivative with respect to the exact Tensor
154 * type.
155 *
156 * Revision 1.5 2003/10/05 21:07:27 e_gourgoulhon
157 * Method std_spectral_base() is now virtual.
158 *
159 * Revision 1.4 2003/10/03 14:08:03 e_gourgoulhon
160 * Added constructor from Tensor.
161 *
162 * Revision 1.3 2003/09/29 13:48:17 j_novak
163 * New class Delta.
164 *
165 * Revision 1.2 2003/09/29 12:52:56 j_novak
166 * Methods for changing the triad are implemented.
167 *
168 * Revision 1.1 2003/09/26 08:07:32 j_novak
169 * New class vector
170 *
171 *
172 * $Header: /cvsroot/Lorene/C++/Include/vector.h,v 1.42 2014/10/13 08:52:37 j_novak Exp $
173 *
174 */
175
176namespace Lorene {
177class Vector_divfree ;
178
179 //-------------------------//
180 // class Vector //
181 //-------------------------//
182
183
188class Vector: public Tensor {
189
190 // Derived data :
191 // ------------
192 protected:
198 mutable Scalar* p_potential[N_MET_MAX] ;
199
205 mutable Vector_divfree* p_div_free[N_MET_MAX] ;
206
219 mutable Scalar* p_eta ;
220
233 mutable Scalar* p_mu ;
234
241 mutable Scalar* p_A ;
242
243 // Constructors - Destructor
244 // -------------------------
245 public:
254 Vector(const Map& map, int tipe, const Base_vect& triad_i) ;
255
264 Vector(const Map& map, int tipe, const Base_vect* triad_i) ;
265
266 Vector(const Vector& a) ;
267
271 Vector(const Tensor& a) ;
272
283 Vector(const Map& map, const Base_vect& triad_i, FILE* fich) ;
284
285 virtual ~Vector() ;
286
287
288 // Memory management
289 // -----------------
290 protected:
291 virtual void del_deriv() const;
292
294 void set_der_0x0() const ;
295
300 virtual void del_derive_met(int) const ;
301
306 void set_der_met_0x0(int) const ;
307
308
309 // Mutators / assignment
310 // ---------------------
311 public:
315 virtual void change_triad(const Base_vect& ) ;
316
318 virtual void operator=(const Vector& a) ;
319
321 virtual void operator=(const Tensor& a) ;
322
332 void set_vr_eta_mu(const Scalar& vr_i, const Scalar& eta_i,
333 const Scalar& mu_i) ;
334
339 void decompose_div(const Metric&) const ;
340
348 const Scalar& potential(const Metric& ) const ;
349
357 const Vector_divfree& div_free(const Metric& ) const;
358
364 virtual void exponential_filter_r(int lzmin, int lzmax, int p,
365 double alpha= -16.) ;
366
372 virtual void exponential_filter_ylm(int lzmin, int lzmax, int p,
373 double alpha= -16.) ;
374
375
376 // Accessors
377 // ---------
378 public:
379 Scalar& set(int ) ;
380
381 const Scalar& operator()(int ) const;
382
392 virtual int position(const Itbl& idx) const {
393 assert (idx.get_ndim() == 1) ;
394 assert (idx.get_dim(0) == 1) ;
395 assert ((idx(0) >= 1) && (idx(0) <= 3)) ;
396
397 return (idx(0) - 1) ;
398
399 } ;
400
411 virtual Itbl indices(int place) const {
412 assert((place>=0) && (place<3)) ;
413
414 Itbl res(1) ;
415 res = place + 1;
416 return res ;
417
418 };
419
424 virtual void std_spectral_base() ;
425
430 virtual void pseudo_spectral_base() ;
431
444 virtual const Scalar& eta() const ;
445
458 virtual const Scalar& mu() const ;
459
460
467 virtual const Scalar& A() const ;
468
469
482 void update_vtvp() ;
483
484
485 // Differential operators/ PDE solvers
486 // -----------------------------------
487 public:
491 const Scalar& divergence(const Metric&) const ;
492
496 const Vector_divfree curl() const ;
497
501 Vector derive_lie(const Vector& v) const ;
502
510 Sym_tensor ope_killing(const Metric& gam) const ;
511
521 Sym_tensor ope_killing_conf(const Metric& gam) const ;
522
537 Vector poisson(double lambda, int method = 6) const ;
538
568 Vector poisson(double lambda, const Metric_flat& met_f, int method = 6) const ;
569
585 Vector poisson(const double lambda, Param& par,
586 int method = 6) const ;
587
600 void poisson_boundary(double lambda, const Mtbl_cf& limit_vr,
601 const Mtbl_cf& limit_eta, const Mtbl_cf& limit_mu,
602 int num_front, double fact_dir, double fact_neu,
603 Vector& resu) const ;
604
605
611 void poisson_boundary2(double lam, Vector& resu, Scalar boundvr,
612 Scalar boundeta, Scalar boundmu, double dir_vr,
613 double neum_vr, double dir_eta, double neum_eta,
614 double dir_mu, double neum_mu ) const ;
615
616
617 Vector poisson_dirichlet(double lambda, const Valeur& limit_vr,
618 const Valeur& limit_vt, const Valeur& limit_vp,
619 int num_front) const ;
620
633 Vector poisson_neumann(double lambda, const Valeur& limit_vr,
634 const Valeur& limit_vt, const Valeur& limit_vp,
635 int num_front) const ;
636
649 Vector poisson_robin(double lambda, const Valeur& limit_vr,
650 const Valeur& limit_vt, const Valeur& limit_vp,
651 double fact_dir, double fact_neu,
652 int num_front) const ;
653
654
668 double flux(double radius, const Metric& met) const ;
669
670 void poisson_block(double lambda, Vector& resu) const ;
671
672 // Graphics
673 // --------
674 public:
694 void visu_arrows(double xmin, double xmax, double ymin, double ymax,
695 double zmin, double zmax, const char* title0 = 0x0,
696 const char* filename0 = 0x0, bool start_dx = true, int nx = 8, int ny = 8,
697 int nz = 8) const ;
698
699 void visu_streamline(double xmin, double xmax, double ymin, double ymax,
700 double zmin, double zmax, const char* title0 = 0x0,
701 const char* filename0 = 0x0, bool start_dx = true, int nx = 8, int ny = 8,
702 int nz = 8) const ;
703
704
705
706};
707
708
709
710 //---------------------------------//
711 // class Vector_divfree //
712 //---------------------------------//
713
714
724class Vector_divfree: public Vector {
725
726 // Data :
727 // -----
728 protected:
730 const Metric* const met_div ;
731
732 // Constructors - Destructor
733 // -------------------------
734 public:
742 Vector_divfree(const Map& map, const Base_vect& triad_i,
743 const Metric& met) ;
744
746
758 Vector_divfree(const Map& map, const Base_vect& triad_i,
759 const Metric& met, FILE* fich) ;
760
761 virtual ~Vector_divfree() ;
762
763
764 // Memory management
765 // -----------------
766 protected:
767 virtual void del_deriv() const;
768
770 void set_der_0x0() const ;
771
772
773 // Mutators / assignment
774 // ---------------------
775
776 public:
778 void operator=(const Vector_divfree& a) ;
779
781 virtual void operator=(const Vector& a) ;
782
784 virtual void operator=(const Tensor& a) ;
785
797 void set_vr_mu(const Scalar& vr_i, const Scalar& mu_i) ;
798
807 void set_vr_eta_mu(const Scalar& vr_i, const Scalar& eta_i, const Scalar& mu_i) ;
808
816 void set_A_mu(const Scalar& A_i, const Scalar& mu_i, const Param* par_bc) ;
817
818 // Computational methods
819 // ---------------------
820 public:
833 virtual const Scalar& eta() const ;
834
844 Vector_divfree poisson() const ;
845
857 void update_etavr() ;
858
868 Vector_divfree poisson(Param& par) const ;
869 protected:
870
882 void sol_Dirac_A(const Scalar& aaa, Scalar& eta, Scalar& vr,
883 const Param* par_bc = 0x0) const ;
884
896 void sol_Dirac_A_tau(const Scalar& aaa, Scalar& eta, Scalar& vr,
897 const Param* par_bc = 0x0) const ;
898
910 void sol_Dirac_A_poisson(const Scalar& aaa, Scalar& eta, Scalar& vr,
911 const Param* par_bc = 0x0) const ;
912
924 void sol_Dirac_A_1z(const Scalar& aaa, Scalar& eta, Scalar& vr,
925 const Param* par_bc = 0x0) const ;
926
927};
928
929
930
931
932
933}
934#endif
Vectorial bases (triads) with respect to which the tensorial components are defined.
Definition base_vect.h:105
Basic integer array class.
Definition itbl.h:122
int get_dim(int i) const
Gives the i th dimension (ie {tt dim.dim[i] )
Definition itbl.h:326
int get_ndim() const
Gives the number of dimensions (ie dim.ndim )
Definition itbl.h:323
Base class for coordinate mappings.
Definition map.h:670
Flat metric for tensor calculation.
Definition metric.h:261
Metric for tensor calculation.
Definition metric.h:90
Coefficients storage for the multi-domain spectral method.
Definition mtbl_cf.h:186
Parameter storage.
Definition param.h:125
Tensor field of valence 0 (or component of a tensorial field).
Definition scalar.h:387
Class intended to describe valence-2 symmetric tensors.
Definition sym_tensor.h:223
Tensor handling.
Definition tensor.h:288
Values and coefficients of a (real-value) function.
Definition valeur.h:287
Divergence-free vectors.
Definition vector.h:724
void set_A_mu(const Scalar &A_i, const Scalar &mu_i, const Param *par_bc)
Defines the components through potentials and .
void set_vr_mu(const Scalar &vr_i, const Scalar &mu_i)
Sets the angular potentials (see member p_mu ), and the component of the vector.
Vector_divfree poisson() const
Computes the solution of a vectorial Poisson equation with *this as a source:
void sol_Dirac_A_1z(const Scalar &aaa, Scalar &eta, Scalar &vr, const Param *par_bc=0x0) const
Solves a one-domain system of two-coupled first-order PDEs obtained from the divergence-free conditio...
void update_etavr()
Computes the components and from the potential A and the divergence-free condition,...
void sol_Dirac_A_poisson(const Scalar &aaa, Scalar &eta, Scalar &vr, const Param *par_bc=0x0) const
Solves via a poisson method a system of two-coupled first-order PDEs obtained from the divergence-fre...
void set_vr_eta_mu(const Scalar &vr_i, const Scalar &eta_i, const Scalar &mu_i)
Defines the components through , and .
virtual void del_deriv() const
Deletes the derived quantities.
void sol_Dirac_A_tau(const Scalar &aaa, Scalar &eta, Scalar &vr, const Param *par_bc=0x0) const
Solves via a tau method a system of two-coupled first-order PDEs obtained from the divergence-free co...
const Metric *const met_div
Metric with respect to which the divergence is defined.
Definition vector.h:730
virtual ~Vector_divfree()
Destructor.
virtual const Scalar & eta() const
Gives the field such that the angular components of the vector are written:
void sol_Dirac_A(const Scalar &aaa, Scalar &eta, Scalar &vr, const Param *par_bc=0x0) const
Solves a system of two-coupled first-order PDEs obtained from the divergence-free condition and the r...
void operator=(const Vector_divfree &a)
Assignment from another Vector_divfree.
void set_der_0x0() const
Sets the pointers on derived quantities to 0x0.
Tensor field of valence 1.
Definition vector.h:188
const Scalar & operator()(int) const
Readonly access to a component.
Definition vector.C:305
Sym_tensor ope_killing(const Metric &gam) const
Computes the Killing operator associated with a given metric.
Definition vector.C:438
Vector poisson(double lambda, int method=6) const
Solves the vector Poisson equation with *this as a source.
virtual void std_spectral_base()
Sets the standard spectal bases of decomposition for each component.
Definition vector.C:316
void decompose_div(const Metric &) const
Makes the Helmholtz decomposition (see documentation of p_potential ) of this with respect to a given...
Definition vector.C:516
const Scalar & divergence(const Metric &) const
The divergence of this with respect to a Metric .
Definition vector.C:381
virtual Itbl indices(int place) const
Returns the index of a component given by its position in the Scalar array cmp .
Definition vector.h:411
void poisson_boundary(double lambda, const Mtbl_cf &limit_vr, const Mtbl_cf &limit_eta, const Mtbl_cf &limit_mu, int num_front, double fact_dir, double fact_neu, Vector &resu) const
Solves the vector Poisson equation with *this as a source with a boundary condition on the excised sp...
Scalar * p_potential[N_MET_MAX]
The potential giving the gradient part in the Helmholtz decomposition of any 3D vector .
Definition vector.h:198
virtual void del_deriv() const
Deletes the derived quantities.
Definition vector.C:219
virtual const Scalar & mu() const
Gives the field such that the angular components of the vector are written:
void update_vtvp()
Computes the components and from the potential and , according to:
virtual int position(const Itbl &idx) const
Returns the position in the Scalar array cmp of a component given by its index.
Definition vector.h:392
Vector poisson_neumann(double lambda, const Valeur &limit_vr, const Valeur &limit_vt, const Valeur &limit_vp, int num_front) const
Solves the vector Poisson equation with *this as a source with a boundary condition on the excised sp...
virtual void del_derive_met(int) const
Logical destructor of the derivatives depending on the i-th element of met_depend in the class Vector...
Definition vector.C:242
double flux(double radius, const Metric &met) const
Computes the flux of the vector accross a sphere r = const.
Definition vector.C:807
void visu_arrows(double xmin, double xmax, double ymin, double ymax, double zmin, double zmax, const char *title0=0x0, const char *filename0=0x0, bool start_dx=true, int nx=8, int ny=8, int nz=8) const
3D visualization via OpenDX.
Definition vector_visu.C:62
Vector_divfree * p_div_free[N_MET_MAX]
The divergence-free vector of the Helmholtz decomposition of any 3D vector .
Definition vector.h:205
const Vector_divfree curl() const
The curl of this with respect to a (flat) Metric .
Definition vector.C:402
Scalar * p_A
Field defined by.
Definition vector.h:241
virtual const Scalar & eta() const
Gives the field such that the angular components of the vector are written:
virtual void exponential_filter_r(int lzmin, int lzmax, int p, double alpha=-16.)
Applies exponential filters to all components (see Scalar::exponential_filter_r ).
Definition vector.C:850
const Vector_divfree & div_free(const Metric &) const
Returns the div-free vector in the Helmholtz decomposition.
Definition vector.C:504
void set_der_0x0() const
Sets the pointers on derived quantities to 0x0.
Definition vector.C:232
virtual ~Vector()
Destructor.
Definition vector.C:208
virtual void pseudo_spectral_base()
Sets the standard spectal bases of decomposition for each component for a pseudo_vector.
Definition vector.C:346
void set_der_met_0x0(int) const
Sets all the i-th components of met_depend in the class Vector (p_potential , etc....
Definition vector.C:258
Vector poisson_robin(double lambda, const Valeur &limit_vr, const Valeur &limit_vt, const Valeur &limit_vp, double fact_dir, double fact_neu, int num_front) const
Solves the vector Poisson equation with *this as a source with a boundary condition on the excised sp...
Scalar * p_mu
Field such that the angular components of the vector are written:
Definition vector.h:233
void poisson_boundary2(double lam, Vector &resu, Scalar boundvr, Scalar boundeta, Scalar boundmu, double dir_vr, double neum_vr, double dir_eta, double neum_eta, double dir_mu, double neum_mu) const
Alternative to previous poisson_boundary method for vectors ; this uses method 6 for vectorial solvin...
Scalar * p_eta
Field such that the angular components of the vector are written:
Definition vector.h:219
Sym_tensor ope_killing_conf(const Metric &gam) const
Computes the conformal Killing operator associated with a given metric.
Definition vector.C:462
virtual void change_triad(const Base_vect &)
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.
virtual void operator=(const Vector &a)
Assignment from a Vector.
Definition vector.C:267
virtual const Scalar & A() const
Gives the field defined by.
const Scalar & potential(const Metric &) const
Returns the potential in the Helmholtz decomposition.
Definition vector.C:492
Scalar & set(int)
Read/write access to a component.
Definition vector.C:296
void set_vr_eta_mu(const Scalar &vr_i, const Scalar &eta_i, const Scalar &mu_i)
Defines the components through potentials and (see members p_eta and p_mu ), as well as the compon...
virtual void exponential_filter_ylm(int lzmin, int lzmax, int p, double alpha=-16.)
Applies exponential filters to all components (see Scalar::exponential_filter_ylm ).
Definition vector.C:867
Vector derive_lie(const Vector &v) const
Computes the Lie derivative of this with respect to some vector field v.
Definition vector.C:392
Lorene prototypes.
Definition app_hor.h:64