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