LORENE
base_val.h
1/*
2 * Definition of Lorene class Base_val
3 *
4 */
5
6/*
7 * Copyright (c) 1999-2000 Jean-Alain Marck
8 * Copyright (c) 1999-2001 Eric Gourgoulhon
9 * Copyright (c) 1999-2001 Philippe Grandclement
10 * Copyright (c) 2001 Jerome Novak
11 *
12 * This file is part of LORENE.
13 *
14 * LORENE is free software; you can redistribute it and/or modify
15 * it under the terms of the GNU General Public License as published by
16 * the Free Software Foundation; either version 2 of the License, or
17 * (at your option) any later version.
18 *
19 * LORENE is distributed in the hope that it will be useful,
20 * but WITHOUT ANY WARRANTY; without even the implied warranty of
21 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
22 * GNU General Public License for more details.
23 *
24 * You should have received a copy of the GNU General Public License
25 * along with LORENE; if not, write to the Free Software
26 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
27 *
28 */
29
30
31#ifndef __BASE_VAL_H_
32#define __BASE_VAL_H_
33
34/*
35 * $Id: base_val.h,v 1.23 2014/10/13 08:52:31 j_novak Exp $
36 * $Log: base_val.h,v $
37 * Revision 1.23 2014/10/13 08:52:31 j_novak
38 * Lorene classes and functions now belong to the namespace Lorene.
39 *
40 * Revision 1.22 2014/10/06 15:09:39 j_novak
41 * Modified #include directives to use c++ syntax.
42 *
43 * Revision 1.21 2013/06/05 15:10:41 j_novak
44 * Suppression of FINJAC sampling in r. This Jacobi(0,2) base is now
45 * available by setting colloc_r to BASE_JAC02 in the Mg3d constructor.
46 *
47 * Revision 1.20 2013/01/11 08:20:10 j_novak
48 * New radial spectral bases with Legendre polynomials (R_LEG, R_LEGP, R_LEGI).
49 *
50 * Revision 1.19 2012/01/17 10:18:17 j_penner
51 * changed nzone from a private member to a protected member
52 *
53 * Revision 1.18 2009/10/23 12:55:46 j_novak
54 * New base T_LEG_MI
55 *
56 * Revision 1.17 2009/10/08 16:19:32 j_novak
57 * Addition of new bases T_COS and T_SIN.
58 *
59 * Revision 1.16 2008/12/03 15:21:20 j_novak
60 * New method mult_cost.
61 *
62 * Revision 1.15 2007/12/11 15:28:05 jl_cornou
63 * Jacobi(0,2) polynomials partially implemented
64 *
65 * Revision 1.14 2005/09/07 13:09:48 j_novak
66 * New method for determining the highest multipole that can be described on a 3D
67 * grid.
68 *
69 * Revision 1.13 2004/11/23 15:03:14 m_forot
70 * Corrected error in comments.
71 *
72 * Revision 1.12 2004/11/04 15:39:45 e_gourgoulhon
73 * Modif documentation: added description of bases without any symmetry
74 * in theta.
75 *
76 * Revision 1.11 2004/08/24 09:14:40 p_grandclement
77 * Addition of some new operators, like Poisson in 2d... It now requieres the
78 * GSL library to work.
79 *
80 * Also, the way a variable change is stored by a Param_elliptic is changed and
81 * no longer uses Change_var but rather 2 Scalars. The codes using that feature
82 * will requiere some modification. (It should concern only the ones about monopoles)
83 *
84 * Revision 1.10 2004/03/22 13:12:40 j_novak
85 * Modification of comments to use doxygen instead of doc++
86 *
87 * Revision 1.9 2004/01/27 14:31:25 j_novak
88 * New method Base_val::mult_sint()
89 *
90 * Revision 1.8 2004/01/27 14:13:58 j_novak
91 * Added method Base_val::mult_x()
92 *
93 * Revision 1.7 2003/10/20 06:41:43 e_gourgoulhon
94 * Corrected documentation.
95 *
96 * Revision 1.6 2003/10/19 19:42:50 e_gourgoulhon
97 * -- member nzone now private
98 * -- introduced new methods get_nzone() and get_b()
99 * -- introduced new methods name_r, name_theta and name_phi.
100 *
101 * Revision 1.5 2003/09/16 08:53:05 j_novak
102 * Addition of the T_LEG_II base (odd in theta, only for odd m) and the
103 * transformation functions from and to T_SIN_P.
104 *
105 * Revision 1.4 2002/10/16 14:36:28 j_novak
106 * Reorganization of #include instructions of standard C++, in order to
107 * use experimental version 3 of gcc.
108 *
109 * Revision 1.3 2002/09/13 09:17:31 j_novak
110 * Modif. commentaires
111 *
112 * Revision 1.2 2002/06/17 14:05:15 j_novak
113 * friend functions are now also declared outside the class definition
114 *
115 * Revision 1.1.1.1 2001/11/20 15:19:27 e_gourgoulhon
116 * LORENE
117 *
118 * Revision 2.20 2001/10/12 14:56:09 novak
119 * Ajout des methodes dsdx(), dsdt(), etc...
120 *
121 * Revision 2.19 2001/05/29 16:08:45 eric
122 * Modif commentaires (mise en conformite Doc++ 3.4.7).
123 *
124 * Revision 2.18 2000/09/28 10:19:50 eric
125 * Modif commentaires (nouvelles bases T_LEG_IP et T_LEG_PI).
126 *
127 * Revision 2.17 2000/09/08 11:42:55 eric
128 * *** empty log message ***
129 *
130 * Revision 2.16 2000/09/08 10:14:11 eric
131 * *** empty log message ***
132 *
133 * Revision 2.15 2000/09/08 10:11:49 eric
134 * *** empty log message ***
135 *
136 * Revision 2.14 2000/09/08 10:10:24 eric
137 * *** empty log message ***
138 *
139 * Revision 2.13 2000/09/08 10:06:19 eric
140 * Ajout des methodes set_base_r, etc... et get_base_r, etc...
141 *
142 * Revision 2.12 1999/12/29 10:49:12 eric
143 * theta_functions et phi_functions declarees const.
144 *
145 * Revision 2.11 1999/12/29 10:36:58 eric
146 * Modif commentaires.
147 *
148 * Revision 2.10 1999/12/28 12:57:44 eric
149 * Introduction des methodes theta_functions et phi_functions.
150 *
151 * Revision 2.9 1999/11/19 14:53:11 phil
152 * *** empty log message ***
153 *
154 * Revision 2.8 1999/11/18 13:48:35 phil
155 * *** empty log message ***
156 *
157 * Revision 2.7 1999/11/18 12:51:12 phil
158 * commentaire de operator*
159 *
160 * Revision 2.6 1999/10/26 13:23:12 phil
161 * *** empty log message ***
162 *
163 * Revision 2.5 1999/10/26 13:18:06 phil
164 * ajout de l'operator*
165 *
166 * Revision 2.4 1999/10/13 15:49:12 eric
167 * *** empty log message ***
168 *
169 * Revision 2.3 1999/10/12 10:02:17 eric
170 * Amelioration des commentaires: description de toutes les bases.
171 *
172 * Revision 2.2 1999/10/01 15:55:58 eric
173 * Depoussierage.
174 * Documentation
175 *
176 * Revision 2.1 1999/09/13 14:38:08 phil
177 * ajout de l'operateur ==
178 *
179 * Revision 2.0 1999/02/22 15:17:47 hyc
180 * *** empty log message ***
181 *
182 *
183 * $Header: /cvsroot/Lorene/C++/Include/base_val.h,v 1.23 2014/10/13 08:52:31 j_novak Exp $
184 *
185 */
186
187#include <cstdio>
188#include <cassert>
189#include "headcpp.h"
190
191
192#include "type_parite.h"
193namespace Lorene {
194class Tbl ;
195class Mg3d ;
196
322class Base_val {
323
324 // Data
325 // ----
326 protected:
327 int nzone ;
328
329 public:
331 int* b ;
332
333 // Constructors - Destructor
334 // -------------------------
335 public:
336 explicit Base_val(int nb_of_domains) ;
337 Base_val(const Base_val& ) ;
338
340 explicit Base_val(FILE* ) ;
341
342 ~Base_val() ;
343
344
345 // Mutators / assignment
346 // ---------------------
347 public:
348 void set_base_nondef() ;
349
359 void set_base_r(int l, int base_r) ;
360
369 void set_base_t(int base_t) ;
370
379 void set_base_p(int base_p) ;
380
381 void operator=(const Base_val& ) ;
382
383 // Accessors
384 // ---------
385 public:
386 bool operator==(const Base_val& ) const ;
387
389 int get_b(int l) const {
390 assert( (l>=0) && (l<nzone) ) ;
391 return b[l] ;
392 }
393
400 int get_base_r(int l) const {
401 assert( (l>=0) && (l<nzone) ) ;
402 return b[l] & MSQ_R ;
403 } ;
404
411 int get_base_t(int l) const {
412 assert( (l>=0) && (l<nzone) ) ;
413 return b[l] & MSQ_T ;
414 } ;
415
422 int get_base_p(int l) const {
423 assert( (l>=0) && (l<nzone) ) ;
424 return b[l] & MSQ_P ;
425 } ;
426
440 void name_r(int l, int k, int j, int i, char* basename) const ;
441
453 void name_theta(int l, int k, int j, char* basename) const ;
454
464 void name_phi(int l, int k, char* basename) const ;
465
466
489 const Tbl& theta_functions(int l, int nt) const ;
490
505 const Tbl& phi_functions(int l, int np) const ;
506
508 int get_nzone() const {return nzone; } ;
509
510 // Manipulation of basis
511 // ----------------------
512 public:
517 void dsdx() ;
518
523 void sx() ;
524
529 void mult_x() ;
530
535 void dsdt() ;
536
541 void ssint() ;
542
547 void mult_sint() ;
548
553 void mult_cost() ;
554
559 void ylm() ;
560
564 void give_quant_numbers (int, int, int,
565 int&, int&, int&) const ;
566
574 int give_lmax(const Mg3d& mgrid, int lz) const ;
575
576 // Outputs
577 // -------
578 public:
579 void sauve(FILE *) const ;
580
581 friend ostream& operator<<(ostream& , const Base_val& ) ;
582 friend Base_val operator*(const Base_val&, const Base_val&) ;
583};
584ostream& operator<<(ostream& , const Base_val& ) ;
602Base_val operator*(const Base_val&, const Base_val&) ;
603
606}
607#endif
Bases of the spectral expansions.
Definition base_val.h:322
void mult_x()
The basis is transformed as with a multiplication by .
void sauve(FILE *) const
Save in a file.
Definition base_val.C:228
int give_lmax(const Mg3d &mgrid, int lz) const
Returns the highest multipole for a given grid.
void name_r(int l, int k, int j, int i, char *basename) const
Name of the basis function in r ( )
void set_base_t(int base_t)
Sets the expansion basis for functions in all domains.
Definition base_val.C:195
void name_phi(int l, int k, char *basename) const
Name of the basis function in .
void dsdt()
The basis is transformed as with a operation.
int get_b(int l) const
Returns the code for the expansion basis in domain no. l.
Definition base_val.h:389
void sx()
The basis is transformed as with a multiplication.
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
void mult_sint()
The basis is transformed as with a multiplication.
void name_theta(int l, int k, int j, char *basename) const
Name of the basis function in .
friend ostream & operator<<(ostream &, const Base_val &)
Display.
Definition base_val.C:238
~Base_val()
Destructor.
Definition base_val.C:177
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
bool operator==(const Base_val &) const
Comparison operator.
Definition base_val.C:335
void give_quant_numbers(int, int, int, int &, int &, int &) const
Computes the various quantum numbers and 1d radial base.
void ssint()
The basis is transformed as with a multiplication.
int * b
Array (size: nzone ) of the spectral basis in each domain.
Definition base_val.h:331
friend Base_val operator*(const Base_val &, const Base_val &)
This operator is used when calling multiplication or division of Valeur .
const Tbl & phi_functions(int l, int np) const
Values of the phi basis functions at the phi collocation points.
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
const Tbl & theta_functions(int l, int nt) const
Values of the theta basis functions at the theta collocation points.
void ylm()
The basis is transformed as with a transformation to basis.
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
int get_nzone() const
Returns the number of domains.
Definition base_val.h:508
void dsdx()
The basis is transformed as with a operation.
void mult_cost()
The basis is transformed as with a multiplication.
Multi-domain grid.
Definition grilles.h:273
Basic array class.
Definition tbl.h:161
Base_val operator*(const Base_val &, const Base_val &)
This operator is used when calling multiplication or division of Valeur .
#define MSQ_R
Extraction de l'info sur R.
#define MSQ_T
Extraction de l'info sur Theta.
#define MSQ_P
Extraction de l'info sur Phi.
Lorene prototypes.
Definition app_hor.h:64