LORENE
base_val_manip.C
1/*
2 * Copyright (c) 2001 Jerome Novak
3 *
4 * This file is part of LORENE.
5 *
6 * LORENE is free software; you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation; either version 2 of the License, or
9 * (at your option) any later version.
10 *
11 * LORENE is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with LORENE; if not, write to the Free Software
18 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
19 *
20 */
21
22
23char base_val_manip_C[] = "$Header: /cvsroot/Lorene/C++/Source/Base_val/base_val_manip.C,v 1.12 2015/03/09 10:32:27 j_novak Exp $" ;
24
25/*
26 * $Id: base_val_manip.C,v 1.12 2015/03/09 10:32:27 j_novak Exp $
27 * $Log: base_val_manip.C,v $
28 * Revision 1.12 2015/03/09 10:32:27 j_novak
29 * Inclusion of r-Legendre bases.
30 *
31 * Revision 1.11 2014/10/13 08:52:38 j_novak
32 * Lorene classes and functions now belong to the namespace Lorene.
33 *
34 * Revision 1.10 2014/10/06 15:12:56 j_novak
35 * Modified #include directives to use c++ syntax.
36 *
37 * Revision 1.9 2009/10/08 16:20:13 j_novak
38 * Addition of new bases T_COS and T_SIN.
39 *
40 * Revision 1.8 2008/12/03 15:21:21 j_novak
41 * New method mult_cost.
42 *
43 * Revision 1.7 2008/02/18 13:53:38 j_novak
44 * Removal of special indentation instructions.
45 *
46 * Revision 1.6 2004/11/23 15:08:00 m_forot
47 * Added the bases for the cases without any equatorial symmetry
48 * (T_COSSIN_C, T_COSSIN_S, T_LEG, R_CHEBPI_P, R_CHEBPI_I).
49 *
50 * Revision 1.5 2004/01/27 14:31:26 j_novak
51 * New method Base_val::mult_sint()
52 *
53 * Revision 1.4 2004/01/27 14:13:59 j_novak
54 * Added method Base_val::mult_x()
55 *
56 * Revision 1.3 2003/09/16 08:54:09 j_novak
57 * Addition of the T_LEG_II base (odd in theta, only for odd m) and the
58 * transformation functions to and from the T_SIN_P base.
59 *
60 * Revision 1.2 2002/10/16 14:36:30 j_novak
61 * Reorganization of #include instructions of standard C++, in order to
62 * use experimental version 3 of gcc.
63 *
64 * Revision 1.1.1.1 2001/11/20 15:19:28 e_gourgoulhon
65 * LORENE
66 *
67 * Revision 1.3 2001/10/12 15:06:50 novak
68 * *** empty log message ***
69 *
70 * Revision 1.2 2001/10/12 15:05:14 novak
71 * *** empty log message ***
72 *
73 * Revision 1.1 2001/10/12 14:57:24 novak
74 * Initial revision
75 *
76 *
77 * $Header: /cvsroot/Lorene/C++/Source/Base_val/base_val_manip.C,v 1.12 2015/03/09 10:32:27 j_novak Exp $
78 *
79 */
80
81// Headers C
82#include <cstdlib>
83#include <cassert>
84
85// Headers Lorene
86#include "headcpp.h"
87#include "type_parite.h"
88#include "base_val.h"
89
90namespace Lorene {
92
93 switch(get_base_r(0)) {
94 case R_CHEBP:
96 break ;
97 case R_CHEBI:
99 break ;
100 case R_CHEBPIM_P:
102 break ;
103 case R_CHEBPIM_I:
105 break ;
106 case R_CHEBPI_P:
108 break ;
109 case R_CHEBPI_I:
111 break ;
112 case R_LEGP:
113 set_base_r(0, R_LEGI) ;
114 break ;
115 case R_LEGI:
116 set_base_r(0, R_LEGP) ;
117 break ;
118 default:
119 break ;
120 }
121 return ;
122}
123
125
126 switch(get_base_r(0)) {
127 case R_CHEBP:
128 set_base_r(0, R_CHEBI) ;
129 break ;
130 case R_CHEBI:
131 set_base_r(0, R_CHEBP) ;
132 break ;
133 case R_CHEBPIM_P:
135 break ;
136 case R_CHEBPIM_I:
138 break ;
139 case R_CHEBPI_P:
141 break ;
142 case R_CHEBPI_I:
144 break ;
145 case R_LEGP:
146 set_base_r(0, R_LEGI) ;
147 break ;
148 case R_LEGI:
149 set_base_r(0, R_LEGP) ;
150 break ;
151 default:
152 break ;
153 }
154 return ;
155}
156
158
159 switch(get_base_r(0)) {
160 case R_CHEBP:
161 set_base_r(0, R_CHEBI) ;
162 break ;
163 case R_CHEBI:
164 set_base_r(0, R_CHEBP) ;
165 break ;
166 case R_CHEBPIM_P:
168 break ;
169 case R_CHEBPIM_I:
171 break ;
172 case R_CHEBPI_P:
174 break ;
175 case R_CHEBPI_I:
177 break ;
178 case R_LEGP:
179 set_base_r(0, R_LEGI) ;
180 break ;
181 case R_LEGI:
182 set_base_r(0, R_LEGP) ;
183 break ;
184 default:
185 break ;
186 }
187 return ;
188}
189
191
192 switch(get_base_t(0)) {
193 case T_COS_P:
195 break ;
196 case T_COS_I:
198 break ;
199 case T_SIN_P:
201 break ;
202 case T_SIN_I:
204 break ;
205 case T_COSSIN_CP:
207 break ;
208 case T_COSSIN_SP:
210 break ;
211 case T_COSSIN_CI:
213 break ;
214 case T_COSSIN_SI:
216 break ;
217 case T_COSSIN_C:
219 break ;
220 case T_COSSIN_S:
222 break ;
223 case T_COS:
225 break ;
226 case T_SIN:
228 break ;
229 default:
230 cout << "Wrong base in Base_val::dsdt()!" << endl ;
231 abort() ;
232 exit(-1) ;
233 break ;
234 }
235 return ;
236}
237
239
240 switch(get_base_t(0)) {
241 case T_COS_P:
243 break ;
244 case T_COS_I:
246 break ;
247 case T_SIN_P:
249 break ;
250 case T_SIN_I:
252 break ;
253 case T_COSSIN_CP:
255 break ;
256 case T_COSSIN_SP:
258 break ;
259 case T_COSSIN_CI:
261 break ;
262 case T_COSSIN_SI:
264 break ;
265 case T_COSSIN_C:
267 break ;
268 case T_COSSIN_S:
270 break ;
271 case T_COS:
273 break ;
274 case T_SIN:
276 break ;
277 default:
278 cout << "Wrong base in Base_val::ssint()!" << endl ;
279 abort() ;
280 exit(-1) ;
281 break ;
282 }
283 return ;
284}
285
287
288 switch(get_base_t(0)) {
289 case T_COS_P:
291 break ;
292 case T_COS_I:
294 break ;
295 case T_SIN_P:
297 break ;
298 case T_SIN_I:
300 break ;
301 case T_COSSIN_CP:
303 break ;
304 case T_COSSIN_SP:
306 break ;
307 case T_COSSIN_CI:
309 break ;
310 case T_COSSIN_SI:
312 break ;
313 case T_COSSIN_C:
315 break ;
316 case T_COSSIN_S:
318 break ;
319 case T_COS:
321 break ;
322 case T_SIN:
324 break ;
325 default:
326 cout << "Wrong base in Base_val::mult_sint()!" << endl ;
327 abort() ;
328 exit(-1) ;
329 break ;
330 }
331 return ;
332}
333
335
336 switch(get_base_t(0)) {
337 case T_COS_P:
339 break ;
340 case T_COS_I:
342 break ;
343 case T_SIN_P:
345 break ;
346 case T_SIN_I:
348 break ;
349 case T_COSSIN_CP:
351 break ;
352 case T_COSSIN_SP:
354 break ;
355 case T_COSSIN_CI:
357 break ;
358 case T_COSSIN_SI:
360 break ;
361 case T_COSSIN_C:
363 break ;
364 case T_COSSIN_S:
366 break ;
367 case T_COS:
369 break ;
370 case T_SIN:
372 break ;
373 default:
374 cout << "Wrong base in Base_val::mult_cost()!" << endl ;
375 abort() ;
376 exit(-1) ;
377 break ;
378 }
379 return ;
380}
381
383
384 switch(get_base_t(0)) {
385 case T_COS_P:
387 break ;
388 case T_COS_I:
390 break ;
391 case T_SIN_I:
393 break ;
394 case T_SIN_P:
396 break ;
397 case T_COSSIN_CP:
399 break ;
400 case T_COSSIN_CI:
402 break ;
403 case T_COSSIN_C:
405 break ;
406 case T_COSSIN_S:
408 break ;
409 case T_COS:
411 break ;
412 default:
413 cout << "Wrong base in Base_val::ylm()!" << endl ;
414 abort() ;
415 exit(-1) ;
416 break ;
417 }
418 return ;
419}
420}
void mult_x()
The basis is transformed as with a multiplication by .
void set_base_t(int base_t)
Sets the expansion basis for functions in all domains.
Definition base_val.C:195
void dsdt()
The basis is transformed as with a operation.
void sx()
The basis is transformed as with a multiplication.
void mult_sint()
The basis is transformed as with a multiplication.
int get_base_r(int l) const
Returns the expansion basis for r ( ) functions in the domain of index l (e.g.
Definition base_val.h:400
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
int get_base_t(int l) const
Returns the expansion basis for functions in the domain of index l (e.g.
Definition base_val.h:411
void ssint()
The basis is transformed as with a multiplication.
void ylm()
The basis is transformed as with a transformation to basis.
void dsdx()
The basis is transformed as with a operation.
void mult_cost()
The basis is transformed as with a multiplication.
#define T_LEG_MP
fct. de Legendre associees avec m pair
#define R_LEGP
base de Legendre paire (rare) seulement
#define T_LEG_PI
fct. de Legendre associees paires avec m impair
#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 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 T_LEG_P
fct. de Legendre associees paires
#define T_LEG_IP
fct. de Legendre associees impaires avec m pair
#define T_COS_P
dev. cos seulement, harmoniques paires
#define T_COSSIN_CI
cos impair-sin pair alternes, cos pour m=0
#define R_CHEBPIM_P
Cheb. pair-impair suivant m, pair pour m=0.
#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 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 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 R_CHEBPI_P
Cheb. pair-impair suivant l pair pour l=0.
Lorene prototypes.
Definition app_hor.h:64