LORENE
tensor_change_triad.C
1 /*
2  * Methods for changing the triad of a Tensor
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 
27 
28 
29 /*
30  * $Id: tensor_change_triad.C,v 1.10 2016/12/05 16:18:17 j_novak Exp $
31  * $Log: tensor_change_triad.C,v $
32  * Revision 1.10 2016/12/05 16:18:17 j_novak
33  * Suppression of some global variables (file names, loch, ...) to prevent redefinitions
34  *
35  * Revision 1.9 2014/10/13 08:53:44 j_novak
36  * Lorene classes and functions now belong to the namespace Lorene.
37  *
38  * Revision 1.8 2014/10/06 15:13:20 j_novak
39  * Modified #include directives to use c++ syntax.
40  *
41  * Revision 1.7 2005/09/15 15:51:26 j_novak
42  * The "rotation" (change of triad) methods take now Scalars as default
43  * arguments.
44  *
45  * Revision 1.6 2005/02/03 14:31:37 f_limousin
46  * Correction of an error in the case Cartesian --> Cartesian for
47  * a Sym_tensor. Now the components of the tensor are modified
48  * using a temporary.
49  *
50  * Revision 1.5 2003/10/28 21:29:08 e_gourgoulhon
51  * -- Read-only access to the components performed via operator()(int, int)
52  * instead of set(int, int).
53  * -- Corrected index range in the case Cartesian -> Cartesian.
54  *
55  * Revision 1.4 2003/10/27 10:50:24 e_gourgoulhon
56  * Added the case of a twice contravariant tensor in the assert.
57  *
58  * Revision 1.3 2003/10/06 14:25:51 j_novak
59  * Added a test #ifndef... to prevent a warning
60  *
61  * Revision 1.2 2003/10/05 21:12:19 e_gourgoulhon
62  * - Modified some assert.
63  * - Corrected bug on index range in line 200.
64  *
65  * Revision 1.1 2003/09/29 12:52:57 j_novak
66  * Methods for changing the triad are implemented.
67  *
68  *
69  * $Header: /cvsroot/Lorene/C++/Source/Tensor/tensor_change_triad.C,v 1.10 2016/12/05 16:18:17 j_novak Exp $
70  *
71  */
72 
73 // C headers
74 #include <cassert>
75 
76 // Lorene headers
77 #include "tensor.h"
78 
79 namespace Lorene {
80 void Tensor::change_triad(const Base_vect& new_triad) {
81 
82  assert (valence == 2) ;
83  assert(triad != 0x0) ;
84 
85  const Base_vect_cart* nbvc = dynamic_cast<const Base_vect_cart*>(&new_triad) ;
86 #ifndef NDEBUG
87  const Base_vect_spher* nbvs
88  = dynamic_cast<const Base_vect_spher*>(&new_triad) ;
89 #endif
90 
91  assert((nbvc != 0x0) || (nbvs != 0x0)) ;
92 
93  const Base_vect_cart* bvc = dynamic_cast<const Base_vect_cart*>(triad) ;
94  const Base_vect_spher* bvs = dynamic_cast<const Base_vect_spher*>(triad) ;
95 
96  assert((bvc != 0x0) || (bvs != 0x0)) ;
97 
98  // ---------------------------------------------
99  // Case where the input triad is a Cartesian one
100  // ---------------------------------------------
101  if (nbvc != 0x0) {
102  assert(nbvs == 0x0) ;
103 
104  // -----------------------------
105  // Case cartesian -> cartesian
106  // -----------------------------
107  if (bvc != 0x0) { // The old triad is a cartesian one
108  assert(bvs == 0x0) ;
109 
110  int ind = nbvc->get_align() * (bvc->get_align()) ;
111 
112  switch (ind) {
113 
114  case 1 : { // the two bases are aligned : nothing to do
115  // -----------------------------------------
116 
117  break ;
118  }
119 
120  case - 1 : { // the two bases are anti-aligned
121  // ------------------------------
122 
123  Tensor copie (*this) ;
124 
125  set(1, 3) = - copie(1, 3) ; // {xz} --> - {xz}
126  set(2, 3) = - copie(2, 3) ; // {yz} --> - {yz}
127  set(3, 1) = - copie(3, 1) ; // {zx} --> - {zx}
128  set(3, 2) = - copie(3, 2) ; // {zy} --> - {zy}
129  // all other components are unchanged
130  break ;
131  }
132  case 0 : { // the two basis have not a special relative orientation
133  // -----------------------------------------------------
134  cout <<
135  "Tensor::change_basis : general value of rot_phi "
136  << " not contemplated yet, sorry !" << endl ;
137  abort() ;
138  break ;
139  }
140 
141  default : { // error
142  cout <<
143  "Tensor::change_basis : unexpected value of ind !" << endl ;
144  cout << " ind = " << ind << endl ;
145  abort() ;
146  break ;
147  }
148  }
149 
150  } // end of the cart -> cart basis case
151 
152 
153  // -----------------------------
154  // Case spherical -> cartesian
155  // -----------------------------
156  if (bvs != 0x0) { // The old triad is a spherical one
157 
158  assert(bvc == 0x0) ;
159 
160  // The triads should be the same as that associated
161  // with the mapping :
162  assert( *nbvc == mp->get_bvect_cart() ) ;
163  assert( *bvs == mp->get_bvect_spher() ) ;
164 
165  // Only for double-covariant tensors or double-contravariant tensors
166  assert( ( (type_indice(0)==COV) && (type_indice(1)==COV) ) ||
167  ( (type_indice(0)==CON) && (type_indice(1)==CON) ) ) ;
168 #ifndef NDEBUG
169  int nz = mp->get_mg()->get_nzone() ;
170  for (int i=0; i<nz; i++) {
171  assert( mp->get_mg()->get_np(i) >= 4) ;
172  assert( mp->get_mg()->get_nt(i) >= 5) ;
173  }
174 #endif
175  // Temporary storage of the components
176  // the Base_vect *this is not valid...
177  Tensor tmp(*mp, 2, COV, *triad) ;
178  for (int i=1; i<=3; i++) {
179  mp->comp_x_from_spherical(operator()(1,i), operator()(2,i),
180  operator()(3,i), tmp.set(1,i)) ;
181  mp->comp_y_from_spherical(operator()(1,i), operator()(2,i),
182  operator()(3,i), tmp.set(2,i)) ;
183  mp->comp_z_from_spherical(operator()(1,i), operator()(2,i), tmp.set(3,i) ) ;
184  }
185  for (int i=1; i<=3; i++) {
186  mp->comp_x_from_spherical(tmp(i,1), tmp(i,2), tmp(i,3), set(i,1)) ;
187  mp->comp_y_from_spherical(tmp(i,1), tmp(i,2), tmp(i,3), set(i,2)) ;
188  mp->comp_z_from_spherical(tmp(i,1), tmp(i,2), set(i,3)) ;
189  }
190 
191  }// End of the spher -> cart case
192  } // End of the case of cartesian new triad
193 
194  // ---------------------------------------------
195  // Case where the new triad is a spherical one
196  // ---------------------------------------------
197  else {
198 
199  assert(nbvc == 0x0) ;
200 
201  // ---------------------------------
202  // Case cartesian -> spherical
203  // ---------------------------------
204  if (bvc != 0x0) { // The old triad is a cartesian one
205  assert(bvs == 0x0) ;
206 
207  // The triads should be the same as that associated
208  // with the mapping :
209  assert( *nbvs == mp->get_bvect_spher() ) ;
210  assert( *bvc == mp->get_bvect_cart() ) ;
211 
212  // Only for double-covariant tensors or double-contravariant tensors
213  assert( ( (type_indice(0)==COV) && (type_indice(1)==COV) ) ||
214  ( (type_indice(0)==CON) && (type_indice(1)==CON) ) ) ;
215 #ifndef NDEBUG
216  int nz = mp->get_mg()->get_nzone() ;
217  for (int i=0; i<nz; i++) {
218  assert( mp->get_mg()->get_np(i) >= 4) ;
219  assert( mp->get_mg()->get_nt(i) >= 5) ;
220  }
221 #endif
222 
223  // Temporary storage of the components
224  Tensor tmp(*mp, 2, COV, *triad) ;
225  for (int i=1; i<=3; i++) {
226  mp->comp_r_from_cartesian(operator()(1,i), operator()(2,i),
227  operator()(3,i), tmp.set(1,i)) ;
228  mp->comp_t_from_cartesian(operator()(1,i), operator()(2,i),
229  operator()(3,i), tmp.set(2,i)) ;
230  mp->comp_p_from_cartesian(operator()(1,i), operator()(2,i), tmp.set(3,i)) ;
231  }
232  for (int i=1; i<=3; i++) {
233  mp->comp_r_from_cartesian(tmp(i,1), tmp(i,2), tmp(i,3), set(i,1)) ;
234  mp->comp_t_from_cartesian(tmp(i,1), tmp(i,2), tmp(i,3), set(i,2)) ;
235  mp->comp_p_from_cartesian(tmp(i,1), tmp(i,2), set(i,3)) ;
236  }
237  } // end of the case cart -> spher
238 
239 
240  // ------------------------------------
241  // Case spherical -> spherical
242  // ------------------------------------
243  if (bvs != 0x0) {
244 
245  assert(bvc == 0x0) ;
246 
247  cout << "Tensor::change_triad : case not treated yet !" << endl ;
248  abort() ;
249  } // end of the spher->spher basis case
250 
251  } // End of the case of spherical new triad
252 
253  triad = &new_triad ;
254 
255 }
256 }
int get_align() const
Returns the indicator of alignment with respect to the absolute frame.
Definition: base_vect.h:285
int get_np(int l) const
Returns the number of points in the azimuthal direction ( ) in domain no. l.
Definition: grilles.h:479
const Base_vect_spher & get_bvect_spher() const
Returns the orthonormal vectorial basis associated with the coordinates of the mapping.
Definition: map.h:801
Lorene prototypes.
Definition: app_hor.h:67
virtual void comp_p_from_cartesian(const Scalar &v_x, const Scalar &v_y, Scalar &v_p) const =0
Computes the Spherical component (with respect to bvect_spher ) of a vector given by its cartesian c...
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:783
const Base_vect * triad
Vectorial basis (triad) with respect to which the tensor components are defined.
Definition: tensor.h:309
Vectorial bases (triads) with respect to which the tensorial components are defined.
Definition: base_vect.h:105
virtual void comp_x_from_spherical(const Scalar &v_r, const Scalar &v_theta, const Scalar &v_phi, Scalar &v_x) const =0
Computes the Cartesian x component (with respect to bvect_cart ) of a vector given by its spherical c...
virtual void comp_t_from_cartesian(const Scalar &v_x, const Scalar &v_y, const Scalar &v_z, Scalar &v_t) const =0
Computes the Spherical component (with respect to bvect_spher ) of a vector given by its cartesian c...
virtual void comp_z_from_spherical(const Scalar &v_r, const Scalar &v_theta, Scalar &v_z) const =0
Computes the Cartesian z component (with respect to bvect_cart ) of a vector given by its spherical c...
Itbl type_indice
1D array of integers (class Itbl ) of size valence containing the type of each index: COV for a cova...
Definition: tensor.h:316
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:465
Tensor handling.
Definition: tensor.h:294
Cartesian vectorial bases (triads).
Definition: base_vect.h:201
Spherical orthonormal vectorial bases (triads).
Definition: base_vect.h:308
int valence
Valence of the tensor (0 = scalar, 1 = vector, etc...)
Definition: tensor.h:304
const Base_vect_cart & get_bvect_cart() const
Returns the Cartesian basis associated with the coordinates (x,y,z) of the mapping, i.e.
Definition: map.h:809
Scalar & set(const Itbl &ind)
Returns the value of a component (read/write version).
Definition: tensor.C:663
virtual void comp_r_from_cartesian(const Scalar &v_x, const Scalar &v_y, const Scalar &v_z, Scalar &v_r) const =0
Computes the Spherical r component (with respect to bvect_spher ) of a vector given by its cartesian ...
virtual void change_triad(const Base_vect &new_triad)
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.
int get_nt(int l) const
Returns the number of points in the co-latitude direction ( ) in domain no. l.
Definition: grilles.h:474
const Map *const mp
Mapping on which the numerical values at the grid points are defined.
Definition: tensor.h:301
virtual void comp_y_from_spherical(const Scalar &v_r, const Scalar &v_theta, const Scalar &v_phi, Scalar &v_y) const =0
Computes the Cartesian y component (with respect to bvect_cart ) of a vector given by its spherical c...