LORENE
valeur_ylm_i.C
1 /*
2  * Copyright (c) 1999-2003 Eric Gourgoulhon
3  * Copyright (c) 1999-2001 Philippe Grandclement
4  *
5  * This file is part of LORENE.
6  *
7  * LORENE is free software; you can redistribute it and/or modify
8  * it under the terms of the GNU General Public License as published by
9  * the Free Software Foundation; either version 2 of the License, or
10  * (at your option) any later version.
11  *
12  * LORENE is distributed in the hope that it will be useful,
13  * but WITHOUT ANY WARRANTY; without even the implied warranty of
14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15  * GNU General Public License for more details.
16  *
17  * You should have received a copy of the GNU General Public License
18  * along with LORENE; if not, write to the Free Software
19  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
20  *
21  */
22 
23 
24 
25 
26 /*
27  * Fonction membre de la classe Valeur qui calcule les coefficients
28  * de la decomposition en cos(l*theta) / sin(l*theta)
29  * a partir des coefficients de la decomp. en harmoniques spheriques
30  *
31  */
32 
33 
34 /*
35  * $Id: valeur_ylm_i.C,v 1.13 2016/12/05 16:18:21 j_novak Exp $
36  * $Log: valeur_ylm_i.C,v $
37  * Revision 1.13 2016/12/05 16:18:21 j_novak
38  * Suppression of some global variables (file names, loch, ...) to prevent redefinitions
39  *
40  * Revision 1.12 2014/10/13 08:53:51 j_novak
41  * Lorene classes and functions now belong to the namespace Lorene.
42  *
43  * Revision 1.11 2014/10/06 15:13:25 j_novak
44  * Modified #include directives to use c++ syntax.
45  *
46  * Revision 1.10 2009/10/23 12:56:29 j_novak
47  * New base T_LEG_MI
48  *
49  * Revision 1.9 2009/10/13 13:49:58 j_novak
50  * New base T_LEG_MP.
51  *
52  * Revision 1.8 2005/11/22 10:13:44 p_grandclement
53  * Correction of the normalization in the case nt=1
54  *
55  * Revision 1.7 2005/02/16 15:33:14 m_forot
56  * Correct the case T_LEG
57  *
58  * Revision 1.6 2004/11/23 15:17:20 m_forot
59  * Added the bases for the cases without any equatorial symmetry
60  * (T_COSSIN_C, T_COSSIN_S, T_LEG, R_CHEBPI_P, R_CHEBPI_I).
61  *
62  * Revision 1.5 2003/10/13 20:52:58 e_gourgoulhon
63  * Loop variables i and l have now local scope.
64  *
65  * Revision 1.4 2003/09/17 12:30:22 j_novak
66  * New checks for changing to T_LEG* bases.
67  *
68  * Revision 1.3 2003/09/16 08:54:09 j_novak
69  * Addition of the T_LEG_II base (odd in theta, only for odd m) and the
70  * transformation functions to and from the T_SIN_P base.
71  *
72  * Revision 1.2 2002/10/16 14:37:16 j_novak
73  * Reorganization of #include instructions of standard C++, in order to
74  * use experimental version 3 of gcc.
75  *
76  * Revision 1.1.1.1 2001/11/20 15:19:27 e_gourgoulhon
77  * LORENE
78  *
79  * Revision 2.11 2000/09/29 16:10:37 eric
80  * Ajout des bases T_LEG_IP et T_LEG_PI.
81  *
82  * Revision 2.10 2000/03/31 15:58:26 phil
83  * changement des bases meme si etat est zero
84  *
85  * Revision 2.9 1999/12/22 16:25:52 eric
86  * Traitement du cas ETATZERO
87  * Test sur c_cf avant d'appeler coef().
88  *
89  * Revision 2.8 1999/12/16 16:41:43 phil
90  * *** empty log message ***
91  *
92  * Revision 2.7 1999/12/16 16:10:07 phil
93  * correction cas nt = 1.
94  *
95  * Revision 2.6 1999/11/30 12:46:40 eric
96  * Valeur::base est desormais du type Base_val et non plus Base_val*.
97  *
98  * Revision 2.5 1999/11/22 11:36:14 eric
99  * Suppression de #include "tenseur.h"
100  *
101  * Revision 2.4 1999/10/18 14:12:39 eric
102  * Les bases sont desormais membres des Mtbl_cf.
103  *
104  * Revision 2.3 1999/04/14 10:19:52 phil
105  * *** empty log message ***
106  *
107  * Revision 2.2 1999/04/14 09:53:21 phil
108  * remplacement de malloc en new
109  *
110  * Revision 2.1 1999/04/14 09:38:47 phil
111  * Changement liberation memoire : free -> delete
112  *
113  * Revision 2.0 1999/04/13 16:47:28 phil
114  * *** empty log message ***
115  *
116  *
117  * $Header: /cvsroot/Lorene/C++/Source/Valeur/valeur_ylm_i.C,v 1.13 2016/12/05 16:18:21 j_novak Exp $
118  *
119  */
120 
121 // headers C
122 #include <cassert>
123 #include <cstdlib>
124 #include <cmath>
125 
126 // headers Lorene
127 #include "type_parite.h"
128 #include "valeur.h"
129 #include "proto.h"
130 
131 namespace Lorene {
132 void ylm_i_pasprevu(const int*, const double*, double*) ;
133 
135 
136  static void (*chbase_t[MAX_BASE])(const int*, const double*,
137  double*) ;
138  static int nouv_base_t[MAX_BASE] ;
139  static int premier_appel = 1 ;
140 
141  int deg[3] ;
142 
143  if (premier_appel==1) {
144  premier_appel = 0 ;
145 
146  for (int i=0; i<MAX_BASE; i++) {
147  chbase_t[i] = ylm_i_pasprevu ;
148  nouv_base_t[i] = NONDEF ;
149  }
150 
151  chbase_t[T_LEG_P >> TRA_T] = chb_legp_cossincp ;
152  nouv_base_t[T_LEG_P >> TRA_T] = T_COSSIN_CP ;
153 
154  chbase_t[T_LEG_I >> TRA_T] = chb_legi_cossinci ;
155  nouv_base_t[T_LEG_I >> TRA_T] = T_COSSIN_CI ;
156 
157  chbase_t[T_LEG >> TRA_T] = chb_leg_cossinc ;
158  nouv_base_t[T_LEG >> TRA_T] = T_COSSIN_C ;
159 
160  chbase_t[T_LEG_MP >> TRA_T] = chb_legmp_cos ;
161  nouv_base_t[T_LEG_MP >> TRA_T] = T_COS ;
162 
163  chbase_t[T_LEG_MI >> TRA_T] = chb_legmi_sin ;
164  nouv_base_t[T_LEG_MI >> TRA_T] = T_SIN ;
165 
166  chbase_t[T_LEG_PP >> TRA_T] = chb_legpp_cosp ;
167  nouv_base_t[T_LEG_PP >> TRA_T] = T_COS_P ;
168 
169  chbase_t[T_LEG_IP >> TRA_T] = chb_legip_cosi ;
170  nouv_base_t[T_LEG_IP >> TRA_T] = T_COS_I ;
171 
172  chbase_t[T_LEG_PI >> TRA_T] = chb_legpi_sini ;
173  nouv_base_t[T_LEG_PI >> TRA_T] = T_SIN_I ;
174 
175  chbase_t[T_LEG_II >> TRA_T] = chb_legii_sinp ;
176  nouv_base_t[T_LEG_II >> TRA_T] = T_SIN_P ;
177  }
178 
179 //---------------------------------------------------------------------------
180 // fin des operation de premier appel
181 //---------------------------------------------------------------------------
182 
183  // Tout null ?
184  int nzone = get_mg()->get_nzone() ;
185  if (etat == ETATZERO) {
186  for (int l=0 ; l<nzone ; l++) {
187  int vbase_r = base.b[l] & MSQ_R ;
188  int vbase_t = base.b[l] & MSQ_T ;
189  int vbase_p = base.b[l] & MSQ_P ;
190 
191  int vbase_t_tra = vbase_t >> TRA_T ;
192  base.b[l] = ( vbase_p | nouv_base_t[vbase_t_tra] ) | vbase_r ;
193  }
194  return ;
195  }
196 
197  // Protection
198  assert(etat != ETATNONDEF) ;
199 
200  if (c_cf == 0x0) {
201  coef() ; // The coefficients are required
202  }
203 
204 // Boucle sur les differentes zones
205 
206 
207  for (int l=0; l<nzone; l++) {
208 
209 // On recupere les anciennes bases en r, phi et theta :
210  int vbase_r = base.b[l] & MSQ_R ;
211  int vbase_t = base.b[l] & MSQ_T ;
212  int vbase_p = base.b[l] & MSQ_P ;
213 
214  if ((vbase_t != T_COSSIN_CP) && (vbase_t != T_COS_I) &&
215  (vbase_t != T_COS_P) && (vbase_t != T_COSSIN_CI) &&
216  (vbase_t != T_SIN_P) && (vbase_t != T_SIN_I) &&
217  (vbase_t != T_COSSIN_C) && (vbase_t != T_COSSIN_S) &&
218  (vbase_t != T_COS) && (vbase_t != T_SIN) )
219  { // cas ou le calcul est necessaire
220 
221  int vbase_t_tra = vbase_t >> TRA_T ;
222  assert(vbase_t_tra < MAX_BASE) ;
223  bool pair = ( (vbase_t == T_LEG_PP) || (vbase_t == T_LEG_IP) ||
224  (vbase_t == T_LEG_MP) ) ;
225  bool impair = ( (vbase_t == T_LEG_PI) || (vbase_t == T_LEG_II) ||
226  (vbase_t == T_LEG_MI) ) ;
227 
228 // Nouvelle base :
229  base.b[l] = ( vbase_p | nouv_base_t[vbase_t_tra] ) | vbase_r ;
230  if (get_mg()->get_nt(l)==1) {
231  *c_cf->t[l] /= sqrt(2.) ;
232  continue ;
233  }
234 //... tbl contenant les coefficients dans la zone l :
235  Tbl* cf = c_cf->t[l] ;
236 
237  if (cf->get_etat() == ETATZERO) continue ; // On ne fait rien si le tbl = 0
238 
239 //... resultat du calcul :
240  double* resu = new double [cf->get_taille()] ;
241 
242  deg[0] = get_mg()->get_np(l) ; // nb. de degres de liberte en phi
243  deg[1] = get_mg()->get_nt(l) ; // nb. de degres de liberte en theta
244  deg[2] = get_mg()->get_nr(l) ; // nb. de degres de liberte en r
245 
246 
247 // Transformation en theta:
248 //-------------------------
249 
250  if ((pair && (vbase_p == P_COSSIN_I)) ||
251  (impair && (vbase_p == P_COSSIN_P)) )
252  ylm_i_pasprevu(deg, (cf->t), resu ) ;
253  else
254  chbase_t[vbase_t_tra](deg, (cf->t), resu ) ;
255 
256  // On branche le tbl contenant les coef sur resu :
257  delete [] cf->t ; // les anciens coef. sont oublies
258  cf->t = resu ; // nouveaux coef.
259 
260  } // fin du cas ou la transformation devait etre effectuee
261 
262  } // fin de la boucle sur les differentes zones
263 
264  // On met les bonnes bases dans c_cf :
265  c_cf->base = base ;
266 
267 }
268 
269 //******************************************************************************
270 
271 void ylm_i_pasprevu(const int*, const double*, double*) {
272 
273  cout <<
274  "Valeur::ylm_i: change of basis not implemented yet !"
275  << endl ;
276  abort() ;
277 }
278 
279 
280 }
#define T_LEG
fct. de Legendre associees
Definition: type_parite.h:236
Mtbl_cf * c_cf
Coefficients of the spectral expansion of the function.
Definition: valeur.h:312
#define T_LEG_MP
fct. de Legendre associees avec m pair
Definition: type_parite.h:238
void ylm_i()
Inverse of ylm()
Definition: valeur_ylm_i.C:134
int get_np(int l) const
Returns the number of points in the azimuthal direction ( ) in domain no. l.
Definition: grilles.h:479
void coef() const
Computes the coeffcients of *this.
Definition: valeur_coef.C:151
Cmp sqrt(const Cmp &)
Square root.
Definition: cmp_math.C:223
#define T_LEG_PI
fct. de Legendre associees paires avec m impair
Definition: type_parite.h:224
#define T_LEG_MI
fct. de Legendre associees avec m impair
Definition: type_parite.h:240
Lorene prototypes.
Definition: app_hor.h:67
#define TRA_T
Translation en Theta, used for a bitwise shift (in hex)
Definition: type_parite.h:160
#define MSQ_P
Extraction de l&#39;info sur Phi.
Definition: type_parite.h:156
#define T_COS
dev. cos seulement
Definition: type_parite.h:196
int get_etat() const
Gives the logical state.
Definition: tbl.h:414
#define T_SIN
dev. sin seulement
Definition: type_parite.h:198
const Mg3d * get_mg() const
Returns the Mg3d on which the this is defined.
Definition: valeur.h:763
#define T_LEG_I
fct. de Legendre associees impaires
Definition: type_parite.h:220
#define T_COS_I
dev. cos seulement, harmoniques impaires
Definition: type_parite.h:204
#define T_COS_P
dev. cos seulement, harmoniques paires
Definition: type_parite.h:200
#define MSQ_T
Extraction de l&#39;info sur Theta.
Definition: type_parite.h:154
#define T_COSSIN_C
dev. cos-sin alternes, cos pour m=0
Definition: type_parite.h:192
Base_val base
Bases on which the spectral expansion is performed.
Definition: valeur.h:315
double * t
The array of double.
Definition: tbl.h:176
int * b
Array (size: nzone ) of the spectral basis in each domain.
Definition: base_val.h:334
#define MSQ_R
Extraction de l&#39;info sur R.
Definition: type_parite.h:152
#define T_SIN_P
dev. sin seulement, harmoniques paires
Definition: type_parite.h:202
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:465
#define T_LEG_IP
fct. de Legendre associees impaires avec m pair
Definition: type_parite.h:222
int get_nr(int l) const
Returns the number of points in the radial direction ( ) in domain no. l.
Definition: grilles.h:469
#define T_LEG_P
fct. de Legendre associees paires
Definition: type_parite.h:216
#define NONDEF
base inconnue
Definition: type_parite.h:150
#define T_COSSIN_CI
cos impair-sin pair alternes, cos pour m=0
Definition: type_parite.h:212
int get_taille() const
Gives the total size (ie dim.taille)
Definition: tbl.h:417
Base_val base
Bases of the spectral expansions.
Definition: mtbl_cf.h:210
int etat
Logical state (ETATNONDEF , ETATQCQ or ETATZERO ).
Definition: valeur.h:305
#define P_COSSIN_I
dev. sur Phi = 2*phi, freq. impaires
Definition: type_parite.h:249
Basic array class.
Definition: tbl.h:164
int get_nt(int l) const
Returns the number of points in the co-latitude direction ( ) in domain no. l.
Definition: grilles.h:474
#define P_COSSIN_P
dev. sur Phi = 2*phi, freq. paires
Definition: type_parite.h:247
#define T_LEG_II
fct. de Legendre associees impaires avec m impair
Definition: type_parite.h:226
#define T_COSSIN_CP
cos pair-sin impair alternes, cos pour m=0
Definition: type_parite.h:208
#define MAX_BASE
Nombre max. de bases differentes.
Definition: type_parite.h:144
#define T_SIN_I
dev. sin seulement, harmoniques impaires
Definition: type_parite.h:206
Tbl ** t
Array (size nzone ) of pointers on the Tbl &#39;s which contain the spectral coefficients in each domain...
Definition: mtbl_cf.h:215
#define T_LEG_PP
fct. de Legendre associees paires avec m pair
Definition: type_parite.h:218
#define T_COSSIN_S
dev. cos-sin alternes, sin pour m=0
Definition: type_parite.h:194