LORENE
connection_fcart.C
1/*
2 * Methods of class Connection_fcart.
3 *
4 * (see file connection.h for documentation)
5 *
6 */
7
8/*
9 * Copyright (c) 2003-2004 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 connection_fcart_C[] = "$Header: /cvsroot/Lorene/C++/Source/Connection/connection_fcart.C,v 1.14 2014/10/13 08:52:50 j_novak Exp $" ;
29
30/*
31 * $Id: connection_fcart.C,v 1.14 2014/10/13 08:52:50 j_novak Exp $
32 * $Log: connection_fcart.C,v $
33 * Revision 1.14 2014/10/13 08:52:50 j_novak
34 * Lorene classes and functions now belong to the namespace Lorene.
35 *
36 * Revision 1.13 2014/10/06 15:13:04 j_novak
37 * Modified #include directives to use c++ syntax.
38 *
39 * Revision 1.12 2004/01/28 13:25:40 j_novak
40 * The ced_mult_r arguments have been suppressed from the Scalar::*dsd* methods.
41 * In the div/mult _r_dzpuis, there is no more default value.
42 *
43 * Revision 1.11 2004/01/04 21:00:50 e_gourgoulhon
44 * Better handling of tensor symmetries in methods p_derive_cov() and
45 * p_divergence() (thanks to the new class Tensor_sym).
46 *
47 * Revision 1.10 2004/01/01 11:24:04 e_gourgoulhon
48 * Full reorganization of method p_derive_cov: the main loop is now
49 * on the indices of the *output* tensor (to take into account
50 * symmetries in the input and output tensors).
51 *
52 * Revision 1.9 2003/12/27 14:59:52 e_gourgoulhon
53 * -- Method derive_cov() suppressed.
54 * -- Change of the position of the derivation index from the first one
55 * to the last one in methods p_derive_cov() and p_divergence().
56 *
57 * Revision 1.8 2003/10/17 13:46:15 j_novak
58 * The argument is now between 1 and 3 (instead of 0->2)
59 *
60 * Revision 1.7 2003/10/16 21:37:08 e_gourgoulhon
61 * Corrected deriv index in divergence.
62 *
63 * Revision 1.6 2003/10/16 15:26:03 e_gourgoulhon
64 * Suppressed unsued variable
65 *
66 * Revision 1.5 2003/10/16 14:21:36 j_novak
67 * The calculation of the divergence of a Tensor is now possible.
68 *
69 * Revision 1.4 2003/10/11 16:45:43 e_gourgoulhon
70 * Suppressed the call to Itbl::set_etat_qcq() after
71 * the construction of the Itbl's.
72 *
73 * Revision 1.3 2003/10/11 14:39:50 e_gourgoulhon
74 * Suppressed declaration of unusued arguments in some methods.
75 *
76 * Revision 1.2 2003/10/06 13:58:46 j_novak
77 * The memory management has been improved.
78 * Implementation of the covariant derivative with respect to the exact Tensor
79 * type.
80 *
81 * Revision 1.1 2003/10/03 14:11:48 e_gourgoulhon
82 * Methods of class Connection_fcart.
83 *
84 *
85 *
86 * $Header: /cvsroot/Lorene/C++/Source/Connection/connection_fcart.C,v 1.14 2014/10/13 08:52:50 j_novak Exp $
87 *
88 */
89
90// C++ headers
91#include "headcpp.h"
92
93// C headers
94#include <cstdlib>
95
96// Lorene headers
97#include "connection.h"
98
99
100 //------------------------------//
101 // Constructors //
102 //------------------------------//
103
104
105
106// Contructor from a Cartesian flat-metric-orthonormal basis
107
108namespace Lorene {
110 : Connection_flat(mpi, bi) {
111
112}
113
114// Copy constructor
119
120
121 //----------------------------//
122 // Destructor //
123 //----------------------------//
124
125
129
130
131 //--------------------------------//
132 // Mutators / assignment //
133 //--------------------------------//
134
135
137
138 cout << "Connection_fcart::operator= : not implemented yet !" << endl ;
139 abort() ;
140
141}
142
143
144
145 //-----------------------------//
146 // Computational methods //
147 //-----------------------------//
148
149// Covariant derivative, returning a pointer.
150//-------------------------------------------
151
153
154 // Notations: suffix 0 in name <=> input tensor
155 // suffix 1 in name <=> output tensor
156
157 int valence0 = uu.get_valence() ;
158 int valence1 = valence0 + 1 ;
159 int valence1m1 = valence1 - 1 ; // same as valence0, but introduced for
160 // the sake of clarity
161
162 // Protections
163 // -----------
164 if (valence0 >= 1) {
165 assert(uu.get_triad() == triad) ;
166 }
167
168 // Creation of the result (pointer)
169 // --------------------------------
170 Tensor* resu ;
171
172 // If uu is a Scalar, the result is a vector
173 if (valence0 == 0)
174 resu = new Vector(*mp, COV, triad) ;
175 else {
176
177 // Type of indices of the result :
178 Itbl tipe(valence1) ;
179 const Itbl& tipeuu = uu.get_index_type() ;
180 for (int id = 0; id<valence0; id++) {
181 tipe.set(id) = tipeuu(id) ; // First indices = same as uu
182 }
183 tipe.set(valence1m1) = COV ; // last index is the derivation index
184
185 // if uu is a Tensor_sym, the result is also a Tensor_sym:
186 const Tensor* puu = &uu ;
187 const Tensor_sym* puus = dynamic_cast<const Tensor_sym*>(puu) ;
188 if ( puus != 0x0 ) { // the input tensor is symmetric
189 resu = new Tensor_sym(*mp, valence1, tipe, *triad,
190 puus->sym_index1(), puus->sym_index2()) ;
191 }
192 else {
193 resu = new Tensor(*mp, valence1, tipe, *triad) ; // no symmetry
194 }
195
196 }
197
198 int ncomp1 = resu->get_n_comp() ;
199
200 Itbl ind1(valence1) ; // working Itbl to store the indices of resu
201 Itbl ind0(valence0) ; // working Itbl to store the indices of uu
202
203 // Loop on all the components of the output tensor
204 // -----------------------------------------------
205 for (int ic=0; ic<ncomp1; ic++) {
206
207 // indices corresponding to the component no. ic in the output tensor
208 ind1 = resu->indices(ic) ;
209
210 // Component no. ic:
211 Scalar& cresu = resu->set(ind1) ;
212
213 // Indices of the input tensor
214 for (int id = 0; id < valence0; id++) {
215 ind0.set(id) = ind1(id) ;
216 }
217
218 // Value of last index (derivation index)
219 int k = ind1(valence1m1) ;
220
221 // Partial derivation with respect to x^k:
222
223 cresu = (uu(ind0)).deriv(k) ;
224
225 }
226
227 // C'est fini !
228 // -----------
229 return resu ;
230
231}
232
233
234
235// Divergence, returning a pointer.
236//---------------------------------
237
239
240 // Notations: suffix 0 in name <=> input tensor
241 // suffix 1 in name <=> output tensor
242
243 int valence0 = uu.get_valence() ;
244 int valence1 = valence0 - 1 ;
245
246 // Protections
247 // -----------
248 assert (valence0 >= 1) ;
249 assert (uu.get_triad() == triad) ;
250
251 // Last index must be contravariant:
252 assert (uu.get_index_type(valence0-1) == CON) ;
253
254
255 // Creation of the pointer on the result tensor
256 // --------------------------------------------
257 Tensor* resu ;
258
259 if (valence0 == 1) // if u is a Vector, the result is a Scalar
260 resu = new Scalar(*mp) ;
261 else {
262
263 // Type of indices of the result :
264 Itbl tipe(valence1) ;
265 const Itbl& tipeuu = uu.get_index_type() ;
266 for (int id = 0; id<valence1; id++) {
267 tipe.set(id) = tipeuu(id) ; // type of remaining indices =
268 } // same as uu indices
269
270 if (valence0 == 2) { // if u is a rank 2 tensor, the result is a Vector
271 resu = new Vector(*mp, tipe(0), *triad) ;
272 }
273 else {
274 // if uu is a Tensor_sym, the result might be also a Tensor_sym:
275 const Tensor* puu = &uu ;
276 const Tensor_sym* puus = dynamic_cast<const Tensor_sym*>(puu) ;
277 if ( puus != 0x0 ) { // the input tensor is symmetric
278
279 if (puus->sym_index2() != valence0 - 1) {
280
281 // the symmetry is preserved:
282
283 if (valence1 == 2) {
284 resu = new Sym_tensor(*mp, tipe, *triad) ;
285 }
286 else {
287 resu = new Tensor_sym(*mp, valence1, tipe, *triad,
288 puus->sym_index1(), puus->sym_index2()) ;
289 }
290 }
291 else { // the symmetry is lost:
292
293 resu = new Tensor(*mp, valence1, tipe, *triad) ;
294 }
295 }
296 else { // no symmetry in the input tensor:
297 resu = new Tensor(*mp, valence1, tipe, *triad) ;
298 }
299 }
300 }
301
302
303 int ncomp1 = resu->get_n_comp() ;
304
305 Itbl ind0(valence0) ; // working Itbl to store the indices of uu
306
307 Itbl ind1(valence1) ; // working Itbl to store the indices of resu
308
309 // Loop on all the components of the output tensor
310 for (int ic=0; ic<ncomp1; ic++) {
311
312 ind1 = resu->indices(ic) ;
313 Scalar& cresu = resu->set(ind1) ;
314 cresu.set_etat_zero() ;
315
316 for (int k=1; k<=3; k++) {
317
318 // indices (ind1,k) in the input tensor
319 for (int id = 0; id < valence1; id++) {
320 ind0.set(id) = ind1(id) ;
321 }
322 ind0.set(valence0-1) = k ;
323
324 cresu += uu(ind0).deriv(k) ; //Addition of dT^i/dx^i
325 }
326
327 }
328
329 // C'est fini !
330 // -----------
331 return resu ;
332
333}
334
335}
Cartesian vectorial bases (triads).
Definition base_vect.h:201
Class Connection_fcart.
Definition connection.h:546
void operator=(const Connection_fcart &)
Assignment to another Connection_fcart.
Connection_fcart(const Map &, const Base_vect_cart &)
Contructor from a Cartesian flat-metric-orthonormal basis.
virtual ~Connection_fcart()
destructor
virtual Tensor * p_divergence(const Tensor &tens) const
Computes the divergence of a tensor (with respect to the current connection).
virtual Tensor * p_derive_cov(const Tensor &tens) const
Computes the covariant derivative of a tensor (with respect to the current connection).
Class Connection_flat.
Definition connection.h:354
const Base_vect *const triad
Triad with respect to which the connection coefficients are defined.
Definition connection.h:124
const Map *const mp
Reference mapping.
Definition connection.h:119
Basic integer array class.
Definition itbl.h:122
int & set(int i)
Read/write of a particular element (index i ) (1D case)
Definition itbl.h:247
Base class for coordinate mappings.
Definition map.h:670
Tensor field of valence 0 (or component of a tensorial field).
Definition scalar.h:387
virtual void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition scalar.C:324
Class intended to describe valence-2 symmetric tensors.
Definition sym_tensor.h:223
Symmetric tensors (with respect to two of their arguments).
Definition tensor.h:1037
Tensor handling.
Definition tensor.h:288
Tensor field of valence 1.
Definition vector.h:188
int sym_index1() const
Number of the first symmetric index (0<= id_sym1 < valence )
Definition tensor.h:1149
int get_index_type(int i) const
Gives the type (covariant or contravariant) of the index number i .
Definition tensor.h:886
int sym_index2() const
Number of the second symmetric index (id_sym1 < id_sym2 < valence )
Definition tensor.h:1154
int get_valence() const
Returns the valence.
Definition tensor.h:869
virtual Itbl indices(int pos) const
Returns the indices of a component given by its position in the array cmp .
Definition tensor.C:539
int get_n_comp() const
Returns the number of stored components.
Definition tensor.h:872
const Base_vect * get_triad() const
Returns the vectorial basis (triad) on which the components are defined.
Definition tensor.h:866
Scalar & set(const Itbl &ind)
Returns the value of a component (read/write version).
Definition tensor.C:654
Lorene prototypes.
Definition app_hor.h:64