LORENE
mtbl_cf.h
1/*
2 * Definition of Lorene class Mtbl_cf
3 *
4 */
5
6/*
7 * Copyright (c) 1999-2000 Jean-Alain Marck
8 * Copyright (c) 1999-2005 Eric Gourgoulhon
9 * Copyright (c) 1999-2001 Philippe Grandclement
10 * Copyright (c) 1999-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 __MTBL_CF_H_
32#define __MTBL_CF_H_
33
34/*
35 * $Id: mtbl_cf.h,v 1.10 2014/10/13 08:52:36 j_novak Exp $
36 * $Log: mtbl_cf.h,v $
37 * Revision 1.10 2014/10/13 08:52:36 j_novak
38 * Lorene classes and functions now belong to the namespace Lorene.
39 *
40 * Revision 1.9 2006/06/06 14:56:59 j_novak
41 * Summation functions for angular coefficients at xi=+/-1.
42 *
43 * Revision 1.8 2005/04/04 21:30:42 e_gourgoulhon
44 * Added argument lambda to method poisson_angu
45 * to treat the generalized angular Poisson equation:
46 * Lap_ang u + lambda u = source.
47 *
48 * Revision 1.7 2004/03/22 13:12:42 j_novak
49 * Modification of comments to use doxygen instead of doc++
50 *
51 * Revision 1.6 2003/11/06 14:43:37 e_gourgoulhon
52 * Gave a name to const arguments in certain method prototypes (e.g.
53 * constructors) to correct a bug of DOC++.
54 *
55 * Revision 1.5 2003/10/19 19:44:41 e_gourgoulhon
56 * Introduced new method display (to replace the old affiche_seuil).
57 *
58 * Revision 1.4 2003/10/15 21:09:22 e_gourgoulhon
59 * Added method poisson_regu.
60 *
61 * Revision 1.3 2002/09/13 09:17:33 j_novak
62 * Modif. commentaires
63 *
64 * Revision 1.2 2002/06/17 14:05:17 j_novak
65 * friend functions are now also declared outside the class definition
66 *
67 * Revision 1.1.1.1 2001/11/20 15:19:27 e_gourgoulhon
68 * LORENE
69 *
70 * Revision 2.28 2000/09/11 13:52:21 eric
71 * Ajout des methodes mult_cp() et mult_sp().
72 *
73 * Revision 2.27 2000/08/16 10:42:54 eric
74 * Suppression du membre dzpuis.
75 *
76 * Revision 2.26 2000/08/04 11:41:52 eric
77 * Ajout de l'operateur (int l) et de la fonction set(int l) pour l'acces
78 * individuel aux Tbl.
79 *
80 * Revision 2.25 2000/03/06 10:26:24 eric
81 * Ajout des fonctions val_point_symy, val_point_asymy.
82 *
83 * Revision 2.24 2000/02/25 13:53:44 eric
84 * Suppression de la fonction nettoie().
85 *
86 * Revision 2.23 1999/12/29 13:11:34 eric
87 * Ajout de la fonction val_point_jk.
88 *
89 * Revision 2.22 1999/12/07 14:51:52 eric
90 * Changement ordre des arguments (phi,theta,xi) --> (xi,theta,phi)
91 * dans la routine val_point.
92 *
93 * Revision 2.21 1999/12/06 16:45:48 eric
94 * Ajout de la fonction val_point.
95 *
96 * Revision 2.20 1999/11/23 14:30:12 novak
97 * Ajout des membres mult_ct et mult_st
98 *
99 * Revision 2.19 1999/11/16 13:06:22 novak
100 * Ajout de mult_x et scost
101 *
102 * Revision 2.18 1999/10/29 15:07:10 eric
103 * Suppression des fonctions membres min() et max():
104 * elles deviennent des fonctions externes.
105 * Ajout de fonctions mathematiques (abs, norme, etc...).
106 *
107 * Revision 2.17 1999/10/21 13:42:00 eric
108 * *** empty log message ***
109 *
110 * Revision 2.16 1999/10/21 12:49:13 eric
111 * Ajout de la fonction membre nettoie().
112 *
113 * Revision 2.15 1999/10/18 15:06:50 eric
114 * La fonction membre annule() est rebaptisee annule_hard().
115 * Introduction de la fonction membre annule(int, int).
116 *
117 * Revision 2.14 1999/10/18 13:39:13 eric
118 * Suppression de l'argument base dans les routines de derivation.
119 *
120 * Revision 2.13 1999/10/13 15:49:22 eric
121 * Ajout du membre base.
122 * Modification des constructeurs (la base doit etre passee en argument).
123 *
124 * Revision 2.12 1999/10/01 14:49:19 eric
125 * Depoussierage.
126 * Documentation.
127 *
128 * Revision 2.11 1999/09/07 14:33:43 phil
129 * ajout de la fonction ssint(int*)
130 *
131 * Revision 2.10 1999/04/26 17:02:40 phil
132 * ajout de sx2(int*)
133 *
134 * Revision 2.9 1999/04/26 16:24:02 phil
135 * ajout de mult2_xm1_zec(int*)
136 *
137 * Revision 2.8 1999/04/26 16:12:27 phil
138 * ajout de mult_xm1_zec(int *)
139 *
140 * Revision 2.7 1999/04/26 15:48:43 phil
141 * ajout de sxm1_zec(int*)
142 *
143 * Revision 2.6 1999/04/26 14:50:30 phil
144 * ajout de sx(int*)
145 *
146 * Revision 2.5 1999/04/26 12:19:38 phil
147 * ajout lapang
148 *
149 * Revision 2.4 1999/03/03 10:32:44 hyc
150 * *** empty log message ***
151 *
152 * Revision 2.3 1999/03/02 18:55:00 eric
153 * Ajout de la fonction affiche_seuil.
154 *
155 *
156 * $Header: /cvsroot/Lorene/C++/Include/mtbl_cf.h,v 1.10 2014/10/13 08:52:36 j_novak Exp $
157 *
158 */
159
160// Headers Lorene
161#include "tbl.h"
162#include "base_val.h"
163#include "grilles.h"
164
165namespace Lorene {
166class Mg3d ;
167
186class Mtbl_cf {
187
188 // Data :
189 // -----
190 private:
192 const Mg3d* mg ;
194 int nzone ;
196 int etat ;
197
198 public:
201
205 Tbl** t ;
206
207 // Constructors - Destructor
208 // -------------------------
209
210 public:
212 Mtbl_cf(const Mg3d& mgrid, const Base_val& basis) ;
214 Mtbl_cf(const Mg3d* p_mgrid, const Base_val& basis) ;
215
217 Mtbl_cf(const Mg3d&, FILE* ) ;
218
220 Mtbl_cf(const Mtbl_cf& ) ;
222 ~Mtbl_cf() ;
223
224 // Assignement
225 // -----------
227 void operator=(const Mtbl_cf& ) ;
229 void operator=(double ) ;
231 void operator=(int ) ;
232
233 // Memory management
234 // -----------------
235 private:
239 void del_t() ;
240
241 public:
242
247 void set_etat_nondef() ;
248
253 void set_etat_zero() ;
254
261 void set_etat_qcq() ;
262
271 void annule_hard() ;
272
283 void annule(int l_min, int l_max) ;
284
285 // Access to individual elements
286 // -----------------------------
287 public:
288
294 Tbl& set(int l) {
295 assert(l < nzone) ;
296 assert(etat == ETATQCQ) ;
297 return *(t[l]) ;
298 };
299
300
306 const Tbl& operator()(int l) const {
307 assert(l < nzone) ;
308 assert(etat == ETATQCQ) ;
309 return *(t[l]) ;
310 };
311
312
319 double& set(int l, int k, int j, int i) {
320 assert(l < nzone) ;
321 assert(etat == ETATQCQ) ;
322 return (t[l])->set(k, j, i) ;
323 };
324
325
332 double operator()(int l, int k, int j, int i) const {
333 assert(etat != ETATNONDEF) ;
334 assert(l < nzone) ;
335 if (etat == ETATZERO) {
336 double zero = 0. ;
337 return zero ;
338 }
339 else return (*t[l])(k, j, i) ;
340 };
341
352 double val_point(int l, double x, double theta, double phi) const ;
353
365 double val_point_symy(int l, double x, double theta, double phi) const ;
366
378 double val_point_asymy(int l, double x, double theta, double phi) const ;
379
392 double val_point_jk(int l, double x, int j, int k) const ;
393
407 double val_point_jk_symy(int l, double x, int j, int k) const ;
408
422 double val_point_jk_asymy(int l, double x, int j, int k) const ;
423
433 double val_out_bound_jk(int l, int j, int k) const ;
434
445 double val_in_bound_jk(int l, int j, int k) const ;
446
447
448
449 // Extraction of information
450 // -------------------------
451 public:
453 const Mg3d* get_mg() const { return mg ; };
454
456 int get_etat() const { return etat ; };
457
459 int get_nzone() const { return nzone ; } ;
460
461
462 // Outputs
463 // -------
464 public:
465
467 void sauve(FILE *) const ;
468
476 void display(double threshold = 1.e-7, int precision = 4,
477 ostream& ostr = cout) const ;
478
485 void affiche_seuil(ostream& ostr, int precision = 4,
486 double threshold = 1.e-7) const ;
488 friend ostream& operator<<(ostream& , const Mtbl_cf& ) ;
489
490 // Member arithmetics
491 // ------------------
492 public:
494 void operator+=(const Mtbl_cf & ) ;
496 void operator-=(const Mtbl_cf & ) ;
498 void operator*=(double ) ;
500 void operator/=(double ) ;
501
502 // Linear operators
503 // ----------------
504 public:
506 void dsdx() ;
507
509 void d2sdx2() ;
510
515 void sx() ;
516
521 void sx2() ;
522
527 void mult_x() ;
528
532 void sxm1_zec() ;
533
537 void mult_xm1_zec() ;
538
542 void mult2_xm1_zec() ;
543
545 void dsdt() ;
546
548 void d2sdt2() ;
549
551 void ssint() ;
552
554 void scost() ;
555
557 void mult_ct() ;
558
560 void mult_st() ;
561
563 void dsdp() ;
564
566 void d2sdp2() ;
567
569 void mult_cp() ;
570
572 void mult_sp() ;
573
575 void lapang() ;
576
577 // PDE resolution
578 //---------------
579 public:
595 void poisson_angu(double lambda = 0) ;
596
597} ;
598ostream& operator<<(ostream& , const Mtbl_cf& ) ;
599
607Mtbl_cf operator+(const Mtbl_cf& ) ;
609Mtbl_cf operator-(const Mtbl_cf& ) ;
611Mtbl_cf operator+(const Mtbl_cf&, const Mtbl_cf& ) ;
613Mtbl_cf operator-(const Mtbl_cf&, const Mtbl_cf& ) ;
615Mtbl_cf operator*(const Mtbl_cf&, double ) ;
617Mtbl_cf operator*(double, const Mtbl_cf& ) ;
619Mtbl_cf operator*(const Mtbl_cf&, int ) ;
621Mtbl_cf operator*(int, const Mtbl_cf& ) ;
623Mtbl_cf operator/(const Mtbl_cf&, double ) ;
625Mtbl_cf operator/(const Mtbl_cf&, int ) ;
626
628Mtbl_cf abs(const Mtbl_cf& ) ;
629
635Tbl max(const Mtbl_cf& ) ;
636
642Tbl min(const Mtbl_cf& ) ;
643
649Tbl norme(const Mtbl_cf& ) ;
650
659Tbl diffrel(const Mtbl_cf& a, const Mtbl_cf& b) ;
660
669Tbl diffrelmax(const Mtbl_cf& a, const Mtbl_cf& b) ;
670
673}
674#endif
Bases of the spectral expansions.
Definition base_val.h:322
Multi-domain grid.
Definition grilles.h:273
Coefficients storage for the multi-domain spectral method.
Definition mtbl_cf.h:186
double val_point_jk(int l, double x, int j, int k) const
Computes the value of the field represented by *this at an arbitrary point in , but collocation point...
Base_val base
Bases of the spectral expansions.
Definition mtbl_cf.h:200
void del_t()
Logical destructor: dellocates the memory occupied by the Tbl array t .
Definition mtbl_cf.C:276
void set_etat_nondef()
Sets the logical state to ETATNONDEF (undefined).
Definition mtbl_cf.C:294
void mult2_xm1_zec()
Id (r sampling = RARE, FIN ) \ (r -sampling = UNSURR )
~Mtbl_cf()
Destructor.
Definition mtbl_cf.C:143
void mult_xm1_zec()
Id (r sampling = RARE, FIN ) \ (r -sampling = UNSURR )
void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition mtbl_cf.C:288
void poisson_angu(double lambda=0)
Resolution of the generalized angular Poisson equation.
Definition mtbl_cf_pde.C:83
double val_point_symy(int l, double x, double theta, double phi) const
Computes the value of the field represented by *this at an arbitrary point, by means of the spectral ...
double val_point_jk_asymy(int l, double x, int j, int k) const
Computes the value of the field represented by *this at an arbitrary point in , but collocation point...
void sx()
(r -sampling = RARE ) \ Id (r sampling = FIN ) \ (r -sampling = UNSURR )
Definition valeur_sx.C:153
int get_etat() const
Returns the logical state.
Definition mtbl_cf.h:456
const Tbl & operator()(int l) const
Read-only of the Tbl containing the coefficients in a given domain.
Definition mtbl_cf.h:306
double val_point_jk_symy(int l, double x, int j, int k) const
Computes the value of the field represented by *this at an arbitrary point in , but collocation point...
double & set(int l, int k, int j, int i)
Read/write of a particular element.
Definition mtbl_cf.h:319
void sxm1_zec()
Id (r sampling = RARE, FIN ) \ (r -sampling = UNSURR )
const Mg3d * mg
Pointer on the multi-grid Mgd3 on which this is defined.
Definition mtbl_cf.h:192
void operator*=(double)
*= double
Definition mtbl_cf.C:429
void annule(int l_min, int l_max)
Sets the Mtbl_cf to zero in some domains.
Definition mtbl_cf.C:332
void sauve(FILE *) const
Save in a file.
Definition mtbl_cf.C:204
Tbl & set(int l)
Read/write of the Tbl containing the coefficients in a given domain.
Definition mtbl_cf.h:294
void operator-=(const Mtbl_cf &)
-= Mtbl_cf
void sx2()
(r -sampling = RARE ) \ Id (r sampling = FIN ) \ (r -sampling = UNSURR )
Definition valeur_sx2.C:157
int get_nzone() const
Returns the number of zones (domains)
Definition mtbl_cf.h:459
void affiche_seuil(ostream &ostr, int precision=4, double threshold=1.e-7) const
Prints only the values greater than a given threshold.
Definition mtbl_cf.C:397
const Mg3d * get_mg() const
Returns the Mg3d on which the Mtbl_cf is defined.
Definition mtbl_cf.h:453
friend ostream & operator<<(ostream &, const Mtbl_cf &)
Display.
Definition mtbl_cf.C:366
void operator=(const Mtbl_cf &)
Assignement to another Mtbl_cf.
Definition mtbl_cf.C:220
int etat
Logical state (ETATNONDEF , ETATQCQ or ETATZERO ).
Definition mtbl_cf.h:196
Tbl ** t
Array (size nzone ) of pointers on the Tbl 's which contain the spectral coefficients in each domain.
Definition mtbl_cf.h:205
double val_point_asymy(int l, double x, double theta, double phi) const
Computes the value of the field represented by *this at an arbitrary point, by means of the spectral ...
void display(double threshold=1.e-7, int precision=4, ostream &ostr=cout) const
Prints the coefficients whose values are greater than a given threshold, as well as the corresponding...
void mult_x()
(r -sampling = RARE ) \ Id (r sampling = FIN ) \ (r -sampling = UNSURR )
double val_point(int l, double x, double theta, double phi) const
Computes the value of the field represented by *this at an arbitrary point, by means of the spectral ...
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition mtbl_cf.C:300
void operator+=(const Mtbl_cf &)
+= Mtbl_cf
void annule_hard()
Sets the Mtbl_cf to zero in a hard way.
Definition mtbl_cf.C:312
void operator/=(double)
/= double
Definition mtbl_cf.C:434
void lapang()
Angular Laplacian.
double val_out_bound_jk(int l, int j, int k) const
Computes the angular coefficient of index j,k of the field represented by *this at by means of the s...
double val_in_bound_jk(int l, int j, int k) const
Computes the angular coefficient of index j,k of the field represented by *this at by means of the s...
int nzone
Number of domains (zones)
Definition mtbl_cf.h:194
double operator()(int l, int k, int j, int i) const
Read-only of a particular element.
Definition mtbl_cf.h:332
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 .
Cmp operator-(const Cmp &)
- Cmp
Definition cmp_arithm.C:108
Cmp operator/(const Cmp &, const Cmp &)
Cmp / Cmp.
Definition cmp_arithm.C:457
Tbl diffrel(const Cmp &a, const Cmp &b)
Relative difference between two Cmp (norme version).
Definition cmp_math.C:504
Tbl norme(const Cmp &)
Sums of the absolute values of all the values of the Cmp in each domain.
Definition cmp_math.C:481
Tbl min(const Cmp &)
Minimum values of a Cmp in each domain.
Definition cmp_math.C:458
Tbl max(const Cmp &)
Maximum values of a Cmp in each domain.
Definition cmp_math.C:435
Cmp abs(const Cmp &)
Absolute value.
Definition cmp_math.C:410
Tbl diffrelmax(const Cmp &a, const Cmp &b)
Relative difference between two Cmp (max version).
Definition cmp_math.C:539
Cmp operator+(const Cmp &)
Definition cmp_arithm.C:104
Lorene prototypes.
Definition app_hor.h:64