LORENE
base_val.C
1/*
2 * Methods of class Base_val
3 *
4 * (see file base_val.h for documentation)
5 *
6 */
7
8/*
9 * Copyright (c) 1999-2000 Jean-Alain Marck
10 * Copyright (c) 1999-2001 Eric Gourgoulhon
11 * Copyright (c) 1999-2001 Philippe Grandclement
12 *
13 * This file is part of LORENE.
14 *
15 * LORENE is free software; you can redistribute it and/or modify
16 * it under the terms of the GNU General Public License as published by
17 * the Free Software Foundation; either version 2 of the License, or
18 * (at your option) any later version.
19 *
20 * LORENE is distributed in the hope that it will be useful,
21 * but WITHOUT ANY WARRANTY; without even the implied warranty of
22 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
23 * GNU General Public License for more details.
24 *
25 * You should have received a copy of the GNU General Public License
26 * along with LORENE; if not, write to the Free Software
27 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
28 *
29 */
30
31
32char base_val_C[] = "$Header: /cvsroot/Lorene/C++/Source/Base_val/base_val.C,v 1.17 2014/10/13 08:52:38 j_novak Exp $" ;
33
34/*
35 * $Id: base_val.C,v 1.17 2014/10/13 08:52:38 j_novak Exp $
36 * $Log: base_val.C,v $
37 * Revision 1.17 2014/10/13 08:52:38 j_novak
38 * Lorene classes and functions now belong to the namespace Lorene.
39 *
40 * Revision 1.16 2014/10/06 15:12:56 j_novak
41 * Modified #include directives to use c++ syntax.
42 *
43 * Revision 1.15 2013/01/11 08:20:11 j_novak
44 * New radial spectral bases with Legendre polynomials (R_LEG, R_LEGP, R_LEGI).
45 *
46 * Revision 1.14 2012/01/17 14:44:19 j_penner
47 * Modified phi variables to only use 16 integers in arrays
48 *
49 * Revision 1.13 2009/10/23 12:55:16 j_novak
50 * New base T_LEG_MI
51 *
52 * Revision 1.12 2009/10/08 16:20:13 j_novak
53 * Addition of new bases T_COS and T_SIN.
54 *
55 * Revision 1.11 2008/08/19 06:41:59 j_novak
56 * Minor modifications to avoid warnings with gcc 4.3. Most of them concern
57 * cast-type operations, and constant strings that must be defined as const char*
58 *
59 * Revision 1.10 2008/02/18 13:53:38 j_novak
60 * Removal of special indentation instructions.
61 *
62 * Revision 1.9 2007/12/11 15:28:09 jl_cornou
63 * Jacobi(0,2) polynomials partially implemented
64 *
65 * Revision 1.8 2004/12/17 13:35:01 m_forot
66 * Add the case T_LEG
67 *
68 * Revision 1.7 2004/11/04 15:19:02 e_gourgoulhon
69 * operator<< : added names R_CHEBPI_P, R_CHEBPI_I, T_COSSIN_C, T_COSSIN_S.
70 *
71 * Revision 1.6 2004/08/24 09:14:41 p_grandclement
72 * Addition of some new operators, like Poisson in 2d... It now requieres the
73 * GSL library to work.
74 *
75 * Also, the way a variable change is stored by a Param_elliptic is changed and
76 * no longer uses Change_var but rather 2 Scalars. The codes using that feature
77 * will requiere some modification. (It should concern only the ones about monopoles)
78 *
79 * Revision 1.5 2003/09/16 08:54:09 j_novak
80 * Addition of the T_LEG_II base (odd in theta, only for odd m) and the
81 * transformation functions to and from the T_SIN_P base.
82 *
83 * Revision 1.4 2002/11/13 15:05:59 j_novak
84 * Affichage de la base T_COS
85 *
86 * Revision 1.3 2002/10/16 14:36:30 j_novak
87 * Reorganization of #include instructions of standard C++, in order to
88 * use experimental version 3 of gcc.
89 *
90 * Revision 1.2 2001/12/04 21:27:52 e_gourgoulhon
91 *
92 * All writing/reading to a binary file are now performed according to
93 * the big endian convention, whatever the system is big endian or
94 * small endian, thanks to the functions fwrite_be and fread_be
95 *
96 * Revision 1.1.1.1 2001/11/20 15:19:28 e_gourgoulhon
97 * LORENE
98 *
99 * Revision 2.8 2000/09/28 10:20:19 eric
100 * Affichage: nouvelles bases T_LEG_IP et T_LEG_PI.
101 *
102 * Revision 2.7 2000/09/08 10:06:45 eric
103 * Ajout des methodes set_base_r, etc...
104 *
105 * Revision 2.6 1999/12/28 12:57:26 eric
106 * Reorganisation des headers.
107 *
108 * Revision 2.5 1999/10/12 10:02:51 eric
109 * Implementation de sauve().
110 * Modif de << : affichage du nom des bases et non plus de leur numero.
111 *
112 * Revision 2.4 1999/10/01 15:56:21 eric
113 * Depoussierage.
114 * Documentation.
115 *
116 * Revision 2.3 1999/09/13 14:57:06 phil
117 * ajout de l'operateur ==
118 *
119 * Revision 2.2 1999/03/02 15:22:23 eric
120 * Affichage des bases.
121 *
122 * Revision 2.1 1999/03/01 14:54:01 eric
123 * *** empty log message ***
124 *
125 * Revision 2.0 1999/02/22 15:19:20 hyc
126 * *** empty log message ***
127 *
128 *
129 * $Header: /cvsroot/Lorene/C++/Source/Base_val/base_val.C,v 1.17 2014/10/13 08:52:38 j_novak Exp $
130 *
131 */
132
133// Headers C
134#include <cstdio>
135#include <cassert>
136
137// Headers Lorene
138#include "headcpp.h"
139#include "type_parite.h"
140#include "base_val.h"
141#include "utilitaires.h"
142
143
144 //---------------//
145 // Constructeurs //
146 //---------------//
147
148// Constructeur
149namespace Lorene {
150Base_val::Base_val(int n) : nzone(n) {
151 b = new int[nzone] ;
152 for (int i=0 ; i<nzone ; i++) { // Boucle sur les zones
153 b[i] = NONDEF ;
154 }
155}
156
157// Copie
158Base_val::Base_val(const Base_val & bi) : nzone(bi.nzone) {
159 b = new int[nzone] ;
160 for (int i=0 ; i<nzone ; i++) { // Boucle sur les zones
161 b[i] = bi.b[i] ;
162 }
163}
164
165// From file
167 fread_be(&nzone, sizeof(int), 1, fd) ; // nzone
168 b = new int[nzone] ;
169 fread_be(b, sizeof(int), nzone, fd) ; // b[]
170}
171
172 //--------------//
173 // Destructeurs //
174 //--------------//
175
176// Destructeur
178 delete [] b ;
179}
180
181 //-------------//
182 // Affectation //
183 //-------------//
184
185void Base_val::set_base_r(int l, int base_r) {
186
187 assert( (l>=0) && (l<nzone) ) ;
188
189 int base_t = b[l] & MSQ_T ;
190 int base_p = b[l] & MSQ_P ;
191 b[l] = base_p | base_t | base_r ;
192
193}
194
195void Base_val::set_base_t(int base_t) {
196
197 for (int l=0; l<nzone; l++) {
198 int base_r = b[l] & MSQ_R ;
199 int base_p = b[l] & MSQ_P ;
200 b[l] = base_p | base_t | base_r ;
201 }
202}
203
204void Base_val::set_base_p(int base_p) {
205
206 for (int l=0; l<nzone; l++) {
207 int base_r = b[l] & MSQ_R ;
208 int base_t = b[l] & MSQ_T ;
209 b[l] = base_p | base_t | base_r ;
210 }
211}
212
213
214// From Base_val
216 // Protection
217 assert(nzone == bi.nzone) ;
218 for (int i=0 ; i<nzone ; i++) { // Boucle sur les zones
219 b[i] = bi.b[i] ;
220 }
221}
222
223 //------------//
224 // Sauvegarde //
225 //------------//
226
227// Save in a file
228void Base_val::sauve(FILE* fd) const {
229 fwrite_be(&nzone, sizeof(int), 1, fd) ; // nzone
230 fwrite_be(b, sizeof(int), nzone, fd) ; // b[]
231}
232
233 //------------//
234 // Impression //
235 //------------//
236
237// Operateurs <<
238ostream& operator<<(ostream& o, const Base_val & bi) {
239
240 static bool premier_appel = true ;
241 static const char* nom_r[MAX_BASE] ;
242 static const char* nom_t[MAX_BASE] ;
243 static const char* nom_p[MAX_BASE_2] ;
244
245 if (premier_appel) { // First call initializations
246
247 premier_appel = false ;
248
249 for (int i=0; i<MAX_BASE; i++) {
250 nom_r[i] = "UNDEFINED" ;
251 nom_t[i] = "UNDEFINED" ;
252 if(i%2==0){
253 nom_p[i/2] = "UNDEFINED" ; // saves a loop
254 }
255 }
256
257 nom_r[R_CHEB >> TRA_R] = "R_CHEB " ;
258 nom_r[R_CHEBU >> TRA_R] = "R_CHEBU " ;
259 nom_r[R_CHEBP >> TRA_R] = "R_CHEBP " ;
260 nom_r[R_CHEBI >> TRA_R] = "R_CHEBI " ;
261 nom_r[R_CHEBPIM_P >> TRA_R] = "R_CHEBPIM_P" ;
262 nom_r[R_CHEBPIM_I >> TRA_R] = "R_CHEBPIM_I" ;
263 nom_r[R_CHEBPI_P >> TRA_R] = "R_CHEBPI_P " ;
264 nom_r[R_CHEBPI_I >> TRA_R] = "R_CHEBPI_I " ;
265 nom_r[R_LEG >> TRA_R] = "R_LEG " ;
266 nom_r[R_LEGP >> TRA_R] = "R_LEGP " ;
267 nom_r[R_LEGI >> TRA_R] = "R_LEGI " ;
268 nom_r[R_JACO02 >> TRA_R] = "R_JACO02 " ;
269
270 nom_t[T_COS >> TRA_T] = "T_COS " ;
271 nom_t[T_SIN >> TRA_T] = "T_SIN " ;
272 nom_t[T_COS_P >> TRA_T] = "T_COS_P " ;
273 nom_t[T_COS_I >> TRA_T] = "T_COS_I " ;
274 nom_t[T_SIN_P >> TRA_T] = "T_SIN_P " ;
275 nom_t[T_SIN_I >> TRA_T] = "T_SIN_I " ;
276 nom_t[T_COSSIN_CP >> TRA_T] = "T_COSSIN_CP" ;
277 nom_t[T_COSSIN_SI >> TRA_T] = "T_COSSIN_SI" ;
278 nom_t[T_COSSIN_SP >> TRA_T] = "T_COSSIN_SP" ;
279 nom_t[T_COSSIN_CI >> TRA_T] = "T_COSSIN_CI" ;
280 nom_t[T_COSSIN_C >> TRA_T] = "T_COSSIN_C " ;
281 nom_t[T_COSSIN_S >> TRA_T] = "T_COSSIN_S " ;
282 nom_t[T_LEG >> TRA_T] = "T_LEG " ;
283 nom_t[T_LEG_MP >> TRA_T] = "T_LEG_MP " ;
284 nom_t[T_LEG_MI >> TRA_T] = "T_LEG_MI " ;
285 nom_t[T_LEG_P >> TRA_T] = "T_LEG_P " ;
286 nom_t[T_LEG_PP >> TRA_T] = "T_LEG_PP " ;
287 nom_t[T_LEG_I >> TRA_T] = "T_LEG_I " ;
288 nom_t[T_LEG_IP >> TRA_T] = "T_LEG_IP " ;
289 nom_t[T_LEG_PI >> TRA_T] = "T_LEG_PI " ;
290 nom_t[T_LEG_II >> TRA_T] = "T_LEG_II " ;
291 nom_t[T_CL_COS_P >> TRA_T] = "T_CL_COS_P " ;
292 nom_t[T_CL_SIN_P >> TRA_T] = "T_CL_SIN_P " ;
293 nom_t[T_CL_COS_I >> TRA_T] = "T_CL_COS_I " ;
294 nom_t[T_CL_SIN_I >> TRA_T] = "T_CL_SIN_I " ;
295
296 nom_p[P_COSSIN >> TRA_P] = "P_COSSIN " ;
297 nom_p[P_COSSIN_P >> TRA_P] = "P_COSSIN_P " ;
298 nom_p[P_COSSIN_I >> TRA_P] = "P_COSSIN_I " ;
299
300
301 } // End of first call operations
302
303
304 // Intro - Nombre de zones
305 int nzone = bi.nzone ;
306 o << "Bases of spectral expansions: " ;
307 for (int l=0 ; l<nzone ; l++) {
308 int base_r = (bi.b[l] & MSQ_R) >> TRA_R ;
309 int base_t = (bi.b[l] & MSQ_T) >> TRA_T ;
310 int base_p = (bi.b[l] & MSQ_P) >> TRA_P ;
311 o << endl ;
312 o << "Domain #" << l << " : r: " << nom_r[base_r]
313 << ", theta: " << nom_t[base_t]
314 << ", phi: " << nom_p[base_p] ;
315 }
316 o << endl ;
317
318 //Termine
319 return o ;
320}
321
322 //----------------------//
323 // Manipulation de base //
324 //----------------------//
325
327 for (int i=0 ; i<nzone ; i++) {
328 b[i] = NONDEF ;
329 }
330}
331 //----------------------//
332 // operateur logique //
333 //----------------------//
334
335bool Base_val::operator== (const Base_val& c2) const {
336
337 return (*b == *c2.b) ;
338}
339}
Bases of the spectral expansions.
Definition base_val.h:322
Base_val(int nb_of_domains)
Standard constructor.
Definition base_val.C:150
void sauve(FILE *) const
Save in a file.
Definition base_val.C:228
void set_base_t(int base_t)
Sets the expansion basis for functions in all domains.
Definition base_val.C:195
~Base_val()
Destructor.
Definition base_val.C:177
void set_base_r(int l, int base_r)
Sets the expansion basis for r ( ) functions in a given domain.
Definition base_val.C:185
bool operator==(const Base_val &) const
Comparison operator.
Definition base_val.C:335
int * b
Array (size: nzone ) of the spectral basis in each domain.
Definition base_val.h:331
int nzone
Number of domains (zones)
Definition base_val.h:327
void set_base_nondef()
Sets the spectral bases to NONDEF.
Definition base_val.C:326
void operator=(const Base_val &)
Assignment.
Definition base_val.C:215
void set_base_p(int base_p)
Sets the expansion basis for functions in all domains.
Definition base_val.C:204
#define MAX_BASE_2
Smaller maximum bases used for phi (and higher dimensions for now)
#define T_CL_COS_I
CL of odd cosines.
#define T_CL_SIN_I
CL of odd sines.
#define TRA_T
Translation en Theta, used for a bitwise shift (in hex)
#define T_LEG_MP
fct. de Legendre associees avec m pair
#define P_COSSIN_P
dev. sur Phi = 2*phi, freq. paires
#define R_LEGP
base de Legendre paire (rare) seulement
#define MAX_BASE
Nombre max. de bases differentes.
#define T_LEG_PI
fct. de Legendre associees paires avec m impair
#define R_CHEBU
base de Chebychev ordinaire (fin), dev. en 1/r
#define R_JACO02
base de Jacobi(0,2) ordinaire (finjac)
#define T_COSSIN_SP
sin pair-cos impair alternes, sin pour m=0
#define R_LEGI
base de Legendre impaire (rare) seulement
#define R_CHEBI
base de Cheb. impaire (rare) seulement
#define T_LEG
fct. de Legendre associees
#define T_SIN_P
dev. sin seulement, harmoniques paires
#define MSQ_R
Extraction de l'info sur R.
#define T_COSSIN_S
dev. cos-sin alternes, sin pour m=0
#define R_CHEBPIM_I
Cheb. pair-impair suivant m, impair pour m=0.
#define T_COSSIN_SI
sin impair-cos pair alternes, sin pour m=0
#define R_CHEBPI_I
Cheb. pair-impair suivant l impair pour l=0.
#define R_LEG
base de Legendre ordinaire (fin)
#define T_CL_SIN_P
CL of even sines.
#define T_LEG_P
fct. de Legendre associees paires
#define T_LEG_IP
fct. de Legendre associees impaires avec m pair
#define NONDEF
base inconnue
#define T_COS_P
dev. cos seulement, harmoniques paires
#define T_COSSIN_CI
cos impair-sin pair alternes, cos pour m=0
#define P_COSSIN
dev. standart
#define P_COSSIN_I
dev. sur Phi = 2*phi, freq. impaires
#define R_CHEBPIM_P
Cheb. pair-impair suivant m, pair pour m=0.
#define TRA_R
Translation en R, used for a bitwise shift (in hex)
#define T_CL_COS_P
CL of even cosines.
#define T_LEG_MI
fct. de Legendre associees avec m impair
#define MSQ_T
Extraction de l'info sur Theta.
#define T_COSSIN_CP
cos pair-sin impair alternes, cos pour m=0
#define T_LEG_II
fct. de Legendre associees impaires avec m impair
#define R_CHEB
base de Chebychev ordinaire (fin)
#define T_LEG_I
fct. de Legendre associees impaires
#define T_SIN_I
dev. sin seulement, harmoniques impaires
#define T_COS
dev. cos seulement
#define R_CHEBP
base de Cheb. paire (rare) seulement
#define MSQ_P
Extraction de l'info sur Phi.
#define T_SIN
dev. sin seulement
#define T_LEG_PP
fct. de Legendre associees paires avec m pair
#define T_COS_I
dev. cos seulement, harmoniques impaires
#define T_COSSIN_C
dev. cos-sin alternes, cos pour m=0
#define TRA_P
Translation en Phi, used for a bitwise shift (in hex)
#define R_CHEBPI_P
Cheb. pair-impair suivant l pair pour l=0.
int fread_be(int *aa, int size, int nb, FILE *fich)
Reads integer(s) from a binary file according to the big endian convention.
Definition fread_be.C:69
int fwrite_be(const int *aa, int size, int nb, FILE *fich)
Writes integer(s) into a binary file according to the big endian convention.
Definition fwrite_be.C:70
Lorene prototypes.
Definition app_hor.h:64