LORENE
tenseur.h
1 /*
2  * Definition of Lorene classes Tenseur
3  * Tenseur_sym
4  *
5  */
6 
7 /*
8  * Copyright (c) 1999-2001 Philippe Grandclement
9  * Copyright (c) 2000-2001 Eric Gourgoulhon
10  * Copyright (c) 2002 Jerome Novak
11  *
12  * This file is part of LORENE.
13  *
14  * LORENE is free software; you can redistribute it and/or modify
15  * it under the terms of the GNU General Public License as published by
16  * the Free Software Foundation; either version 2 of the License, or
17  * (at your option) any later version.
18  *
19  * LORENE is distributed in the hope that it will be useful,
20  * but WITHOUT ANY WARRANTY; without even the implied warranty of
21  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
22  * GNU General Public License for more details.
23  *
24  * You should have received a copy of the GNU General Public License
25  * along with LORENE; if not, write to the Free Software
26  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
27  *
28  */
29 
30 
31 #ifndef __TENSEUR_H_
32 #define __TENSEUR_H_
33 
34 
35 /*
36  * $Id: tenseur.h,v 1.19 2016/09/19 15:26:22 j_novak Exp $
37  * $Log: tenseur.h,v $
38  * Revision 1.19 2016/09/19 15:26:22 j_novak
39  * Correction of several bugs preventing the shared library compilation.
40  *
41  * Revision 1.18 2014/10/13 08:52:37 j_novak
42  * Lorene classes and functions now belong to the namespace Lorene.
43  *
44  * Revision 1.17 2010/02/02 13:34:12 e_gourgoulhon
45  * Marked DEPRECATED (in the documentation).
46  *
47  * Revision 1.16 2005/08/30 08:35:10 p_grandclement
48  * Addition of the Tau version of the vectorial Poisson equation for the Tensors
49  *
50  * Revision 1.15 2004/12/29 16:24:03 k_taniguchi
51  * Addition of the method for vectorial Poisson equations with a multipole
52  * falloff condition at the outer boundary.
53  *
54  * Revision 1.14 2004/12/22 18:25:12 k_taniguchi
55  * Change an argument of poisson_vect_falloff
56  *
57  * Revision 1.13 2004/11/30 20:43:32 k_taniguchi
58  * Addition of the method for vectorial Poisson equations with falloff
59  * condition at the outer boundary.
60  *
61  * Revision 1.12 2004/03/22 13:12:43 j_novak
62  * Modification of comments to use doxygen instead of doc++
63  *
64  * Revision 1.11 2003/11/06 12:17:31 r_prix
65  * fixed mini-bug in documentation: without explicit argument in function-prototype,
66  * doc++ seemed to merge docu of Tensor::operator=(Cmp&) and operator=(Tenseur&)
67  *
68  * Revision 1.10 2003/06/20 14:23:38 f_limousin
69  * Add the functions compare().
70  *
71  * Revision 1.9 2002/10/16 14:36:29 j_novak
72  * Reorganization of #include instructions of standard C++, in order to
73  * use experimental version 3 of gcc.
74  *
75  * Revision 1.8 2002/09/19 14:12:37 e_gourgoulhon
76  * Modif documentation for LaTeX compliance.
77  *
78  * Revision 1.7 2002/09/10 13:44:17 j_novak
79  * The method "manipule" of one indice has been removed for Tenseur_sym objects
80  * (the result cannot be a Tenseur_sym).
81  * The method "sans_trace" now computes the traceless part of a Tenseur (or
82  * Tenseur_sym) of valence 2.
83  *
84  * Revision 1.6 2002/09/06 14:49:25 j_novak
85  * Added method lie_derive for Tenseur and Tenseur_sym.
86  * Corrected various errors for derive_cov and arithmetic.
87  *
88  * Revision 1.5 2002/08/14 13:46:14 j_novak
89  * Derived quantities of a Tenseur can now depend on several Metrique's
90  *
91  * Revision 1.4 2002/08/08 15:10:44 j_novak
92  * The flag "plat" has been added to the class Metrique to show flat metrics.
93  *
94  * Revision 1.3 2002/08/07 16:14:11 j_novak
95  * class Tenseur can now also handle tensor densities, this should be transparent to older codes
96  *
97  * Revision 1.2 2002/06/17 14:05:17 j_novak
98  * friend functions are now also declared outside the class definition
99  *
100  * Revision 1.1.1.1 2001/11/20 15:19:27 e_gourgoulhon
101  * LORENE
102  *
103  * Revision 2.46 2001/08/27 10:03:56 eric
104  * Ajout de l'operator% (produit tensoriel avec desaliasing)
105  *
106  * Revision 2.45 2001/06/19 15:38:00 eric
107  * Modif commentaires: mise en conformite Doc++ 3.4.8.
108  *
109  * Revision 2.44 2001/06/18 13:56:07 novak
110  * Ajout de la fonction abs()
111  *
112  * Revision 2.43 2001/05/29 16:10:43 eric
113  * Modif commentaires (mise en conformite Doc++ 3.4.7).
114  *
115  * Revision 2.42 2001/05/26 15:48:17 eric
116  * *** empty log message ***
117  *
118  * Revision 2.41 2001/05/26 15:42:54 eric
119  * Ajout de la fonction flat_scalar_prod_desal (desaliasage)
120  *
121  * Revision 2.40 2000/10/19 10:37:51 phil
122  * *** empty log message ***
123  *
124  * Revision 2.39 2000/10/19 09:47:15 phil
125  * *** empty log message ***
126  *
127  * Revision 2.38 2000/10/19 09:41:48 phil
128  * ajout de inverse_poisson_vect
129  *
130  * Revision 2.37 2000/10/06 12:37:09 keisuke
131  * Add a vectorial Poisson equation Tenseur::poisson_vect_regu.
132  *
133  * Revision 2.36 2000/09/27 08:52:34 eric
134  * Modif commentaires.
135  *
136  * Revision 2.35 2000/09/13 12:21:36 eric
137  * Modif commentaires.
138  *
139  * Revision 2.34 2000/09/13 12:11:26 eric
140  * Ajout de la fonction allocate_all().
141  *
142  * Revision 2.33 2000/05/22 14:39:11 phil
143  * ajout de inc_dzpuis et dec_dzpuis
144  *
145  * Revision 2.32 2000/04/03 15:18:54 phil
146  * suppression de poisson_vect_dirichlet
147  *
148  * Revision 2.31 2000/03/31 13:29:43 phil
149  * poisson_vect_neumann devient poisson_vect_dirichlet
150  *
151  * Revision 2.30 2000/03/30 16:10:46 phil
152  * *** empty log message ***
153  *
154  * Revision 2.29 2000/02/18 10:45:05 eric
155  * Modif commentaires.
156  *
157  * Revision 2.28 2000/02/11 19:11:55 phil
158  * commentaires
159  *
160  * Revision 2.27 2000/02/10 16:26:48 eric
161  * Modif commentaires.
162  *
163  * Revision 2.26 2000/02/10 16:10:49 eric
164  * Ajout de la fonction change_triad.
165  *
166  * Revision 2.25 2000/02/09 19:29:25 eric
167  * MODIF IMPORTANTE: la triade de decomposition est desormais passee en
168  * argument des constructeurs.
169  *
170  * Revision 2.24 2000/02/09 09:53:56 phil
171  * ajout de poisson_vect_oohara
172  * ,
173  *
174  * Revision 2.23 2000/02/08 19:04:25 eric
175  * Les fonctions arithmetiques ne sont plus amies.
176  * Ajout de nouvelles fonctions arithmetiques.
177  *
178  * Revision 2.22 2000/02/01 15:40:19 eric
179  * Ajout de la fonction sqrt
180  *
181  * Revision 2.21 2000/02/01 14:13:56 eric
182  * Modif commentaires.
183  * Ajout de la fonction amie flat_scalar_prod.
184  *
185  * Revision 2.20 2000/01/21 12:48:55 phil
186  * changement prototypage de Tenseur::poisson_vect
187  *
188  * Revision 2.19 2000/01/20 16:02:20 eric
189  * Ajout des operator=(double ) et operator=(int ).
190  *
191  * Revision 2.18 2000/01/20 14:52:08 phil
192  * *** empty log message ***
193  *
194  * Revision 2.17 2000/01/20 14:50:56 phil
195  * *** empty log message ***
196  *
197  * Revision 2.16 2000/01/20 14:39:31 phil
198  * *** empty log message ***
199  *
200  * Revision 2.15 2000/01/20 13:10:56 phil
201  * Ajout de Tenseur::poisson_vect (double)
202  *
203  * Revision 2.14 2000/01/20 11:20:17 phil
204  * changement prototypage
205  *
206  * Revision 2.13 2000/01/20 10:31:30 phil
207  * ajout de xksk
208  *
209  * Revision 2.12 2000/01/14 14:17:47 eric
210  * Modif commentaires.
211  *
212  * Revision 2.11 2000/01/14 14:04:07 eric
213  * Ajout de la fonction annule.
214  * classe Tenseur: constructeurs pour les classes derivees et fonctions
215  * de gestion memoire (del_t(), ...) declarees protected et non plus
216  * private.
217  *
218  * Revision 2.10 2000/01/13 14:15:30 eric
219  * Modif commentaires.
220  *
221  * Revision 2.9 2000/01/13 14:10:18 eric
222  * Ajout du constructeur par copie d'un Cmp (pour un scalaire)
223  * ainsi que l'affectation a un Cmp.
224  *
225  * Revision 2.8 2000/01/13 13:45:56 eric
226  * Ajout du membre p_gradient_spher et des fonctions fait_gradient_spher(),
227  * gradient_spher() pour le calcul du gradient d'un scalaire en
228  * coordonnees spheriques sur la triade spherique associee.
229  *
230  * Revision 2.7 2000/01/12 13:17:49 eric
231  * Les operator::(...) renvoie desormais une reference const sur le c[...]
232  * correspondant et non plus un Cmp copie de c[...].
233  * (ceci grace a Map::cmp_zero()).
234  *
235  * Revision 2.6 2000/01/11 11:13:16 eric
236  * Changement de nom pour la base vectorielle : base --> triad
237  *
238  * Revision 2.5 2000/01/10 17:22:26 eric
239  * Modif des #include
240  *
241  * Revision 2.4 2000/01/10 15:14:58 eric
242  * Ajout du membre base (base vectorielle sur laquelle sont definies
243  * les composantes).
244  *
245  * Revision 2.3 1999/12/09 12:39:39 phil
246  * changement prototypage des derivees
247  *
248  * Revision 2.2 1999/12/07 15:24:17 phil
249  * ajout include
250  *
251  * Revision 2.1 1999/12/03 09:37:11 phil
252  * *** empty log message ***
253  *
254  * Revision 2.0 1999/12/02 17:15:32 phil
255  * *** empty log message ***
256  *
257  * Revision 1.1 1999/12/02 17:13:29 phil
258  * Initial revision
259  *
260  *
261  * $Header: /cvsroot/Lorene/C++/Include/tenseur.h,v 1.19 2016/09/19 15:26:22 j_novak Exp $
262  *
263  */
264 
265 #define COV -1
266 #define CON +1
267 
268 #define N_MET_MAX 5
269 
270 // Headers Lorene
271 #include "cmp.h"
272 #include "itbl.h"
273 #include "base_vect.h"
274 
275 namespace Lorene {
276 class Metrique ;
277 class Tenseur_sym ;
278 
279  //---------------------------------//
280  // class Tenseur //
281  //---------------------------------//
282 
283 
301 class Tenseur {
302 
303  // Data :
304  // -----
305  protected:
306  const Map* const mp ;
307  int valence ;
308 
312  const Base_vect* triad ;
313 
319 
320  int n_comp ;
321  int etat ;
322  Cmp** c ;
323  double poids ;
324  const Metrique* metric ;
326 
327  // Derived data :
328  // ------------
329  protected:
337  const Metrique** met_depend ;
338 
343  mutable Tenseur* p_gradient ;
344 
351 
359 
366 
373 
374  // Constructors - Destructor :
375  // -------------------------
376  protected:
378  bool verif() const ;
379 
385  void new_der_met() ;
386 
387  public:
388  explicit Tenseur (const Map& map, const Metrique* met = 0x0,
389  double weight = 0) ;
390 
392  explicit Tenseur (const Cmp& cmp, const Metrique* met = 0x0,
393  double weight = 0) ;
394 
412  Tenseur (const Map& map, int val, const Itbl& tipe,
413  const Base_vect& triad_i, const Metrique* met = 0x0,
414  double weight = 0) ;
415 
434  Tenseur (const Map& map, int val, const Itbl& tipe,
435  const Base_vect* triad_i, const Metrique* met = 0x0,
436  double weight = 0) ;
437 
450  Tenseur (const Map& map, int val, int tipe, const
451  Base_vect& triad_i, const Metrique* met = 0x0,
452  double weight = 0) ;
453 
454  Tenseur (const Tenseur&) ;
455 
457  explicit Tenseur (const Tenseur_sym&) ;
458 
471  Tenseur (const Map& map, const Base_vect& triad_i, FILE* fich,
472  const Metrique* met = 0x0) ;
473 
483  Tenseur (const Map& map, FILE* fich, const Metrique* met = 0x0) ;
484 
485 
486  protected:
506  Tenseur (const Map& map, int val, const Itbl& tipe, int n_comp,
507  const Base_vect& triad_i, const Metrique* met = 0x0,
508  double weight = 0) ;
509 
524  Tenseur (const Map&, int val, int tipe, int n_comp,
525  const Base_vect& triad_i, const Metrique* met = 0x0,
526  double weight = 0) ;
527 
528 
529  public:
530 
531  virtual ~Tenseur() ;
532 
533  // Memory management
534  // -----------------
535  protected:
536  void del_t() ;
537 
542  void del_derive_met(int i) const ;
543 
547  void del_derive() const ;
548 
554  void set_der_met_0x0(int i) const ;
555 
560  void set_der_0x0() const ;
561 
562  // Mutators / assignment
563  // ---------------------
564  public:
569  void set_etat_nondef() ;
570 
575  void set_etat_zero() ;
576 
581  void set_etat_qcq() ;
582 
592  void allocate_all() ;
593 
597  void change_triad(const Base_vect& new_triad) ;
598 
605  void set_triad(const Base_vect& new_triad) ;
606  void set_poids(double weight) ;
607  void set_metric(const Metrique& met) ;
609 
611  virtual void operator=(const Tenseur& tens) ;
612 
614  void operator=(const Cmp& field) ;
615 
617  void operator=(double ) ;
618 
620  void operator=(int ) ;
621 
623  Cmp& set () ;
624  Cmp& set (int) ;
625  Cmp& set (int, int) ;
626  Cmp& set (int, int, int) ;
627  Cmp& set (const Itbl&) ;
628 
634  void annule(int l) ;
635 
647  void annule(int l_min, int l_max) ;
648 
654  void set_std_base() ;
655 
656  void dec_dzpuis() ;
657  void inc_dzpuis() ;
658  void dec2_dzpuis() ;
659  void inc2_dzpuis() ;
660  void mult_r_zec() ;
661 
666  Tenseur inverse_poisson_vect (double lambda) const ;
667 
668  // Accessors
669  // ---------
670  public:
681  virtual int donne_place (const Itbl& idx) const ;
682 
696  virtual Itbl donne_indices (int place) const ;
697 
699  const Map* get_mp() const {return mp ;} ;
700 
704  const Base_vect* get_triad() const {return triad;} ;
705 
707  int get_etat() const {return etat ;} ;
708 
710  int get_valence() const {return valence ; } ;
711 
713  int get_n_comp() const {return n_comp ;} ;
714 
726  int get_type_indice (int i) const {return type_indice(i) ;};
727 
735  Itbl get_type_indice () const {return type_indice ; } ;
736 
738  double get_poids() const {return poids ; } ;
739 
745  const Metrique* get_metric() const {return metric ; } ;
746 
747  const Cmp& operator()() const ;
748  const Cmp& operator()(int) const ;
749  const Cmp& operator()(int, int) const ;
750  const Cmp& operator()(int, int, int) const ;
751  const Cmp& operator()(const Itbl&) const ;
752 
753  // Outputs
754  // -------
755  public:
756  void sauve(FILE *) const ;
757  friend ostream& operator<<(ostream& , const Tenseur & ) ;
758 
759  // Computation of derived members
760  // ------------------------------
761  protected:
766  virtual void fait_gradient () const ;
767 
773  void fait_gradient_spher () const ;
774 
780  virtual void fait_derive_cov (const Metrique& met, int i) const ;
781 
787  virtual void fait_derive_con (const Metrique&, int i) const ;
788 
794  void fait_carre_scal (const Metrique&, int i) const ;
795 
805  void set_dependance (const Metrique& met) const ;
806 
812  int get_place_met(const Metrique& metre) const ;
813 
814  // Differential operators
815  // ----------------------
816  public:
818  const Tenseur& gradient() const ;
819 
823  const Tenseur& gradient_spher() const ;
824 
829  const Tenseur& derive_cov (const Metrique& met) const ;
830 
835  const Tenseur& derive_con (const Metrique&) const ;
836 
841  const Tenseur& carre_scal (const Metrique&) const ;
842 
843  // Resolution d'EDP :
869  void poisson_vect(double lambda, Param& par, Tenseur& shift, Tenseur& vect
870  , Tenseur& scal) const ;
871 
872  /*
873  * Same as poisson_vect with a Tau method
874  **/
875  void poisson_vect_tau(double lambda, Param& par, Tenseur& shift, Tenseur& vect
876  , Tenseur& scal) const ;
877 
878  void poisson_vect_falloff(double lambda, Param& par, Tenseur& shift,
879  Tenseur& vect, Tenseur& scal, int* k_falloff) const ;
880 
881  void poisson_vect_ylm(double lambda, Param& para, Tenseur& shift,
882  Tenseur& vecteur, Tenseur& scalaire, int nylm,
883  double* intvec) const ;
884 
912  Tenseur poisson_vect(double lambda, Tenseur& vect , Tenseur& scal ) const ;
913 
914  /*
915  * Same as poisson_vect with a Tau method
916  **/
917  Tenseur poisson_vect_tau(double lambda, Tenseur& vect , Tenseur& scal ) const ;
918 
919  Tenseur poisson_vect_falloff(double lambda, Tenseur& vect ,
920  Tenseur& scal, int* k_falloff ) const ;
921 
922  Tenseur poisson_vect_ylm(double lambda, Tenseur& vecteur,
923  Tenseur& scalaire, int nylm, double* intvec) const ;
924 
950  void poisson_vect_oohara(double lambda, Param& par, Tenseur& shift,
951  Tenseur& scal) const ;
952 
953  /*
954  * Same as poisson_vect_oohara with a Tau method
955  **/
956  void poisson_vect_oohara_tau(double lambda, Param& par, Tenseur& shift,
957  Tenseur& scal) const ;
958 
987  Tenseur poisson_vect_oohara(double lambda, Tenseur& scal) const ;
988  /*
989  * Same as poisson_vect_oohara with a Tau method
990  **/
991  Tenseur poisson_vect_oohara_tau(double lambda, Tenseur& scal) const ;
992 
1021  void poisson_vect_regu(int k_div, int nzet, double unsgam1,
1022  double lambda, Param& par, Tenseur& shift,
1023  Tenseur& vect, Tenseur& scal) const ;
1024 
1025 
1026  // Friend classes
1027  // ---------------
1028  friend class Tenseur_sym ;
1029  friend class Metrique ;
1030 
1031  // Mathematical operators
1032  // ----------------------
1033 
1034  friend Tenseur operator* (const Tenseur&, const Tenseur&) ;
1035  friend Tenseur operator% (const Tenseur&, const Tenseur&) ;
1036  friend Tenseur contract(const Tenseur&, int id1, int id2) ;
1037  friend Tenseur contract(const Tenseur&, int id1, const Tenseur&,
1038  int id2) ;
1039  friend Tenseur contract_desal(const Tenseur&, int id1, const Tenseur&,
1040  int id2) ;
1041  friend Tenseur flat_scalar_prod(const Tenseur& t1, const Tenseur& t2) ;
1042  friend Tenseur flat_scalar_prod_desal(const Tenseur& t1,
1043  const Tenseur& t2) ;
1044  friend Tenseur manipule(const Tenseur&, const Metrique&, int idx) ;
1045  friend Tenseur manipule(const Tenseur&, const Metrique&) ;
1046  friend Tenseur skxk (const Tenseur&) ;
1047  friend Tenseur lie_derive(const Tenseur& , const Tenseur& ,
1048  const Metrique* ) ;
1049 
1050 };
1051 
1052 
1059 Tenseur operator*(const Tenseur&, const Tenseur&) ;
1061 
1063 Tenseur operator%(const Tenseur&, const Tenseur&) ;
1064 
1087 Tenseur contract(const Tenseur&, int id1, int id2) ;
1088 
1111 Tenseur contract(const Tenseur&, int id1, const Tenseur&, int id2) ;
1112 
1119 Tenseur flat_scalar_prod(const Tenseur& t1, const Tenseur& t2) ;
1120 
1124 Tenseur flat_scalar_prod_desal(const Tenseur& t1, const Tenseur& t2) ;
1125 
1130 Tenseur manipule(const Tenseur&, const Metrique&, int idx) ;
1131 
1136 Tenseur manipule(const Tenseur&, const Metrique&) ;
1137 
1145 Tenseur skxk (const Tenseur&) ;
1146 
1153 Tenseur lie_derive (const Tenseur& t, const Tenseur& x, const Metrique* = 0x0);
1154 
1163 Tenseur sans_trace(const Tenseur& tens, const Metrique& metre) ;
1164 
1165 
1174 Tenseur operator+(const Tenseur& ) ;
1175 Tenseur operator-(const Tenseur& ) ;
1176 Tenseur operator+(const Tenseur&, const Tenseur &) ;
1177 
1179 Tenseur operator+(const Tenseur&, double ) ;
1180 
1182 Tenseur operator+(double, const Tenseur& ) ;
1183 
1185 Tenseur operator+(const Tenseur&, int ) ;
1186 
1188 Tenseur operator+(int, const Tenseur& ) ;
1189 
1190 Tenseur operator-(const Tenseur &, const Tenseur &) ;
1191 
1193 Tenseur operator-(const Tenseur&, double ) ;
1194 
1196 Tenseur operator-(double, const Tenseur& ) ;
1197 
1199 Tenseur operator-(const Tenseur&, int ) ;
1200 
1202 Tenseur operator-(int, const Tenseur& ) ;
1203 
1205 Tenseur operator*(const Tenseur&, double ) ;
1206 
1208 Tenseur operator*(double, const Tenseur& ) ;
1209 
1211 Tenseur operator*(const Tenseur&, int ) ;
1212 
1214 Tenseur operator*(int, const Tenseur& ) ;
1215 
1217 Tenseur operator/(const Tenseur& a, const Tenseur& b) ;
1218 
1219 Tenseur operator/(const Tenseur&, double ) ;
1220 
1222 Tenseur operator/(double, const Tenseur &) ;
1223 
1224 Tenseur operator/(const Tenseur&, int ) ;
1225 
1227 Tenseur operator/(int, const Tenseur &) ;
1228 
1229 Tenseur exp(const Tenseur& ) ;
1230 Tenseur log(const Tenseur& ) ;
1231 Tenseur sqrt(const Tenseur& ) ;
1232 Tenseur abs(const Tenseur& ) ;
1233 Tenseur pow(const Tenseur&, int ) ;
1234 Tenseur pow(const Tenseur&, double ) ;
1235 
1242  //---------------------------------//
1243  // class Tenseur_sym //
1244  //---------------------------------//
1245 
1253 class Tenseur_sym : public Tenseur {
1254 
1255  // Constructors - Destructor :
1256  // -------------------------
1257 
1258  public:
1273  Tenseur_sym (const Map& map, int val, const Itbl& tipe,
1274  const Base_vect& triad_i, const Metrique* met = 0x0,
1275  double weight = 0) ;
1276 
1286  Tenseur_sym (const Map& map, int val, int tipe,
1287  const Base_vect& triad_i, const Metrique* met = 0x0,
1288  double weight = 0) ;
1289 
1290  Tenseur_sym (const Tenseur_sym&) ;
1291 
1295  explicit Tenseur_sym (const Tenseur&) ;
1296 
1307  Tenseur_sym (const Map& map, const Base_vect& triad_i, FILE* fich,
1308  const Metrique* met = 0x0) ;
1309 
1310  virtual ~Tenseur_sym() ;
1311 
1312  // Mutators / assignment
1313  // ---------------------
1314  public:
1320  virtual void operator= (const Tenseur&) ;
1321 
1322 
1323  // Accessors
1324  // ---------
1325  public:
1336  virtual int donne_place (const Itbl& idx) const ;
1337 
1349  virtual Itbl donne_indices (int place) const ;
1350 
1351  // Computation of derived members
1352  // ------------------------------
1353  protected:
1358  virtual void fait_gradient () const ;
1359 
1365  virtual void fait_derive_cov (const Metrique& met, int i) const ;
1366 
1372  virtual void fait_derive_con (const Metrique&, int i) const ;
1373 
1374  // Mathematical operators
1375  // ----------------------
1376  friend Tenseur_sym operator* (const Tenseur&, const Tenseur_sym&) ;
1377  friend Tenseur_sym manipule(const Tenseur_sym&, const Metrique&) ;
1378  friend Tenseur lie_derive (const Tenseur& , const Tenseur& ,
1379  const Metrique* );
1380 
1381 } ;
1388 Tenseur_sym operator* (const Tenseur&, const Tenseur_sym&) ;
1390 
1395 Tenseur_sym manipule(const Tenseur_sym&, const Metrique&) ;
1396 
1404 Tenseur_sym lie_derive (const Tenseur_sym& t, const Tenseur& x,
1405  const Metrique* = 0x0);
1406 
1415 Tenseur_sym sans_trace(const Tenseur_sym& tens, const Metrique& metre) ;
1416 
1427 Tenseur_sym operator+(const Tenseur_sym& ) ;
1428 Tenseur_sym operator-(const Tenseur_sym& ) ;
1429 
1431 Tenseur_sym operator+(const Tenseur_sym&, const Tenseur_sym &) ;
1432 
1434 Tenseur_sym operator-(const Tenseur_sym &, const Tenseur_sym &) ;
1435 
1437 Tenseur_sym operator*(const Tenseur_sym&, double ) ;
1438 
1440 Tenseur_sym operator*(double, const Tenseur_sym& ) ;
1441 
1443 Tenseur_sym operator*(const Tenseur_sym&, int ) ;
1444 
1446 Tenseur_sym operator*(int, const Tenseur_sym& ) ;
1447 
1449 Tenseur_sym operator/(const Tenseur_sym& a, const Tenseur& b) ;
1450 
1451 Tenseur_sym operator/(const Tenseur_sym&, double ) ;
1452 
1453 Tenseur_sym operator/(const Tenseur_sym&, int ) ;
1454 
1461 }
1462 #endif
Lorene::Tenseur::operator%
friend Tenseur operator%(const Tenseur &, const Tenseur &)
Tensorial product with desaliasing.
Definition: tenseur_operateur.C:199
Lorene::Tenseur::carre_scal
const Tenseur & carre_scal(const Metrique &) const
Returns the scalar square of *this , with respect to met .
Definition: tenseur.C:1577
Lorene::Tenseur::etat
int etat
Logical state ETATZERO , ETATQCQ or ETATNONDEF.
Definition: tenseur.h:321
Lorene::Tenseur::set_etat_zero
void set_etat_zero()
Sets the logical state to ETATZERO (zero state).
Definition: tenseur.C:645
Lorene::Tenseur::poisson_vect_regu
void poisson_vect_regu(int k_div, int nzet, double unsgam1, double lambda, Param &par, Tenseur &shift, Tenseur &vect, Tenseur &scal) const
Solves the vectorial Poisson equation : .
Definition: tenseur_pde_regu.C:71
Lorene::Tenseur
Tensor handling *** DEPRECATED : use class Tensor instead ***.
Definition: tenseur.h:301
Lorene::Tenseur::new_der_met
void new_der_met()
Builds the arrays met_depend , p_derive_cov , p_derive_con and p_carre_scal and fills them with null ...
Definition: tenseur.C:196
Lorene::Tenseur::p_derive_cov
Tenseur ** p_derive_cov
Array of pointers on the covariant derivatives of *this with respect to the corresponding metric in...
Definition: tenseur.h:358
Lorene::Tenseur::set_metric
void set_metric(const Metrique &met)
Sets the pointer on the metric for a tensor density.
Definition: tenseur.C:685
Lorene::Tenseur::p_gradient_spher
Tenseur * p_gradient_spher
Pointer on the gradient of *this in a spherical orthonormal basis (scalar field only).
Definition: tenseur.h:350
Lorene::Tenseur::del_derive_met
void del_derive_met(int i) const
Logical destructor of the derivatives depending on the i-th element of *met_depend .
Definition: tenseur.C:554
Lorene::Tenseur::fait_gradient
virtual void fait_gradient() const
Calculates, if needed, the gradient of *this .
Definition: tenseur.C:1345
Lorene::Tenseur::inc_dzpuis
void inc_dzpuis()
dzpuis += 1 ;
Definition: tenseur.C:1117
Lorene::Tenseur::set
Cmp & set()
Read/write for a scalar (see also operator=(const Cmp&) ).
Definition: tenseur.C:824
Lorene::log
Cmp log(const Cmp &)
Neperian logarithm.
Definition: cmp_math.C:296
Lorene::Tenseur_sym::manipule
friend Tenseur_sym manipule(const Tenseur_sym &, const Metrique &)
Raise or lower all the indices, depending on their type, using the given Metrique .
Definition: tenseur_sym_operateur.C:137
Lorene::Tenseur::inc2_dzpuis
void inc2_dzpuis()
dzpuis += 2 ;
Definition: tenseur.C:1143
Lorene::exp
Cmp exp(const Cmp &)
Exponential.
Definition: cmp_math.C:270
Lorene
Lorene prototypes.
Definition: app_hor.h:64
Lorene::Tenseur::fait_gradient_spher
void fait_gradient_spher() const
Calculates, if needed, the gradient of *this in a spherical orthonormal basis (scalar field only).
Definition: tenseur.C:1392
Lorene::Tenseur::get_type_indice
int get_type_indice(int i) const
Returns the type of the index number i .
Definition: tenseur.h:726
Lorene::Tenseur::mult_r_zec
void mult_r_zec()
Multiplication by r in the external zone.
Definition: tenseur.C:1156
Lorene::Tenseur::get_metric
const Metrique * get_metric() const
Returns a pointer on the metric defining the conformal factor for tensor densities.
Definition: tenseur.h:745
Lorene::Tenseur::valence
int valence
Valence.
Definition: tenseur.h:307
Lorene::operator*
Base_val operator*(const Base_val &, const Base_val &)
This operator is used when calling multiplication or division of Valeur .
Definition: base_val_mult.C:102
Lorene::Tenseur::gradient
const Tenseur & gradient() const
Returns the gradient of *this (Cartesian coordinates)
Definition: tenseur.C:1542
Lorene::Tenseur::poisson_vect
void poisson_vect(double lambda, Param &par, Tenseur &shift, Tenseur &vect, Tenseur &scal) const
Solves the vectorial Poisson equation : .
Definition: tenseur_pde.C:118
Lorene::Tenseur::gradient_spher
const Tenseur & gradient_spher() const
Returns the gradient of *this (Spherical coordinates) (scalar field only).
Definition: tenseur.C:1548
Lorene::Tenseur::poids
double poids
For tensor densities: the weight.
Definition: tenseur.h:323
Lorene::Tenseur::type_indice
Itbl type_indice
Array of size valence contening the type of each index, COV for a covariant one and CON for a contrav...
Definition: tenseur.h:318
Lorene::Tenseur::sauve
void sauve(FILE *) const
Save in a file.
Definition: tenseur.C:1325
Lorene::Tenseur::p_derive_con
Tenseur ** p_derive_con
Array of pointers on the contravariant derivatives of *this with respect to the corresponding metri...
Definition: tenseur.h:365
Lorene::abs
Cmp abs(const Cmp &)
Absolute value.
Definition: cmp_math.C:410
Lorene::Tenseur::get_triad
const Base_vect * get_triad() const
Returns the vectorial basis (triad) on which the components are defined.
Definition: tenseur.h:704
Lorene::Itbl
Basic integer array class.
Definition: itbl.h:122
Lorene::Tenseur::set_triad
void set_triad(const Base_vect &new_triad)
Assigns a new vectorial basis (triad) of decomposition.
Definition: tenseur.C:674
Lorene::pow
Cmp pow(const Cmp &, int)
Power .
Definition: cmp_math.C:348
Lorene::Tenseur::skxk
friend Tenseur skxk(const Tenseur &)
Contraction of the last index of (*this) with or , depending on the type of S .
Definition: tenseur_operateur.C:580
Lorene::Tenseur::set_poids
void set_poids(double weight)
Sets the weight for a tensor density.
Definition: tenseur.C:680
Lorene::skxk
Tenseur skxk(const Tenseur &)
Contraction of the last index of (*this) with or , depending on the type of S .
Definition: tenseur_operateur.C:580
Lorene::manipule
Tenseur manipule(const Tenseur &, const Metrique &, int idx)
Raise or lower the index idx depending on its type, using the given Metrique .
Definition: tenseur_operateur.C:509
Lorene::Tenseur::get_etat
int get_etat() const
Returns the logical state.
Definition: tenseur.h:707
Lorene::Tenseur::contract
friend Tenseur contract(const Tenseur &, int id1, int id2)
Self contraction of two indices of a Tenseur .
Definition: tenseur_operateur.C:279
Lorene::flat_scalar_prod_desal
Tenseur flat_scalar_prod_desal(const Tenseur &t1, const Tenseur &t2)
Same as flat_scalar_prod but with desaliasing.
Definition: tenseur_operateur.C:735
Lorene::Tenseur::Tenseur
Tenseur(const Map &map, const Metrique *met=0x0, double weight=0)
Constructor for a scalar field.
Definition: tenseur.C:209
Lorene::Tenseur_sym::fait_derive_con
virtual void fait_derive_con(const Metrique &, int i) const
Calculates, if needed, the contravariant derivative of *this , with respect to met .
Definition: tenseur_sym.C:441
Lorene::Tenseur::set_dependance
void set_dependance(const Metrique &met) const
To be used to describe the fact that the derivatives members have been calculated with met .
Definition: tenseur.C:608
Lorene::Tenseur::inverse_poisson_vect
Tenseur inverse_poisson_vect(double lambda) const
Compute of *this , *this being of valence 1.
Definition: tenseur_inv_pois_vect.C:59
Lorene::Tenseur::get_valence
int get_valence() const
Returns the valence.
Definition: tenseur.h:710
Lorene::Tenseur::operator=
virtual void operator=(const Tenseur &tens)
Assignment to another Tenseur.
Definition: tenseur.C:718
Lorene::Tenseur::del_t
void del_t()
Logical destructor.
Definition: tenseur.C:545
Lorene::Tenseur::derive_con
const Tenseur & derive_con(const Metrique &) const
Returns the contravariant derivative of *this , with respect to met .
Definition: tenseur.C:1568
Lorene::lie_derive
Tenseur lie_derive(const Tenseur &t, const Tenseur &x, const Metrique *=0x0)
Lie Derivative of t with respect to x .
Definition: tenseur_operateur.C:816
Lorene::operator/
Cmp operator/(const Cmp &, const Cmp &)
Cmp / Cmp.
Definition: cmp_arithm.C:457
Lorene::Tenseur_sym
Class intended to describe tensors with a symmetry on the two last indices *** DEPRECATED : use class...
Definition: tenseur.h:1253
Lorene::Cmp
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition: cmp.h:446
Lorene::Tenseur::operator()
const Cmp & operator()() const
Read only for a scalar.
Definition: tenseur.C:943
Lorene::Tenseur::change_triad
void change_triad(const Base_vect &new_triad)
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.
Definition: tenseur.C:668
Lorene::operator-
Cmp operator-(const Cmp &)
- Cmp
Definition: cmp_arithm.C:108
Lorene::Tenseur::set_der_met_0x0
void set_der_met_0x0(int i) const
Sets the pointers of the derivatives depending on the i-th element of *met_depend to zero (as well as...
Definition: tenseur.C:583
Lorene::flat_scalar_prod
Tenseur flat_scalar_prod(const Tenseur &t1, const Tenseur &t2)
Scalar product of two Tenseur when the metric is : performs the contraction of the last index of t1 w...
Definition: tenseur_operateur.C:653
Lorene::Tenseur::p_gradient
Tenseur * p_gradient
Pointer on the gradient of *this .
Definition: tenseur.h:343
Lorene::Tenseur::derive_cov
const Tenseur & derive_cov(const Metrique &met) const
Returns the covariant derivative of *this , with respect to met .
Definition: tenseur.C:1554
Lorene::Tenseur::operator*
friend Tenseur operator*(const Tenseur &, const Tenseur &)
Tensorial product.
Definition: tenseur_operateur.C:119
Lorene::sans_trace
Tenseur sans_trace(const Tenseur &tens, const Metrique &metre)
Computes the traceless part of a Tenseur of valence 2.
Definition: tenseur_operateur.C:897
Lorene::Tenseur::verif
bool verif() const
Returns false for a tensor density without a defined metric.
Definition: tenseur.C:192
Lorene::sqrt
Cmp sqrt(const Cmp &)
Square root.
Definition: cmp_math.C:220
Lorene::Tenseur::fait_derive_cov
virtual void fait_derive_cov(const Metrique &met, int i) const
Calculates, if needed, the covariant derivative of *this , with respect to met .
Definition: tenseur.C:1429
Lorene::Tenseur_sym::~Tenseur_sym
virtual ~Tenseur_sym()
Destructor.
Definition: tenseur_sym.C:201
Lorene::Tenseur::p_carre_scal
Tenseur ** p_carre_scal
Array of pointers on the scalar squares of *this with respect to the corresponding metric in *met_d...
Definition: tenseur.h:372
Lorene::Tenseur::get_mp
const Map * get_mp() const
Returns pointer on the mapping.
Definition: tenseur.h:699
Lorene::operator%
Cmp operator%(const Cmp &, const Cmp &)
Cmp * Cmp with desaliasing.
Definition: cmp_arithm.C:364
Lorene::Tenseur::set_etat_qcq
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: tenseur.C:636
Lorene::Tenseur_sym::donne_indices
virtual Itbl donne_indices(int place) const
Returns the indices of a component given by its position in the Cmp 1-D array c .
Definition: tenseur_sym.C:252
Lorene::Tenseur::flat_scalar_prod_desal
friend Tenseur flat_scalar_prod_desal(const Tenseur &t1, const Tenseur &t2)
Same as flat_scalar_prod but with desaliasing.
Definition: tenseur_operateur.C:735
Lorene::Tenseur::del_derive
void del_derive() const
Logical destructor of all the derivatives.
Definition: tenseur.C:573
Lorene::Tenseur::fait_derive_con
virtual void fait_derive_con(const Metrique &, int i) const
Calculates, if needed, the contravariant derivative of *this , with respect to met .
Definition: tenseur.C:1499
Lorene::Tenseur::fait_carre_scal
void fait_carre_scal(const Metrique &, int i) const
Calculates, if needed, the scalar square of *this , with respect to met .
Definition: tenseur.C:1517
Lorene::Tenseur::mp
const Map *const mp
Reference mapping.
Definition: tenseur.h:306
Lorene::Tenseur::get_poids
double get_poids() const
Returns the weight.
Definition: tenseur.h:738
Lorene::Tenseur_sym::operator*
friend Tenseur_sym operator*(const Tenseur &, const Tenseur_sym &)
Tensorial product.
Definition: tenseur_sym_operateur.C:84
Lorene::Tenseur::~Tenseur
virtual ~Tenseur()
Destructor.
Definition: tenseur.C:533
Lorene::operator+
Cmp operator+(const Cmp &)
Definition: cmp_arithm.C:104
Lorene::Param
Parameter storage.
Definition: param.h:125
Lorene::Tenseur::annule
void annule(int l)
Sets the Tenseur to zero in a given domain.
Definition: tenseur.C:900
Lorene::Tenseur::get_place_met
int get_place_met(const Metrique &metre) const
Returns the position of the pointer on metre in the array met_depend .
Definition: tenseur.C:598
Lorene::Tenseur::set_etat_nondef
void set_etat_nondef()
Sets the logical state to ETATNONDEF (undefined state).
Definition: tenseur.C:650
Lorene::Tenseur_sym::operator=
virtual void operator=(const Tenseur &)
Assignment from a Tenseur .
Definition: tenseur_sym.C:283
Lorene::Tenseur::c
Cmp ** c
The components.
Definition: tenseur.h:322
Lorene::Tenseur_sym::fait_gradient
virtual void fait_gradient() const
Calculates, if needed, the gradient of *this .
Definition: tenseur_sym.C:325
Lorene::Tenseur_sym::fait_derive_cov
virtual void fait_derive_cov(const Metrique &met, int i) const
Calculates, if needed, the covariant derivative of *this , with respect to met .
Definition: tenseur_sym.C:372
Lorene::Tenseur::poisson_vect_oohara
void poisson_vect_oohara(double lambda, Param &par, Tenseur &shift, Tenseur &scal) const
Solves the vectorial Poisson equation .
Definition: tenseur_pde.C:218
Lorene::Tenseur::triad
const Base_vect * triad
Vectorial basis (triad) with respect to which the tensor components are defined.
Definition: tenseur.h:312
Lorene::Tenseur::met_depend
const Metrique ** met_depend
Array of pointers on the Metrique 's used to calculate derivatives members.
Definition: tenseur.h:337
Lorene::Tenseur::lie_derive
friend Tenseur lie_derive(const Tenseur &, const Tenseur &, const Metrique *)
Lie Derivative of t with respect to x .
Definition: tenseur_operateur.C:816
Lorene::Tenseur::donne_place
virtual int donne_place(const Itbl &idx) const
Returns the position in the Cmp 1-D array c of a component given by its indices.
Definition: tenseur.C:690
Lorene::Map
Base class for coordinate mappings.
Definition: map.h:670
Lorene::Tenseur_sym::donne_place
virtual int donne_place(const Itbl &idx) const
Returns the position in the Cmp 1-D array c of a component given by its indices.
Definition: tenseur_sym.C:207
Lorene::Tenseur::dec2_dzpuis
void dec2_dzpuis()
dzpuis -= 2 ;
Definition: tenseur.C:1130
Lorene::Tenseur::get_n_comp
int get_n_comp() const
Returns the number of components.
Definition: tenseur.h:713
Lorene::Tenseur::donne_indices
virtual Itbl donne_indices(int place) const
Returns the indices of a component given by its position in the Cmp 1-D array c .
Definition: tenseur.C:704
Lorene::Tenseur::set_der_0x0
void set_der_0x0() const
Sets the pointers of all the derivatives to zero.
Definition: tenseur.C:591
Lorene::Tenseur::metric
const Metrique * metric
For tensor densities: the metric defining the conformal factor.
Definition: tenseur.h:325
Lorene::Base_vect
Vectorial bases (triads) with respect to which the tensorial components are defined.
Definition: base_vect.h:105
Lorene::Tenseur::n_comp
int n_comp
Number of components, depending on the symmetry.
Definition: tenseur.h:320
Lorene::Tenseur::manipule
friend Tenseur manipule(const Tenseur &, const Metrique &, int idx)
Raise or lower the index idx depending on its type, using the given Metrique .
Definition: tenseur_operateur.C:509
Lorene::Tenseur::get_type_indice
Itbl get_type_indice() const
Returns the types of all the indices.
Definition: tenseur.h:735
Lorene::Tenseur_sym::lie_derive
friend Tenseur lie_derive(const Tenseur &, const Tenseur &, const Metrique *)
Lie Derivative of t with respect to x .
Definition: tenseur_operateur.C:816
Lorene::Tenseur::set_std_base
void set_std_base()
Set the standard spectal basis of decomposition for each component.
Definition: tenseur.C:1170
Lorene::Tenseur::allocate_all
void allocate_all()
Sets the logical state to ETATQCQ (ordinary state) and performs the memory allocation of all the elem...
Definition: tenseur.C:657
Lorene::contract
Tenseur contract(const Tenseur &, int id1, int id2)
Self contraction of two indices of a Tenseur .
Definition: tenseur_operateur.C:279
Lorene::Tenseur::dec_dzpuis
void dec_dzpuis()
dzpuis -= 1 ;
Definition: tenseur.C:1104
Lorene::Tenseur::flat_scalar_prod
friend Tenseur flat_scalar_prod(const Tenseur &t1, const Tenseur &t2)
Scalar product of two Tenseur when the metric is : performs the contraction of the last index of t1 w...
Definition: tenseur_operateur.C:653