LORENE
base_val_quantum.C
1 /*
2  * Copyright (c) 2004 Philippe Grandclement
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 
23 char base_val_quantum_C[] = "$Header: /cvsroot/Lorene/C++/Source/Base_val/base_val_quantum.C,v 1.10 2014/10/13 08:52:39 j_novak Exp $" ;
24 
25 /*
26  * $Id: base_val_quantum.C,v 1.10 2014/10/13 08:52:39 j_novak Exp $
27  * $Log: base_val_quantum.C,v $
28  * Revision 1.10 2014/10/13 08:52:39 j_novak
29  * Lorene classes and functions now belong to the namespace Lorene.
30  *
31  * Revision 1.9 2014/10/06 15:12:57 j_novak
32  * Modified #include directives to use c++ syntax.
33  *
34  * Revision 1.8 2013/01/11 08:20:11 j_novak
35  * New radial spectral bases with Legendre polynomials (R_LEG, R_LEGP, R_LEGI).
36  *
37  * Revision 1.7 2009/10/23 12:55:16 j_novak
38  * New base T_LEG_MI
39  *
40  * Revision 1.6 2009/10/08 16:20:13 j_novak
41  * Addition of new bases T_COS and T_SIN.
42  *
43  * Revision 1.5 2007/12/11 15:28:09 jl_cornou
44  * Jacobi(0,2) polynomials partially implemented
45  *
46  * Revision 1.4 2005/09/07 13:09:50 j_novak
47  * New method for determining the highest multipole that can be described on a 3D
48  * grid.
49  *
50  * Revision 1.3 2004/11/23 15:08:01 m_forot
51  * Added the bases for the cases without any equatorial symmetry
52  * (T_COSSIN_C, T_COSSIN_S, T_LEG, R_CHEBPI_P, R_CHEBPI_I).
53  *
54  * Revision 1.2 2004/08/30 16:27:59 r_prix
55  * added #include <stdlib.h> (got ERROR 'abort' is undefined without this...)
56  *
57  * Revision 1.1 2004/08/24 09:14:41 p_grandclement
58  * Addition of some new operators, like Poisson in 2d... It now requieres the
59  * GSL library to work.
60  *
61  * Also, the way a variable change is stored by a Param_elliptic is changed and
62  * no longer uses Change_var but rather 2 Scalars. The codes using that feature
63  * will requiere some modification. (It should concern only the ones about monopoles)
64  *
65  *
66  * $Header: /cvsroot/Lorene/C++/Source/Base_val/base_val_quantum.C,v 1.10 2014/10/13 08:52:39 j_novak Exp $
67  *
68  */
69 
70 // Headers C
71 #include <cstdio>
72 #include <cassert>
73 #include <cstdlib>
74 
75 // Headers Lorene
76 #include "grilles.h"
77 #include "base_val.h"
78 #include "utilitaires.h"
79 
80 namespace Lorene {
81 void Base_val::give_quant_numbers (int l, int k, int j,
82  int& m_quant, int& l_quant, int& base_r_1d) const {
83 
84  int base_p = get_base_p(l) ;
85  int base_t = get_base_t(l) ;
86  int base_r = get_base_r(l) ;
87 
88  switch (base_p) {
89  case P_COSSIN :
90  m_quant = k/2 ;
91  break;
92 
93  case P_COSSIN_P :
94  if (k%2 == 0)
95  m_quant = k ;
96  else
97  m_quant = (k-1) ;
98  break;
99 
100  case P_COSSIN_I :
101  m_quant = 2*( (k-1) / 2) + 1 ;
102  break;
103  default:
104  cout << "Unknown basis in phi in give_quant_numbers ..." << endl ;
105  abort() ;
106  break ;
107  }
108 
109  switch (base_t) {
110  case T_COS_P :
111  l_quant = 2*j ;
112  break;
113 
114  case T_SIN_P :
115  l_quant = 2*j ;
116  break;
117 
118  case T_COS_I :
119  l_quant = 2*j+1 ;
120  break;
121 
122  case T_SIN_I :
123  l_quant = 2*j+1 ;
124  break;
125 
126  case T_COSSIN_CP :
127  if (m_quant%2 == 0)
128  l_quant = 2*j ;
129  else
130  l_quant = 2*j+1 ;
131  break ;
132 
133  case T_COSSIN_SP :
134  if (m_quant%2 == 0)
135  l_quant = 2*j ;
136  else
137  l_quant = 2*j+1 ;
138  break ;
139 
140  case T_COSSIN_CI :
141  if (m_quant%2 == 0)
142  l_quant = 2*j+1 ;
143  else
144  l_quant = 2*j ;
145  break ;
146 
147  case T_COSSIN_SI :
148  if (m_quant%2 == 0)
149  l_quant = 2*j+1 ;
150  else
151  l_quant = 2*j ;
152  break ;
153 
154  case T_COSSIN_C :
155  l_quant = j ;
156  break ;
157 
158  case T_COSSIN_S :
159  l_quant = j ;
160  break ;
161 
162  case T_COS :
163  l_quant = j ;
164  break ;
165 
166  case T_LEG_P :
167  if (m_quant%2 == 0)
168  l_quant = 2*j ;
169  else
170  l_quant = 2*j+1 ;
171  break ;
172 
173  case T_LEG_PP :
174  l_quant = 2*j ;
175  break ;
176 
177  case T_LEG_I :
178  if (m_quant%2 == 0)
179  l_quant = 2*j+1 ;
180  else
181  l_quant = 2*j ;
182  break ;
183 
184  case T_LEG_IP :
185  l_quant = 2*j+1 ;
186  break ;
187 
188  case T_LEG_PI :
189  l_quant = 2*j+1 ;
190  break ;
191 
192  case T_LEG_II :
193  l_quant = 2*j ;
194  break ;
195 
196  case T_LEG :
197  l_quant = j ;
198  break ;
199 
200  case T_LEG_MP :
201  l_quant = j ;
202  break ;
203 
204  case T_LEG_MI :
205  l_quant = j ;
206  break ;
207 
208  case T_CL_COS_P:
209  l_quant = 2*j ;
210  break ;
211 
212  case T_CL_SIN_P:
213  l_quant = 2*j ;
214  break ;
215 
216  case T_CL_COS_I:
217  l_quant = 2*j+1 ;
218  break ;
219 
220  case T_CL_SIN_I:
221  l_quant = 2*j+1 ;
222  break ;
223 
224  default:
225  cout << "Unknown basis in theta in give_quant_numbers ..." << endl ;
226  abort() ;
227  break ;
228  }
229 
230  switch (base_r) {
231  case R_CHEB :
232  base_r_1d = R_CHEB ;
233  break ;
234 
235  case R_CHEBP :
236  base_r_1d = R_CHEBP ;
237  break ;
238 
239  case R_CHEBI :
240  base_r_1d = R_CHEBI ;
241  break ;
242 
243  case R_LEG :
244  base_r_1d = R_LEG ;
245  break ;
246 
247  case R_LEGP :
248  base_r_1d = R_LEGP ;
249  break ;
250 
251  case R_LEGI :
252  base_r_1d = R_LEGI ;
253  break ;
254 
255  case R_JACO02 :
256  base_r_1d = R_JACO02 ;
257  break ;
258 
259  case R_CHEBPIM_P :
260  if (m_quant%2 == 0)
261  base_r_1d = R_CHEBP ;
262  else
263  base_r_1d = R_CHEBI ;
264  break ;
265  case R_CHEBPIM_I :
266  if (m_quant%2 == 0)
267  base_r_1d = R_CHEBI ;
268  else
269  base_r_1d = R_CHEBP ;
270  break ;
271  case R_CHEBPI_P :
272  if (l_quant%2 == 0)
273  base_r_1d = R_CHEBP ;
274  else
275  base_r_1d = R_CHEBI ;
276  break ;
277  case R_CHEBPI_I :
278  if (l_quant%2 == 0)
279  base_r_1d = R_CHEBI ;
280  else
281  base_r_1d = R_CHEBP ;
282  break ;
283  case R_CHEBU :
284  base_r_1d = R_CHEBU ;
285  break ;
286 
287  default:
288  cout << "Unknown basis in r in give_quant_numbers ..." << endl ;
289  abort() ;
290  break ;
291  }
292 }
293 
294 int Base_val::give_lmax(const Mg3d& mgrid, int lz) const {
295 
296 #ifndef NDEBUG
297  int nz = mgrid.get_nzone() ;
298  assert (lz < nz) ;
299 #endif
300 
301  int ntm1 = mgrid.get_nt(lz) - 1;
302  int base_t = get_base_t(lz) ;
303  bool m_odd = (mgrid.get_np(lz) > 2) ;
304 
305  int l_max = 0 ;
306 
307  switch (base_t) {
308  case T_COS_P :
309  l_max = 2*ntm1 ;
310  break;
311 
312  case T_SIN_P :
313  l_max = 2*ntm1 ;
314  break;
315 
316  case T_COS_I :
317  l_max = 2*ntm1+1 ;
318  break;
319 
320  case T_SIN_I :
321  l_max = 2*ntm1+1 ;
322  break;
323 
324  case T_COSSIN_CP :
325  if (!m_odd)
326  l_max = 2*ntm1 ;
327  else
328  l_max = 2*ntm1+1 ;
329  break ;
330 
331  case T_COSSIN_SP :
332  if (!m_odd)
333  l_max = 2*ntm1 ;
334  else
335  l_max = 2*ntm1+1 ;
336  break ;
337 
338  case T_COSSIN_CI :
339  if (!m_odd)
340  l_max = 2*ntm1+1 ;
341  else
342  l_max = 2*ntm1 ;
343  break ;
344 
345  case T_COSSIN_SI :
346  if (!m_odd)
347  l_max = 2*ntm1+1 ;
348  else
349  l_max = 2*ntm1 ;
350  break ;
351 
352  case T_COSSIN_C :
353  l_max = ntm1 ;
354  break ;
355 
356  case T_COSSIN_S :
357  l_max = ntm1 ;
358  break ;
359 
360  case T_COS :
361  l_max = ntm1 ;
362  break ;
363 
364  case T_LEG_P :
365  if (!m_odd)
366  l_max = 2*ntm1 ;
367  else
368  l_max = 2*ntm1+1 ;
369  break ;
370 
371  case T_LEG_PP :
372  l_max = 2*ntm1 ;
373  break ;
374 
375  case T_LEG_I :
376  if (!m_odd)
377  l_max = 2*ntm1+1 ;
378  else
379  l_max = 2*ntm1 ;
380  break ;
381 
382  case T_LEG_IP :
383  l_max = 2*ntm1+1 ;
384  break ;
385 
386  case T_LEG_PI :
387  l_max = 2*ntm1+1 ;
388  break ;
389 
390  case T_LEG_II :
391  l_max = 2*ntm1 ;
392  break ;
393 
394  case T_LEG :
395  l_max = ntm1 ;
396  break ;
397 
398  case T_LEG_MP :
399  l_max = ntm1 ;
400  break ;
401 
402  case T_LEG_MI :
403  l_max = ntm1 ;
404  break ;
405 
406  case T_CL_COS_P:
407  l_max = 2*ntm1 ;
408  break ;
409 
410  case T_CL_SIN_P:
411  l_max = 2*ntm1 ;
412  break ;
413 
414  case T_CL_COS_I:
415  l_max = 2*ntm1+1 ;
416  break ;
417 
418  case T_CL_SIN_I:
419  l_max = 2*ntm1+1 ;
420  break ;
421 
422  default:
423  cout << "Unknown basis in theta in Base_val::get_lmax ..."
424  << endl ;
425  abort() ;
426  break ;
427  }
428  return l_max ;
429 }
430 }
P_COSSIN
#define P_COSSIN
dev. standart
Definition: type_parite.h:245
R_CHEB
#define R_CHEB
base de Chebychev ordinaire (fin)
Definition: type_parite.h:166
Lorene::Mg3d::get_np
int get_np(int l) const
Returns the number of points in the azimuthal direction ( ) in domain no. l.
Definition: grilles.h:462
T_LEG_MI
#define T_LEG_MI
fct. de Legendre associees avec m impair
Definition: type_parite.h:240
T_LEG_MP
#define T_LEG_MP
fct. de Legendre associees avec m pair
Definition: type_parite.h:238
Lorene::Base_val::get_base_p
int get_base_p(int l) const
Returns the expansion basis for functions in the domain of index l (e.g.
Definition: base_val.h:422
Lorene::Mg3d
Multi-domain grid.
Definition: grilles.h:273
T_LEG_PI
#define T_LEG_PI
fct. de Legendre associees paires avec m impair
Definition: type_parite.h:224
T_COS
#define T_COS
dev. cos seulement
Definition: type_parite.h:196
T_COSSIN_SP
#define T_COSSIN_SP
sin pair-cos impair alternes, sin pour m=0
Definition: type_parite.h:210
R_CHEBPI_I
#define R_CHEBPI_I
Cheb. pair-impair suivant l impair pour l=0.
Definition: type_parite.h:174
Lorene::Mg3d::get_nt
int get_nt(int l) const
Returns the number of points in the co-latitude direction ( ) in domain no. l.
Definition: grilles.h:457
Lorene
Lorene prototypes.
Definition: app_hor.h:64
R_JACO02
#define R_JACO02
base de Jacobi(0,2) ordinaire (finjac)
Definition: type_parite.h:188
T_CL_COS_P
#define T_CL_COS_P
CL of even cosines.
Definition: type_parite.h:228
R_LEGP
#define R_LEGP
base de Legendre paire (rare) seulement
Definition: type_parite.h:184
T_COS_I
#define T_COS_I
dev. cos seulement, harmoniques impaires
Definition: type_parite.h:204
R_LEGI
#define R_LEGI
base de Legendre impaire (rare) seulement
Definition: type_parite.h:186
T_LEG_I
#define T_LEG_I
fct. de Legendre associees impaires
Definition: type_parite.h:220
R_CHEBP
#define R_CHEBP
base de Cheb. paire (rare) seulement
Definition: type_parite.h:168
R_CHEBI
#define R_CHEBI
base de Cheb. impaire (rare) seulement
Definition: type_parite.h:170
T_CL_SIN_I
#define T_CL_SIN_I
CL of odd sines.
Definition: type_parite.h:234
Lorene::Base_val::give_quant_numbers
void give_quant_numbers(int, int, int, int &, int &, int &) const
Computes the various quantum numbers and 1d radial base.
Definition: base_val_quantum.C:81
T_COSSIN_C
#define T_COSSIN_C
dev. cos-sin alternes, cos pour m=0
Definition: type_parite.h:192
T_CL_SIN_P
#define T_CL_SIN_P
CL of even sines.
Definition: type_parite.h:230
T_SIN_P
#define T_SIN_P
dev. sin seulement, harmoniques paires
Definition: type_parite.h:202
T_COS_P
#define T_COS_P
dev. cos seulement, harmoniques paires
Definition: type_parite.h:200
T_COSSIN_SI
#define T_COSSIN_SI
sin impair-cos pair alternes, sin pour m=0
Definition: type_parite.h:214
R_CHEBPIM_I
#define R_CHEBPIM_I
Cheb. pair-impair suivant m, impair pour m=0.
Definition: type_parite.h:178
Lorene::Mg3d::get_nzone
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:448
T_LEG_P
#define T_LEG_P
fct. de Legendre associees paires
Definition: type_parite.h:216
T_LEG_IP
#define T_LEG_IP
fct. de Legendre associees impaires avec m pair
Definition: type_parite.h:222
R_CHEBPIM_P
#define R_CHEBPIM_P
Cheb. pair-impair suivant m, pair pour m=0.
Definition: type_parite.h:176
R_CHEBPI_P
#define R_CHEBPI_P
Cheb. pair-impair suivant l pair pour l=0.
Definition: type_parite.h:172
Lorene::Base_val::get_base_r
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
P_COSSIN_I
#define P_COSSIN_I
dev. sur Phi = 2*phi, freq. impaires
Definition: type_parite.h:249
T_COSSIN_CI
#define T_COSSIN_CI
cos impair-sin pair alternes, cos pour m=0
Definition: type_parite.h:212
T_COSSIN_CP
#define T_COSSIN_CP
cos pair-sin impair alternes, cos pour m=0
Definition: type_parite.h:208
P_COSSIN_P
#define P_COSSIN_P
dev. sur Phi = 2*phi, freq. paires
Definition: type_parite.h:247
T_COSSIN_S
#define T_COSSIN_S
dev. cos-sin alternes, sin pour m=0
Definition: type_parite.h:194
T_LEG_PP
#define T_LEG_PP
fct. de Legendre associees paires avec m pair
Definition: type_parite.h:218
R_CHEBU
#define R_CHEBU
base de Chebychev ordinaire (fin), dev. en 1/r
Definition: type_parite.h:180
T_LEG_II
#define T_LEG_II
fct. de Legendre associees impaires avec m impair
Definition: type_parite.h:226
R_LEG
#define R_LEG
base de Legendre ordinaire (fin)
Definition: type_parite.h:182
T_LEG
#define T_LEG
fct. de Legendre associees
Definition: type_parite.h:236
T_CL_COS_I
#define T_CL_COS_I
CL of odd cosines.
Definition: type_parite.h:232
Lorene::Base_val::get_base_t
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
T_SIN_I
#define T_SIN_I
dev. sin seulement, harmoniques impaires
Definition: type_parite.h:206
Lorene::Base_val::give_lmax
int give_lmax(const Mg3d &mgrid, int lz) const
Returns the highest multipole for a given grid.
Definition: base_val_quantum.C:294