LORENE
tensor_arithm.C
1/*
2 * Arithmetics functions for the Tensor class.
3 *
4 * These functions are not member functions of the Tensor class.
5 *
6 * (see file tensor.h for documentation).
7 *
8 */
9
10/*
11 * Copyright (c) 2003 Eric Gourgoulhon & Jerome Novak
12 * Copyright (c) 1999-2001 Philippe Grandclement (Tenseur version)
13 * Copyright (c) 2000-2001 Eric Gourgoulhon (Tenseur version)
14 *
15 * This file is part of LORENE.
16 *
17 * LORENE is free software; you can redistribute it and/or modify
18 * it under the terms of the GNU General Public License as published by
19 * the Free Software Foundation; either version 2 of the License, or
20 * (at your option) any later version.
21 *
22 * LORENE is distributed in the hope that it will be useful,
23 * but WITHOUT ANY WARRANTY; without even the implied warranty of
24 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
25 * GNU General Public License for more details.
26 *
27 * You should have received a copy of the GNU General Public License
28 * along with LORENE; if not, write to the Free Software
29 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
30 *
31 */
32
33char tensor_arithm_C[] = "$Header: /cvsroot/Lorene/C++/Source/Tensor/tensor_arithm.C,v 1.5 2014/10/13 08:53:44 j_novak Exp $" ;
34
35/*
36 * $Id: tensor_arithm.C,v 1.5 2014/10/13 08:53:44 j_novak Exp $
37 * $Log: tensor_arithm.C,v $
38 * Revision 1.5 2014/10/13 08:53:44 j_novak
39 * Lorene classes and functions now belong to the namespace Lorene.
40 *
41 * Revision 1.4 2014/10/06 15:13:20 j_novak
42 * Modified #include directives to use c++ syntax.
43 *
44 * Revision 1.3 2004/01/08 09:24:11 e_gourgoulhon
45 * Added arithmetics common to Scalar and Tensor.
46 * Corrected treatment ETATUN in Tensor / Scalar.
47 *
48 * Revision 1.2 2003/10/01 13:04:44 e_gourgoulhon
49 * The method Tensor::get_mp() returns now a reference (and not
50 * a pointer) onto a mapping.
51 *
52 * Revision 1.1 2003/09/26 14:33:53 j_novak
53 * Arithmetic functions for the class Tensor
54 *
55 *
56 *
57 * $Header: /cvsroot/Lorene/C++/Source/Tensor/tensor_arithm.C,v 1.5 2014/10/13 08:53:44 j_novak Exp $
58 *
59 */
60
61// Headers C
62#include <cstdlib>
63#include <cassert>
64#include <cmath>
65
66// Headers Lorene
67#include "tensor.h"
68
69 //********************//
70 // OPERATEURS UNAIRES //
71 //********************//
72
73namespace Lorene {
75
76 return t ;
77
78}
79
81
83 t.get_triad()) ;
84
85 for (int i=0 ; i<res.get_n_comp() ; i++) {
86 Itbl ind (res.indices(i)) ;
87 res.set(ind) = -t(ind) ;
88 }
89 return res ;
90
91}
92
93 //**********//
94 // ADDITION //
95 //**********//
96
97Tensor operator+(const Tensor & t1, const Tensor & t2) {
98
99 assert (t1.get_valence() == t2.get_valence()) ;
100 assert (t1.get_mp() == t2.get_mp()) ;
101 if (t1.get_valence() != 0) {
102 assert ( *(t1.get_triad()) == *(t2.get_triad()) ) ;
103 }
104
105 for (int i=0 ; i<t1.get_valence() ; i++)
106 assert(t1.get_index_type(i) == t2.get_index_type(i)) ;
107
108 Tensor res(t1.get_mp(), t1.get_valence(), t1.get_index_type(),
109 t1.get_triad()) ;
110
111 for (int i=0 ; i<res.get_n_comp() ; i++) {
112 Itbl ind (res.indices(i)) ;
113 res.set(ind) = t1(ind) + t2(ind) ;
114 }
115 return res ;
116
117}
118
119
121
122 assert (t1.get_valence() == 0) ;
123 assert (t1.get_mp() == t2.get_mp()) ;
124
125 return *(t1.cmp[0]) + t2 ;
126
127}
128
130
131 assert (t2.get_valence() == 0) ;
132 assert (t1.get_mp() == t2.get_mp()) ;
133
134 return t1 + *(t2.cmp[0]) ;
135
136}
137
138 //**************//
139 // SOUSTRACTION //
140 //**************//
141
142Tensor operator-(const Tensor & t1, const Tensor & t2) {
143
144 assert (t1.get_valence() == t2.get_valence()) ;
145 assert (t1.get_mp() == t2.get_mp()) ;
146 if (t1.get_valence() != 0) {
147 assert ( *(t1.get_triad()) == *(t2.get_triad()) ) ;
148 }
149
150 for (int i=0 ; i<t1.get_valence() ; i++)
151 assert(t1.get_index_type(i) == t2.get_index_type(i)) ;
152
153 Tensor res(t1.get_mp(), t1.get_valence(), t1.get_index_type(),
154 t1.get_triad()) ;
155
156 for (int i=0 ; i<res.get_n_comp() ; i++) {
157 Itbl ind (res.indices(i)) ;
158 res.set(ind) = t1(ind) - t2(ind) ;
159 }
160 return res ;
161
162}
163
165
166 assert (t1.get_valence() == 0) ;
167 assert (t1.get_mp() == t2.get_mp()) ;
168
169 return *(t1.cmp[0]) - t2 ;
170
171}
172
174
175 assert (t2.get_valence() == 0) ;
176 assert (t1.get_mp() == t2.get_mp()) ;
177
178 return t1 - *(t2.cmp[0]) ;
179
180}
181
182
183
184 //****************//
185 // MULTIPLICATION //
186 //****************//
187
189
190 assert (&(t1.get_mp()) == &(t2.get_mp())) ;
191
192 if (t1.get_etat() == ETATUN) return t2 ;
193
194 Tensor res(t2.get_mp(), t2.get_valence(), t2.get_index_type(),
195 t2.get_triad()) ;
196
197 for (int ic=0 ; ic<res.get_n_comp() ; ic++) {
198 Itbl ind = res.indices(ic) ;
199 res.set(ind) = t1 * t2(ind) ;
200 }
201
202 return res ;
203}
204
205
207
208 return t1*t2 ;
209}
210
211
212
213Tensor operator*(double x, const Tensor& t) {
214
216 t.get_triad()) ;
217
218 for (int i=0 ; i<res.get_n_comp() ; i++) {
219 Itbl ind (res.indices(i)) ;
220 res.set(ind) = x*t(ind) ;
221 }
222
223 return res ;
224
225}
226
227
228Tensor operator* (const Tensor& t, double x) {
229 return x * t ;
230}
231
232Tensor operator*(int m, const Tensor& t) {
233 return double(m) * t ;
234}
235
236
237Tensor operator* (const Tensor& t, int m) {
238 return double(m) * t ;
239}
240
241
242 //**********//
243 // DIVISION //
244 //**********//
245
247
248 // Protections
249 assert(s2.get_etat() != ETATNONDEF) ;
250 assert(t1.get_mp() == s2.get_mp()) ;
251
252 // Cas particuliers
253 if (s2.get_etat() == ETATZERO) {
254 cout << "Division by 0 in Tensor / Scalar !" << endl ;
255 abort() ;
256 }
257
258 if (s2.get_etat() == ETATUN) return t1 ;
259
260 Tensor res(t1.get_mp(), t1.get_valence(), t1.get_index_type(),
261 t1.get_triad()) ;
262
263 for (int i=0 ; i<res.get_n_comp() ; i++) {
264 Itbl ind (res.indices(i)) ;
265 res.set(ind) = t1(ind) / s2 ; // Scalar / Scalar
266 }
267 return res ;
268
269}
270
271
272Tensor operator/ (const Tensor& t, double x) {
273
274 if ( x == double(0) ) {
275 cout << "Division by 0 in Tensor / double !" << endl ;
276 abort() ;
277 }
278
279 if (x == double(1))
280 return t ;
281 else {
283 t.get_triad()) ;
284
285 for (int i=0 ; i<res.get_n_comp() ; i++) {
286 Itbl ind (res.indices(i)) ;
287 res.set(ind) = t(ind) / x ; // Scalar / double
288 }
289 return res ;
290 }
291
292}
293
294Tensor operator/ (const Tensor& t, int m) {
295
296 return t / double(m) ;
297}
298
299
300
301
302
303
304
305}
Time evolution with partial storage (*** under development ***).
Definition evolution.h:371
Basic integer array class.
Definition itbl.h:122
Tensor field of valence 0 (or component of a tensorial field).
Definition scalar.h:387
Tensor handling.
Definition tensor.h:288
Base_val operator*(const Base_val &, const Base_val &)
This operator is used when calling multiplication or division of Valeur .
Cmp operator-(const Cmp &)
- Cmp
Definition cmp_arithm.C:108
Cmp operator/(const Cmp &, const Cmp &)
Cmp / Cmp.
Definition cmp_arithm.C:457
Cmp operator+(const Cmp &)
Definition cmp_arithm.C:104
const Map & get_mp() const
Returns the mapping.
Definition tensor.h:861
int get_index_type(int i) const
Gives the type (covariant or contravariant) of the index number i .
Definition tensor.h:886
int get_valence() const
Returns the valence.
Definition tensor.h:869
const Base_vect * get_triad() const
Returns the vectorial basis (triad) on which the components are defined.
Definition tensor.h:866
Lorene prototypes.
Definition app_hor.h:64