LORENE
tenseur_sym_arithm.C
1 /*
2  * Arithmetics functions for the Tenseur_sym class.
3  *
4  * These functions are not member functions of the Tenseur_sym class.
5  *
6  * (see file tenseur.h for documentation).
7  *
8  */
9 
10 /*
11  * Copyright (c) 1999-2001 Philippe Grandclement
12  * Copyright (c) 2000-2001 Eric Gourgoulhon
13  *
14  * This file is part of LORENE.
15  *
16  * LORENE is free software; you can redistribute it and/or modify
17  * it under the terms of the GNU General Public License as published by
18  * the Free Software Foundation; either version 2 of the License, or
19  * (at your option) any later version.
20  *
21  * LORENE is distributed in the hope that it will be useful,
22  * but WITHOUT ANY WARRANTY; without even the implied warranty of
23  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
24  * GNU General Public License for more details.
25  *
26  * You should have received a copy of the GNU General Public License
27  * along with LORENE; if not, write to the Free Software
28  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
29  *
30  */
31 
32 
33 
34 
35 /*
36  * $Id: tenseur_sym_arithm.C,v 1.8 2016/12/05 16:18:17 j_novak Exp $
37  * $Log: tenseur_sym_arithm.C,v $
38  * Revision 1.8 2016/12/05 16:18:17 j_novak
39  * Suppression of some global variables (file names, loch, ...) to prevent redefinitions
40  *
41  * Revision 1.7 2014/10/13 08:53:42 j_novak
42  * Lorene classes and functions now belong to the namespace Lorene.
43  *
44  * Revision 1.6 2014/10/06 15:13:19 j_novak
45  * Modified #include directives to use c++ syntax.
46  *
47  * Revision 1.5 2003/06/20 14:54:17 f_limousin
48  * Put an assert on "poids" into comments
49  *
50  * Revision 1.4 2002/10/16 14:37:15 j_novak
51  * Reorganization of #include instructions of standard C++, in order to
52  * use experimental version 3 of gcc.
53  *
54  * Revision 1.3 2002/09/06 14:49:25 j_novak
55  * Added method lie_derive for Tenseur and Tenseur_sym.
56  * Corrected various errors for derive_cov and arithmetic.
57  *
58  * Revision 1.2 2002/08/07 16:14:11 j_novak
59  * class Tenseur can now also handle tensor densities, this should be transparent to older codes
60  *
61  * Revision 1.1.1.1 2001/11/20 15:19:30 e_gourgoulhon
62  * LORENE
63  *
64  * Revision 2.3 2000/02/09 19:30:36 eric
65  * MODIF IMPORTANTE: la triade de decomposition est desormais passee en
66  * argument des constructeurs.
67  *
68  * Revision 2.2 2000/02/08 19:06:40 eric
69  * Les fonctions arithmetiques ne sont plus amies.
70  * Modif de diverses operations (notament division avec double)
71  * Ajout de nouvelles operations (par ex. Tenseur + double, etc...)
72  *
73  * Revision 2.1 2000/01/11 11:15:00 eric
74  * Gestion de la base vectorielle (triad).
75  *
76  * Revision 2.0 1999/12/02 17:18:52 phil
77  * *** empty log message ***
78  *
79  *
80  * $Header: /cvsroot/Lorene/C++/Source/Tenseur/tenseur_sym_arithm.C,v 1.8 2016/12/05 16:18:17 j_novak Exp $
81  *
82  */
83 
84 // Headers C
85 #include <cstdlib>
86 #include <cassert>
87 #include <cmath>
88 
89 // Headers Lorene
90 #include "tenseur.h"
91 
92  //********************//
93  // OPERATEURS UNAIRES //
94  //********************//
95 
96 namespace Lorene {
98 
99  return t ;
100 
101 }
102 
103 
105 
106  assert (t.get_etat() != ETATNONDEF) ;
107  if (t.get_etat() == ETATZERO)
108  return t ;
109  else {
110  Tenseur_sym res(*(t.get_mp()), t.get_valence(), t.get_type_indice(),
111  *(t.get_triad()), t.get_metric(), t.get_poids() ) ;
112 
113  res.set_etat_qcq();
114  for (int i=0 ; i<res.get_n_comp() ; i++) {
115  Itbl indices (res.donne_indices(i)) ;
116  res.set(indices) = -t(indices) ;
117  }
118  return res ;
119  }
120 }
121 
122 
123  //**********//
124  // ADDITION //
125  //**********//
126 
127 Tenseur_sym operator+(const Tenseur_sym & t1, const Tenseur_sym & t2) {
128 
129  assert ((t1.get_etat() != ETATNONDEF) && (t2.get_etat() != ETATNONDEF)) ;
130  assert (t1.get_valence() == t2.get_valence()) ;
131  assert (t1.get_mp() == t2.get_mp()) ;
132  if (t1.get_valence() != 0) {
133  assert ( *(t1.get_triad()) == *(t2.get_triad()) ) ;
134  }
135 
136  for (int i=0 ; i<t1.get_valence() ; i++)
137  assert(t1.get_type_indice(i) == t2.get_type_indice(i)) ;
138  assert (t1.get_metric() == t2.get_metric()) ;
139  assert (fabs(t1.get_poids() - t2.get_poids())<1.e-10) ;
140 
141  if (t1.get_etat() == ETATZERO)
142  return t2 ;
143  else if (t2.get_etat() == ETATZERO)
144  return t1 ;
145  else {
146  Tenseur_sym res(*(t1.get_mp()), t1.get_valence(),
147  t1.get_type_indice(), *(t1.get_triad()),
148  t1.get_metric(), t1.get_poids() ) ;
149 
150  res.set_etat_qcq() ;
151  for (int i=0 ; i<res.get_n_comp() ; i++) {
152  Itbl indices (res.donne_indices(i)) ;
153  res.set(indices) = t1(indices) + t2(indices) ;
154  }
155  return res ;
156  }
157 }
158 
159 
160 
161  //**************//
162  // SOUSTRACTION //
163  //**************//
164 
165 
166 Tenseur_sym operator-(const Tenseur_sym & t1, const Tenseur_sym & t2) {
167 
168  return (t1 + (-t2)) ;
169 
170 }
171 
172 
173  //****************//
174  // MULTIPLICATION //
175  //****************//
176 
177 Tenseur_sym operator*(double x, const Tenseur_sym& t) {
178 
179  assert (t.get_etat() != ETATNONDEF) ;
180  if ( (t.get_etat() == ETATZERO) || (x == double(1)) )
181  return t ;
182  else {
183  Tenseur_sym res(*(t.get_mp()), t.get_valence(), t.get_type_indice(),
184  *(t.get_triad()), t.get_metric(), t.get_poids() ) ;
185 
186  if ( x == double(0) )
187  res.set_etat_zero() ;
188  else {
189  res.set_etat_qcq() ;
190  for (int i=0 ; i<res.get_n_comp() ; i++) {
191  Itbl indices (res.donne_indices(i)) ;
192  res.set(indices) = x*t(indices) ;
193  }
194  }
195  return res ;
196  }
197 }
198 
199 
200 Tenseur_sym operator* (const Tenseur_sym& t, double x) {
201  return x * t ;
202 }
203 
205  return double(m) * t ;
206 }
207 
208 
210  return double(m) * t ;
211 }
212 
213 
214 
215  //**********//
216  // DIVISION //
217  //**********//
218 
219 Tenseur_sym operator/ (const Tenseur_sym& t1, const Tenseur& t2) {
220 
221  // Protections
222  assert(t1.get_etat() != ETATNONDEF) ;
223  assert(t2.get_etat() != ETATNONDEF) ;
224  assert(t2.get_valence() == 0) ; // t2 doit etre un scalaire !
225  assert(t1.get_mp() == t2.get_mp()) ;
226 
227  double poids_res = t1.get_poids() - t2.get_poids() ;
228  poids_res = (fabs(poids_res) < 1.e-10 ? 0. : poids_res) ;
229  const Metrique* met_res = 0x0 ;
230  if (poids_res != 0.) {
231  // assert((t1.get_metric() != 0x0) || (t2.get_metric() != 0x0)) ;
232  if (t1.get_metric() != 0x0) met_res = t1.get_metric() ;
233  else met_res = t2.get_metric() ;
234  }
235 
236  // Cas particuliers
237  if (t2.get_etat() == ETATZERO) {
238  cout << "Division by 0 in Tenseur_sym / Tenseur !" << endl ;
239  abort() ;
240  }
241  if (t1.get_etat() == ETATZERO) {
242  Tenseur_sym resu(t1) ;
243  resu.set_poids(poids_res) ;
244  resu.set_metric(*met_res) ;
245  return resu ;
246  }
247 
248  // Cas general
249 
250  assert(t1.get_etat() == ETATQCQ) ; // sinon...
251  assert(t2.get_etat() == ETATQCQ) ; // sinon...
252 
253  Tenseur_sym res(*(t1.get_mp()), t1.get_valence(), t1.get_type_indice(),
254  *(t1.get_triad()), met_res, poids_res ) ;
255 
256  res.set_etat_qcq() ;
257  for (int i=0 ; i<res.get_n_comp() ; i++) {
258  Itbl indices (res.donne_indices(i)) ;
259  res.set(indices) = t1(indices) / t2() ; // Cmp / Cmp
260  }
261  return res ;
262 
263 }
264 
265 
266 Tenseur_sym operator/ (const Tenseur_sym& t, double x) {
267 
268  assert (t.get_etat() != ETATNONDEF) ;
269 
270  if ( x == double(0) ) {
271  cout << "Division by 0 in Tenseur_sym / double !" << endl ;
272  abort() ;
273  }
274 
275  if ( (t.get_etat() == ETATZERO) || (x == double(1)) )
276  return t ;
277  else {
278  Tenseur_sym res(*(t.get_mp()), t.get_valence(), t.get_type_indice(),
279  *(t.get_triad()), t.get_metric(), t.get_poids() ) ;
280 
281  res.set_etat_qcq() ;
282  for (int i=0 ; i<res.get_n_comp() ; i++) {
283  Itbl indices (res.donne_indices(i)) ;
284  res.set(indices) = t(indices) / x ; // Cmp / double
285  }
286  return res ;
287  }
288 }
289 
290 
291 
293 
294  return t / double(m) ;
295 }
296 
297 
298 }
double get_poids() const
Returns the weight.
Definition: tenseur.h:741
int get_type_indice(int i) const
Returns the type of the index number i .
Definition: tenseur.h:729
int & set(int i)
Read/write of a particular element (index i ) (1D case)
Definition: itbl.h:247
void set_poids(double weight)
Sets the weight for a tensor density.
Definition: tenseur.C:696
Class intended to describe tensors with a symmetry on the two last indices *** DEPRECATED : use class...
Definition: tenseur.h:1256
Lorene prototypes.
Definition: app_hor.h:67
Base_val operator*(const Base_val &, const Base_val &)
This operator is used when calling multiplication or division of Valeur .
Basic integer array class.
Definition: itbl.h:122
Cmp operator/(const Cmp &, const Cmp &)
Cmp / Cmp.
Definition: cmp_arithm.C:460
int get_valence() const
Returns the valence.
Definition: tenseur.h:713
void set_metric(const Metrique &met)
Sets the pointer on the metric for a tensor density.
Definition: tenseur.C:701
const Map * get_mp() const
Returns pointer on the mapping.
Definition: tenseur.h:702
const Metrique * get_metric() const
Returns a pointer on the metric defining the conformal factor for tensor densities.
Definition: tenseur.h:748
Cmp operator+(const Cmp &)
Definition: cmp_arithm.C:107
const Base_vect * get_triad() const
Returns the vectorial basis (triad) on which the components are defined.
Definition: tenseur.h:707
int get_etat() const
Returns the logical state.
Definition: tenseur.h:710
Cmp operator-(const Cmp &)
- Cmp
Definition: cmp_arithm.C:111
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: tenseur.C:652
void set_etat_zero()
Sets the logical state to ETATZERO (zero state).
Definition: tenseur.C:661
Tensor handling *** DEPRECATED : use class Tensor instead ***.
Definition: tenseur.h:304