LORENE
scalar.h
1/*
2 * Definition of Lorene class Scalar
3 *
4 */
5
6/*
7 * Copyright (c) 2003-2004 Eric Gourgoulhon & Jerome Novak
8 *
9 * Copyright (c) 1999-2000 Jean-Alain Marck (for previous class Cmp)
10 * Copyright (c) 1999-2002 Eric Gourgoulhon (for previous class Cmp)
11 * Copyright (c) 1999-2001 Philippe Grandclement (for previous class Cmp)
12 * Copyright (c) 2000-2002 Jerome Novak (for previous class Cmp)
13 * Copyright (c) 2000-2001 Keisuke Taniguchi (for previous class Cmp)
14 *
15 * This file is part of LORENE.
16 *
17 * LORENE is free software; you can redistribute it and/or modify
18 * it under the terms of the GNU General Public License as published by
19 * the Free Software Foundation; either version 2 of the License, or
20 * (at your option) any later version.
21 *
22 * LORENE is distributed in the hope that it will be useful,
23 * but WITHOUT ANY WARRANTY; without even the implied warranty of
24 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
25 * GNU General Public License for more details.
26 *
27 * You should have received a copy of the GNU General Public License
28 * along with LORENE; if not, write to the Free Software
29 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
30 *
31 */
32
33
34#ifndef __SCALAR_H_
35#define __SCALAR_H_
36
37
38/*
39 * $Id: scalar.h,v 1.94 2015/12/18 15:52:51 j_novak Exp $
40 * $Log: scalar.h,v $
41 * Revision 1.94 2015/12/18 15:52:51 j_novak
42 * New method is_nan() for class Scalar.
43 *
44 * Revision 1.93 2015/09/10 13:28:05 j_novak
45 * New methods for the class Hot_Eos
46 *
47 * Revision 1.92 2014/10/13 08:52:36 j_novak
48 * Lorene classes and functions now belong to the namespace Lorene.
49 *
50 * Revision 1.91 2013/06/05 15:06:10 j_novak
51 * Legendre bases are treated as standard bases, when the multi-grid
52 * (Mg3d) is built with BASE_LEG.
53 *
54 * Revision 1.90 2013/01/11 15:44:54 j_novak
55 * Addition of Legendre bases (part 2).
56 *
57 * Revision 1.89 2012/12/19 13:59:56 j_penner
58 * added a few lines to the documentation for scalar_match_tau function
59 *
60 * Revision 1.88 2012/01/17 15:05:46 j_penner
61 * *** empty log message ***
62 *
63 * Revision 1.87 2012/01/17 10:16:27 j_penner
64 * functions added: sarra_filter_r, sarra_filter_r_all_domains, Heaviside
65 *
66 * Revision 1.86 2011/04/08 13:13:09 e_gourgoulhon
67 * Changed the comment of function val_point to indicate specifically the
68 * division by r^dzpuis in the compactified external domain.
69 *
70 * Revision 1.85 2008/09/29 13:23:51 j_novak
71 * Implementation of the angular mapping associated with an affine
72 * mapping. Things must be improved to take into account the domain index.
73 *
74 * Revision 1.84 2008/09/22 19:08:01 j_novak
75 * New methods to deal with boundary conditions
76 *
77 * Revision 1.83 2008/05/24 15:05:22 j_novak
78 * New method Scalar::match_tau to match the output of an explicit time-marching scheme with the tau method.
79 *
80 * Revision 1.82 2007/12/21 16:06:16 j_novak
81 * Methods to filter Tensor, Vector and Sym_tensor objects.
82 *
83 * Revision 1.81 2007/10/31 10:33:11 j_novak
84 * Added exponential filters to smooth Gibbs-type phenomena.
85 *
86 * Revision 1.80 2007/06/21 19:56:36 k_taniguchi
87 * Introduction of another filtre_r.
88 *
89 * Revision 1.79 2007/05/06 10:48:08 p_grandclement
90 * Modification of a few operators for the vorton project
91 *
92 * Revision 1.78 2007/01/16 15:05:59 n_vasset
93 * New constructor (taking a Scalar in mono-domain angular grid for
94 * boundary) for function sol_elliptic_boundary
95 *
96 * Revision 1.77 2006/05/26 09:00:09 j_novak
97 * New members for multiplication or division by cos(theta).
98 *
99 * Revision 1.76 2005/11/30 13:48:06 e_gourgoulhon
100 * Replaced M_PI/2 by 1.57... in argument list of sol_elliptic_sin_zec
101 * (in order not to require the definition of M_PI).
102 *
103 * Revision 1.75 2005/11/30 11:09:03 p_grandclement
104 * Changes for the Bin_ns_bh project
105 *
106 * Revision 1.74 2005/11/17 15:29:46 e_gourgoulhon
107 * Added arithmetics with Mtbl.
108 *
109 * Revision 1.73 2005/10/25 08:56:34 p_grandclement
110 * addition of std_spectral_base in the case of odd functions near the origin
111 *
112 * Revision 1.72 2005/09/07 13:10:47 j_novak
113 * Added a filter setting to zero all mulitpoles in a given range.
114 *
115 * Revision 1.71 2005/08/26 14:02:38 p_grandclement
116 * Modification of the elliptic solver that matches with an oscillatory exterior solution
117 * small correction in Poisson tau also...
118 *
119 * Revision 1.70 2005/08/25 12:14:07 p_grandclement
120 * Addition of a new method to solve the scalar Poisson equation, based on a multi-domain Tau-method
121 *
122 * Revision 1.69 2005/06/09 07:56:25 f_limousin
123 * Implement a new function sol_elliptic_boundary() and
124 * Vector::poisson_boundary(...) which solve the vectorial poisson
125 * equation (method 6) with an inner boundary condition.
126 *
127 * Revision 1.68 2005/06/08 12:35:18 j_novak
128 * New method for solving divergence-like ODEs.
129 *
130 * Revision 1.67 2005/05/20 14:42:30 j_novak
131 * Added the method Scalar::get_spectral_base().
132 *
133 * Revision 1.66 2005/04/04 21:28:57 e_gourgoulhon
134 * Added argument lambda to method poisson_angu
135 * to deal with the generalized angular Poisson equation:
136 * Lap_ang u + lambda u = source.
137 *
138 * Revision 1.65 2004/12/14 09:09:39 f_limousin
139 * Modif. comments.
140 *
141 * Revision 1.64 2004/11/23 12:41:53 f_limousin
142 * Intoduce function poisson_dir_neu(...) to solve a scalar poisson
143 * equation with a mixed boundary condition (Dirichlet + Neumann).
144 *
145 * Revision 1.63 2004/10/11 15:09:00 j_novak
146 * The radial manipulation functions take Scalar as arguments, instead of Cmp.
147 * Added a conversion operator from Scalar to Cmp.
148 * The Cmp radial manipulation function make conversion to Scalar, call to the
149 * Map_radial version with a Scalar argument and back.
150 *
151 * Revision 1.62 2004/08/24 09:14:40 p_grandclement
152 * Addition of some new operators, like Poisson in 2d... It now requieres the
153 * GSL library to work.
154 *
155 * Also, the way a variable change is stored by a Param_elliptic is changed and
156 * no longer uses Change_var but rather 2 Scalars. The codes using that feature
157 * will requiere some modification. (It should concern only the ones about monopoles)
158 *
159 * Revision 1.61 2004/07/27 08:24:26 j_novak
160 * Modif. comments
161 *
162 * Revision 1.60 2004/07/26 16:02:21 j_novak
163 * Added a flag to specify whether the primitive should be zero either at r=0
164 * or at r going to infinity.
165 *
166 * Revision 1.59 2004/07/06 13:36:27 j_novak
167 * Added methods for desaliased product (operator |) only in r direction.
168 *
169 * Revision 1.58 2004/06/22 08:49:57 p_grandclement
170 * Addition of everything needed for using the logarithmic mapping
171 *
172 * Revision 1.57 2004/06/14 15:24:23 e_gourgoulhon
173 * Added method primr (radial primitive).
174 *
175 * Revision 1.56 2004/06/11 14:29:56 j_novak
176 * Scalar::multipole_spectrum() is now a const method.
177 *
178 * Revision 1.55 2004/05/24 14:07:31 e_gourgoulhon
179 * Method set_domain now includes a call to del_deriv() for safety.
180 *
181 * Revision 1.54 2004/05/07 11:26:10 f_limousin
182 * New method filtre_r(int*)
183 *
184 * Revision 1.53 2004/03/22 13:12:43 j_novak
185 * Modification of comments to use doxygen instead of doc++
186 *
187 * Revision 1.52 2004/03/17 15:58:47 p_grandclement
188 * Slight modification of sol_elliptic_no_zec
189 *
190 * Revision 1.51 2004/03/11 12:07:30 e_gourgoulhon
191 * Added method visu_section_anim.
192 *
193 * Revision 1.50 2004/03/08 15:45:38 j_novak
194 * Modif. comment
195 *
196 * Revision 1.49 2004/03/05 15:09:40 e_gourgoulhon
197 * Added method smooth_decay.
198 *
199 * Revision 1.48 2004/03/01 09:57:02 j_novak
200 * the wave equation is solved with Scalars. It now accepts a grid with a
201 * compactified external domain, which the solver ignores and where it copies
202 * the values of the field from one time-step to the next.
203 *
204 * Revision 1.47 2004/02/27 09:43:58 f_limousin
205 * New methods filtre_phi(int) and filtre_theta(int).
206 *
207 * Revision 1.46 2004/02/26 22:46:26 e_gourgoulhon
208 * Added methods derive_cov, derive_con and derive_lie.
209 *
210 * Revision 1.45 2004/02/21 17:03:49 e_gourgoulhon
211 * -- Method "point" renamed "val_grid_point".
212 * -- Method "set_point" renamed "set_grid_point".
213 *
214 * Revision 1.44 2004/02/19 22:07:35 e_gourgoulhon
215 * Added argument "comment" in method spectral_display.
216 *
217 * Revision 1.43 2004/02/11 09:47:44 p_grandclement
218 * Addition of a new elliptic solver, matching with the homogeneous solution
219 * at the outer shell and not solving in the external domain (more details
220 * coming soon ; check your local Lorene dealer...)
221 *
222 * Revision 1.42 2004/01/28 16:46:22 p_grandclement
223 * Addition of the sol_elliptic_fixe_der_zero stuff
224 *
225 * Revision 1.41 2004/01/28 13:25:40 j_novak
226 * The ced_mult_r arguments have been suppressed from the Scalar::*dsd* methods.
227 * In the div/mult _r_dzpuis, there is no more default value.
228 *
229 * Revision 1.40 2004/01/28 10:39:17 j_novak
230 * Comments modified.
231 *
232 * Revision 1.39 2004/01/27 15:10:01 j_novak
233 * New methods Scalar::div_r_dzpuis(int) and Scalar_mult_r_dzpuis(int)
234 * which replace div_r_inc*. Tried to clean the dzpuis handling.
235 * WARNING: no testing at this point!!
236 *
237 * Revision 1.38 2004/01/23 13:25:44 e_gourgoulhon
238 * Added methods set_inner_boundary and set_outer_boundary.
239 * Methods set_val_inf and set_val_hor, which are particular cases of
240 * the above, have been suppressed.
241 *
242 * Revision 1.37 2004/01/22 16:10:09 e_gourgoulhon
243 * Added (provisory) method div_r_inc().
244 *
245 * Revision 1.36 2003/12/16 06:32:20 e_gourgoulhon
246 * Added method visu_box.
247 *
248 * Revision 1.35 2003/12/14 21:46:35 e_gourgoulhon
249 * Added argument start_dx in visu_section.
250 *
251 * Revision 1.34 2003/12/11 16:19:38 e_gourgoulhon
252 * Added method visu_section for visualization with OpenDX.
253 *
254 * Revision 1.33 2003/12/11 14:48:47 p_grandclement
255 * Addition of ALL (and that is a lot !) the files needed for the general elliptic solver ... UNDER DEVELOPEMENT...
256 *
257 * Revision 1.32 2003/11/13 13:43:53 p_grandclement
258 * Addition of things needed for Bhole::update_metric (const Etoile_bin&, double, double)
259 *
260 * Revision 1.31 2003/11/06 14:43:37 e_gourgoulhon
261 * Gave a name to const arguments in certain method prototypes (e.g.
262 * constructors) to correct a bug of DOC++.
263 *
264 * Revision 1.30 2003/11/04 22:55:50 e_gourgoulhon
265 * Added new methods mult_cost(), mult_sint() and div_sint().
266 *
267 * Revision 1.29 2003/10/29 13:09:11 e_gourgoulhon
268 * -- Added integer argument to derivative functions dsdr, etc...
269 * so that one can choose the dzpuis of the result (default=2).
270 * -- Change of method name: laplacien --> laplacian.
271 *
272 * Revision 1.28 2003/10/29 11:00:42 e_gourgoulhon
273 * Virtual functions dec_dzpuis and inc_dzpuis have now an integer argument to
274 * specify by which amount dzpuis is to be increased.
275 * Accordingly virtual methods dec2_dzpuis and inc2_dzpuis have been suppressed.
276 *
277 * Revision 1.27 2003/10/20 14:26:02 j_novak
278 * New assignement operators.
279 *
280 * Revision 1.26 2003/10/19 19:46:33 e_gourgoulhon
281 * -- Method spectral_display now virtual (from Tensor), list of argument
282 * changed.
283 *
284 * Revision 1.25 2003/10/17 13:46:14 j_novak
285 * The argument is now between 1 and 3 (instead of 0->2)
286 *
287 * Revision 1.24 2003/10/16 15:23:41 e_gourgoulhon
288 * Name of method div_r_ced() changed to div_r_inc2().
289 * Name of method div_rsint_ced() changed to div_rsint_inc2().
290 *
291 * Revision 1.23 2003/10/15 21:10:11 e_gourgoulhon
292 * Added method poisson_angu().
293 *
294 * Revision 1.22 2003/10/15 16:03:35 j_novak
295 * Added the angular Laplace operator for Scalar.
296 *
297 * Revision 1.21 2003/10/15 10:29:05 e_gourgoulhon
298 * Added new members p_dsdt and p_stdsdp.
299 * Added new methods dsdt(), stdsdp() and div_tant().
300 *
301 * Revision 1.20 2003/10/13 13:52:39 j_novak
302 * Better managment of derived quantities.
303 *
304 * Revision 1.19 2003/10/10 15:57:27 j_novak
305 * Added the state one (ETATUN) to the class Scalar
306 *
307 * Revision 1.18 2003/10/08 14:24:08 j_novak
308 * replaced mult_r_zec with mult_r_ced
309 *
310 * Revision 1.17 2003/10/06 16:16:02 j_novak
311 * New constructor from a Tensor.
312 *
313 * Revision 1.16 2003/10/06 13:58:45 j_novak
314 * The memory management has been improved.
315 * Implementation of the covariant derivative with respect to the exact Tensor
316 * type.
317 *
318 * Revision 1.15 2003/10/05 21:06:31 e_gourgoulhon
319 * - Added new methods div_r_ced() and div_rsint_ced().
320 * - Added new virtual method std_spectral_base()
321 * - Removed method std_spectral_base_scal()
322 *
323 * Revision 1.14 2003/10/01 13:02:58 e_gourgoulhon
324 * Suppressed the constructor from Map* .
325 *
326 * Revision 1.13 2003/09/29 12:52:56 j_novak
327 * Methods for changing the triad are implemented.
328 *
329 * Revision 1.12 2003/09/25 09:33:36 j_novak
330 * Added methods for integral calculation and various manipulations
331 *
332 * Revision 1.11 2003/09/25 09:11:21 e_gourgoulhon
333 * Added functions for radial operations (divr, etc...)
334 *
335 * Revision 1.10 2003/09/25 08:55:23 e_gourgoulhon
336 * Added members raccord*.
337 *
338 * Revision 1.9 2003/09/25 08:50:11 j_novak
339 * Added the members import
340 *
341 * Revision 1.8 2003/09/25 08:13:51 j_novak
342 * Added method for calculating derivatives
343 *
344 * Revision 1.7 2003/09/25 07:59:26 e_gourgoulhon
345 * Added prototypes for PDE resolutions.
346 *
347 * Revision 1.6 2003/09/25 07:17:58 j_novak
348 * Method asymptot implemented.
349 *
350 * Revision 1.5 2003/09/24 20:53:38 e_gourgoulhon
351 * Added -- constructor by conversion from a Cmp
352 * -- assignment from Cmp
353 *
354 * Revision 1.4 2003/09/24 15:10:54 j_novak
355 * Suppression of the etat flag in class Tensor (still present in Scalar)
356 *
357 * Revision 1.3 2003/09/24 12:01:44 j_novak
358 * Added friend functions for math.
359 *
360 * Revision 1.2 2003/09/24 10:22:01 e_gourgoulhon
361 * still in progress...
362 *
363 * Revision 1.1 2003/09/22 12:50:47 e_gourgoulhon
364 * First version: not ready yet!
365 *
366 *
367 * $Header: /cvsroot/Lorene/C++/Include/scalar.h,v 1.94 2015/12/18 15:52:51 j_novak Exp $
368 *
369 */
370
371// Headers Lorene
372#include "valeur.h"
373#include "tensor.h"
374
375namespace Lorene {
376class Param ;
377class Cmp ;
378class Param_elliptic ;
379
387class Scalar : public Tensor {
388
389 // Data :
390 // -----
391 protected:
392
396 int etat ;
397
403 int dzpuis ;
404
406
407 // Derived data :
408 // ------------
409 protected:
411 mutable Scalar* p_dsdr ;
412
416 mutable Scalar* p_srdsdt ;
417
422
424 mutable Scalar* p_dsdt ;
425
429 mutable Scalar* p_stdsdp ;
430
434 mutable Scalar* p_dsdx ;
435
439 mutable Scalar* p_dsdy ;
440
444 mutable Scalar* p_dsdz ;
445
448 mutable Scalar* p_lap ;
449
452 mutable Scalar* p_lapang ;
453
455 mutable Scalar* p_dsdradial ;
456
458 mutable Scalar* p_dsdrho ;
459
463 mutable int ind_lap ;
464
468 mutable Tbl* p_integ ;
469
470 // Constructors - Destructor
471 // -------------------------
472
473 public:
474
475 explicit Scalar(const Map& mpi) ;
476
478 Scalar(const Tensor& a) ;
479
480 Scalar(const Scalar& a) ;
481
482 Scalar(const Cmp& a) ;
483
485 Scalar(const Map&, const Mg3d&, FILE* ) ;
486
487 virtual ~Scalar() ;
488
489
490 // Memory management
491 // -----------------
492 protected:
493 void del_t() ;
494 virtual void del_deriv() const;
495 void set_der_0x0() const;
496
497 public:
498
504 virtual void set_etat_nondef() ;
505
511 virtual void set_etat_zero() ;
512
519 virtual void set_etat_qcq() ;
520
526 void set_etat_one() ;
527
536 virtual void allocate_all() ;
537
546 void annule_hard() ;
547
548 // Extraction of information
549 // -------------------------
550 public:
554 int get_etat() const {return etat;} ;
555
557 int get_dzpuis() const {return dzpuis;} ;
558
562 bool dz_nonzero() const ;
563
568 bool check_dzpuis(int dzi) const ;
569
575 bool is_nan(bool verbose=false) const ;
576
577 // Assignment
578 // -----------
579 public:
581 void operator=(const Scalar& a) ;
582
584 virtual void operator=(const Tensor& a) ;
585
586 void operator=(const Cmp& a) ;
587 void operator=(const Valeur& a) ;
588 void operator=(const Mtbl& a) ;
589 void operator=(double ) ;
590 void operator=(int ) ;
591
592 // Conversion oprators
593 //---------------------
594 operator Cmp() const ;
595
596 // Access to individual elements
597 // -----------------------------
598 public:
599
601 const Valeur& get_spectral_va() const {return va;} ;
602
604 Valeur& set_spectral_va() {return va;} ;
605
615 Tbl& set_domain(int l) {
616 assert(etat == ETATQCQ) ;
617 del_deriv() ;
618 return va.set(l) ;
619 };
620
625 const Tbl& domain(int l) const {
626 assert( (etat == ETATQCQ) || (etat == ETATUN) ) ;
627 return va(l) ;
628 };
629
630
637 double val_grid_point(int l, int k, int j, int i) const {
638 assert(etat != ETATNONDEF) ;
639 if (etat == ETATZERO) {
640 double zero = 0. ;
641 return zero ;
642 }
643 else {
644 if (etat == ETATUN) {
645 double one = 1. ;
646 return one ;
647 }
648 else{
649 return va(l, k, j, i) ;
650 }
651 }
652 };
653
668 double val_point(double r, double theta, double phi) const ;
669
670
684 double& set_grid_point(int l, int k, int j, int i) {
685 assert(etat == ETATQCQ) ;
686 return va.set(l, k, j, i) ;
687 };
688
689
700 virtual void annule(int l_min, int l_max) ;
701
707 void set_inner_boundary(int l, double x) ;
708
714 void set_outer_boundary(int l, double x) ;
715
724 Tbl multipole_spectrum () const ;
725
732 Tbl tbl_out_bound(int l_dom, bool leave_ylm = false) ;
733
740 Tbl tbl_in_bound(int n, bool leave_ylm = false) ;
741
748 Scalar scalar_out_bound(int n, bool leave_ylm = false) ;
749
750 // Differential operators and others
751 // ---------------------------------
752 public:
757 const Scalar& dsdr() const ;
758
763 const Scalar& srdsdt() const ;
764
769 const Scalar& srstdsdp() const ;
770
773 const Scalar& dsdt() const ;
774
782 const Scalar& dsdradial() const ;
783
788 const Scalar& dsdrho() const ;
789
792 const Scalar& stdsdp() const ;
793
799 const Scalar& dsdx() const ;
800
806 const Scalar& dsdy() const ;
807
813 const Scalar& dsdz() const ;
814
821 const Scalar& deriv(int i) const ;
822
827 const Vector& derive_cov(const Metric& gam) const ;
828
829
835 const Vector& derive_con(const Metric& gam) const ;
836
838 Scalar derive_lie(const Vector& v) const ;
839
840
849 const Scalar& laplacian(int ced_mult_r = 4) const ;
850
858 const Scalar& lapang() const ;
859
861 void div_r() ;
862
867 void div_r_dzpuis(int ced_mult_r) ;
868
872 void div_r_ced() ;
873
875 void mult_r() ;
876
881 void mult_r_dzpuis(int ced_mult_r) ;
882
886 void mult_r_ced() ;
887
889 void mult_rsint() ;
890
895 void mult_rsint_dzpuis(int ced_mult_r) ;
896
898 void div_rsint() ;
899
904 void div_rsint_dzpuis(int ced_mult_r) ;
905
906 void mult_cost() ;
907
908 void div_cost() ;
909
910 void mult_sint() ;
911
912 void div_sint() ;
913
914 void div_tant() ;
915
926 Scalar primr(bool null_infty = true) const ;
927
935 double integrale() const ;
936
946 const Tbl& integrale_domains() const ;
947
952 virtual void dec_dzpuis(int dec = 1) ;
953
958 virtual void inc_dzpuis(int inc = 1) ;
959
963 virtual void change_triad(const Base_vect& new_triad) ;
964
969 void filtre (int n) ;
970
975 void filtre_r (int* nn) ;
976
980 void filtre_r (int n, int nzone) ;
981
992// virtual void exponential_filter_r(int lzmin, int lzmax, int p,
993// double alpha= -16.) ;
994 virtual void exponential_filter_r(int lzmin, int lzmax, int p,
995 double alpha= -16.) ;
996
1007 void sarra_filter_r(int lzmin, int lzmax, double p,
1008 double alpha= -1E-16) ;
1009
1015 void exp_filter_r_all_domains(Scalar &ss, int p, double alpha=-16.) ;
1016
1023 void sarra_filter_r_all_domains(double p, double alpha=1E-16) ;
1024
1032 virtual void exponential_filter_ylm(int lzmin, int lzmax, int p,
1033 double alpha= -16.) ;
1034
1042 void annule_l (int l_min, int l_max, bool ylm_output= false ) ;
1043
1048 void filtre_phi (int n, int zone) ;
1049
1054 void filtre_tp(int nn, int nz1, int nz2) ;
1055
1056
1062 void fixe_decroissance (int puis) ;
1063
1069 void smooth_decay(int k, int n) ;
1070
1075 void raccord(int n) ;
1076
1083 void raccord_c1_zec(int puis, int nbre, int lmax) ;
1084
1088 void raccord_externe(int puis, int nbre, int lmax) ;
1089
1106 void match_tau(Param& par_bc, Param* par_mat=0x0) ;
1107
1116 void match_tau_shell(Param& par_bc, Param* par_mat=0x0) ;
1117
1126 void match_collocation(Param& par_bc, Param* par_mat=0x0) ;
1127
1128 // Outputs
1129 // -------
1130 public:
1131 virtual void sauve(FILE *) const ;
1132
1143 virtual void spectral_display(const char* comment = 0x0,
1144 double threshold = 1.e-7, int precision = 4,
1145 ostream& ostr = cout) const ;
1146
1148 friend ostream& operator<<(ostream& , const Scalar & ) ;
1149
1173 void visu_section(const char section_type, double aa, double umin, double umax, double vmin,
1174 double vmax, const char* title = 0x0, const char* filename = 0x0,
1175 bool start_dx = true, int nu = 200, int nv = 200) const ;
1176
1205 void visu_section(const Tbl& plane, double umin, double umax, double vmin,
1206 double vmax, const char* title = 0x0, const char* filename = 0x0,
1207 bool start_dx = true, int nu = 200, int nv = 200) const ;
1208
1237 void visu_section_anim(const char section_type, double aa, double umin,
1238 double umax, double vmin, double vmax, int jtime, double ttime,
1239 int jgraph = 1, const char* title = 0x0, const char* filename_root = 0x0,
1240 bool start_dx = false, int nu = 200, int nv = 200) const ;
1241
1263 void visu_box(double xmin, double xmax, double ymin, double ymax,
1264 double zmin, double zmax, const char* title0 = 0x0,
1265 const char* filename0 = 0x0, bool start_dx = true, int nx = 40, int ny = 40,
1266 int nz = 40) const ;
1267
1268
1269
1270 // Member arithmetics
1271 // ------------------
1272 public:
1273 void operator+=(const Scalar &) ;
1274 void operator-=(const Scalar &) ;
1275 void operator*=(const Scalar &) ;
1276
1277 // Manipulation of spectral bases
1278 // ------------------------------
1282 virtual void std_spectral_base() ;
1283
1287 virtual void std_spectral_base_odd() ;
1288
1291 void set_spectral_base(const Base_val& ) ;
1292
1294 const Base_val& get_spectral_base( ) const {return va.base ;} ;
1295
1301 void set_dzpuis(int ) ;
1302
1318 Valeur** asymptot(int n, const int flag = 0) const ;
1319
1320
1321 // PDE resolution
1322 // --------------
1323 public:
1333 Scalar poisson() const ;
1334
1346 void poisson(Param& par, Scalar& uu) const ;
1347
1357 Scalar poisson_tau() const ;
1358
1369 void poisson_tau(Param& par, Scalar& uu) const ;
1370
1385 Scalar poisson_dirichlet (const Valeur& limite, int num) const ;
1386
1391 Scalar poisson_neumann (const Valeur&, int) const ;
1392
1393
1412 Scalar poisson_dir_neu (const Valeur& limite , int num,
1413 double fact_dir, double fact_neu) const ;
1414
1421 Scalar poisson_frontiere_double (const Valeur&, const Valeur&, int) const ;
1422
1447 void poisson_regular(int k_div, int nzet, double unsgam1, Param& par,
1448 Scalar& uu, Scalar& uu_regu, Scalar& uu_div,
1449 Tensor& duu_div,
1450 Scalar& source_regu, Scalar& source_div) const ;
1451
1486 Tbl test_poisson(const Scalar& uu, ostream& ostr,
1487 bool detail = false) const ;
1488
1502 Scalar poisson_angu(double lambda =0) const ;
1503
1536 Scalar avance_dalembert(Param& par, const Scalar& fJm1, const Scalar& source)
1537 const ;
1538
1543 Scalar sol_elliptic(Param_elliptic& params) const ;
1544
1554 Scalar sol_elliptic_boundary(Param_elliptic& params, const Mtbl_cf& bound,
1555 double fact_dir, double fact_neu) const ;
1556
1561 Scalar sol_elliptic_boundary(Param_elliptic& params, const Scalar& bound,
1562 double fact_dir, double fact_neu) const ;
1563
1564
1572
1578
1586 Scalar sol_elliptic_no_zec(Param_elliptic& params, double val = 0) const ;
1587
1594 Scalar sol_elliptic_only_zec(Param_elliptic& params, double val) const ;
1595
1605 Scalar sol_elliptic_sin_zec(Param_elliptic& params, double* coefs, double* phases) const ;
1606
1607
1616 Param_elliptic& params) const ;
1617
1618
1629 Scalar sol_divergence(int n) const ;
1630
1631 // Import from other mapping
1632 // -------------------------
1633
1639 void import(const Scalar& ci) ;
1640
1647 void import_symy(const Scalar& ci) ;
1648
1656 void import_asymy(const Scalar& ci) ;
1657
1669 void import(int nzet, const Scalar& ci) ;
1670
1683 void import_symy(int nzet, const Scalar& ci) ;
1684
1698 void import_asymy(int nzet, const Scalar& ci) ;
1699
1700 protected:
1714 void import_gal(int nzet, const Scalar& ci) ;
1715
1729 void import_align(int nzet, const Scalar& ci) ;
1730
1745 void import_anti(int nzet, const Scalar& ci) ;
1746
1761 void import_align_symy(int nzet, const Scalar& ci) ;
1762
1778 void import_anti_symy(int nzet, const Scalar& ci) ;
1779
1795 void import_align_asymy(int nzet, const Scalar& ci) ;
1796
1813 void import_anti_asymy(int nzet, const Scalar& ci) ;
1814
1815
1816 friend Scalar operator-(const Scalar& ) ;
1817 friend Scalar operator+(const Scalar&, const Scalar &) ;
1818 friend Scalar operator+(const Scalar&, const Mtbl&) ;
1819 friend Scalar operator+(const Scalar&, double ) ;
1820 friend Scalar operator-(const Scalar &, const Scalar &) ;
1821 friend Scalar operator-(const Scalar&, const Mtbl&) ;
1822 friend Scalar operator-(const Scalar&, double ) ;
1823 friend Scalar operator*(const Scalar &, const Scalar &) ;
1824 friend Scalar operator%(const Scalar &, const Scalar &) ;
1825 friend Scalar operator|(const Scalar &, const Scalar &) ;
1826 friend Scalar operator*(const Mtbl&, const Scalar &) ;
1827 friend Scalar operator*(double, const Scalar &) ;
1828 friend Scalar operator/(const Scalar &, const Scalar &) ;
1829 friend Scalar operator/(const Scalar &, const Mtbl &) ;
1830 friend Scalar operator/(const Mtbl &, const Scalar &) ;
1831 friend Scalar operator/(const Scalar&, double ) ;
1832 friend Scalar operator/(double, const Scalar &) ;
1833
1834 friend Scalar sin(const Scalar& ) ;
1835 friend Scalar cos(const Scalar& ) ;
1836 friend Scalar tan(const Scalar& ) ;
1837 friend Scalar asin(const Scalar& ) ;
1838 friend Scalar acos(const Scalar& ) ;
1839 friend Scalar atan(const Scalar& ) ;
1840 friend Scalar exp(const Scalar& ) ;
1841 friend Scalar Heaviside(const Scalar& ) ;
1842 friend Scalar log(const Scalar& ) ;
1843 friend Scalar log10(const Scalar& ) ;
1844 friend Scalar sqrt(const Scalar& ) ;
1845 friend Scalar racine_cubique (const Scalar& ) ;
1846 friend Scalar pow(const Scalar& , int ) ;
1847 friend Scalar pow(const Scalar& , double ) ;
1848 friend Scalar abs(const Scalar& ) ;
1849
1850 friend double totalmax(const Scalar& ) ;
1851 friend double totalmin(const Scalar& ) ;
1852 friend Tbl max(const Scalar& ) ;
1853 friend Tbl min(const Scalar& ) ;
1854 friend Tbl norme(const Scalar& ) ;
1855 friend Tbl diffrel(const Scalar& a, const Scalar& b) ;
1856 friend Tbl diffrelmax(const Scalar& a, const Scalar& b) ;
1857
1858};
1859
1860ostream& operator<<(ostream& , const Scalar & ) ;
1861
1862// Prototypage de l'arithmetique
1869Scalar operator+(const Scalar& ) ;
1870Scalar operator-(const Scalar& ) ;
1871Scalar operator+(const Scalar&, const Scalar &) ;
1872Scalar operator+(const Scalar&, const Mtbl&) ;
1873Scalar operator+(const Mtbl&, const Scalar&) ;
1874Scalar operator+(const Scalar&, double ) ;
1875Scalar operator+(double, const Scalar& ) ;
1876Scalar operator+(const Scalar&, int ) ;
1877Scalar operator+(int, const Scalar& ) ;
1878Scalar operator-(const Scalar &, const Scalar &) ;
1879Scalar operator-(const Scalar&, const Mtbl&) ;
1880Scalar operator-(const Mtbl&, const Scalar&) ;
1881Scalar operator-(const Scalar&, double ) ;
1882Scalar operator-(double, const Scalar& ) ;
1883Scalar operator-(const Scalar&, int ) ;
1884Scalar operator-(int, const Scalar& ) ;
1885Scalar operator*(const Scalar &, const Scalar &) ;
1886
1888Scalar operator%(const Scalar &, const Scalar &) ;
1889
1891Scalar operator|(const Scalar &, const Scalar &) ;
1892
1893Scalar operator*(const Mtbl&, const Scalar&) ;
1894Scalar operator*(const Scalar&, const Mtbl&) ;
1895
1896Scalar operator*(const Scalar&, double ) ;
1897Scalar operator*(double, const Scalar &) ;
1898Scalar operator*(const Scalar&, int ) ;
1899Scalar operator*(int, const Scalar& ) ;
1900Scalar operator/(const Scalar &, const Scalar &) ;
1901Scalar operator/(const Scalar&, double ) ;
1902Scalar operator/(double, const Scalar &) ;
1903Scalar operator/(const Scalar&, int ) ;
1904Scalar operator/(int, const Scalar &) ;
1905Scalar operator/(const Scalar &, const Mtbl&) ;
1906Scalar operator/(const Mtbl&, const Scalar &) ;
1907
1908
1909Scalar sin(const Scalar& ) ;
1910Scalar cos(const Scalar& ) ;
1911Scalar tan(const Scalar& ) ;
1912Scalar asin(const Scalar& ) ;
1913Scalar acos(const Scalar& ) ;
1914Scalar atan(const Scalar& ) ;
1915Scalar exp(const Scalar& ) ;
1916Scalar Heaviside(const Scalar& ) ;
1917Scalar log(const Scalar& ) ;
1918Scalar log10(const Scalar& ) ;
1919Scalar sqrt(const Scalar& ) ;
1920Scalar racine_cubique (const Scalar& ) ;
1921Scalar pow(const Scalar& , int ) ;
1922Scalar pow(const Scalar& , double ) ;
1923Scalar abs(const Scalar& ) ;
1924
1930double totalmax(const Scalar& ) ;
1931
1937double totalmin(const Scalar& ) ;
1938
1944Tbl max(const Scalar& ) ;
1945
1951Tbl min(const Scalar& ) ;
1952
1959Tbl norme(const Scalar& ) ;
1960
1969Tbl diffrel(const Scalar& a, const Scalar& b) ;
1970
1979Tbl diffrelmax(const Scalar& a, const Scalar& b) ;
1980
1985void exp_filter_ylm_all_domains(Scalar& ss, int p, double alpha=-16.) ;
1986
1988}
1989#endif
Bases of the spectral expansions.
Definition base_val.h:322
Vectorial bases (triads) with respect to which the tensorial components are defined.
Definition base_vect.h:105
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition cmp.h:446
Base class for coordinate mappings.
Definition map.h:670
Metric for tensor calculation.
Definition metric.h:90
Multi-domain grid.
Definition grilles.h:273
Coefficients storage for the multi-domain spectral method.
Definition mtbl_cf.h:186
Multi-domain array.
Definition mtbl.h:118
This class contains the parameters needed to call the general elliptic solver.
Parameter storage.
Definition param.h:125
Tensor field of valence 0 (or component of a tensorial field).
Definition scalar.h:387
Scalar poisson_neumann(const Valeur &, int) const
Idem as Scalar::poisson_dirichlet , the boundary condition being on the radial derivative of the solu...
virtual void sauve(FILE *) const
Save in a file.
Definition scalar.C:686
void annule_l(int l_min, int l_max, bool ylm_output=false)
Sets all the multipolar components between l_min and l_max to zero.
friend Scalar exp(const Scalar &)
Exponential.
void div_cost()
Division by .
Scalar sol_divergence(int n) const
Resolution of a divergence-like equation.
const Scalar & lapang() const
Returns the angular Laplacian of *this , where .
void set_etat_one()
Sets the logical state to ETATUN (one).
Definition scalar.C:334
const Scalar & srdsdt() const
Returns of *this .
void raccord(int n)
Performs the matching of the nucleus with respect to the first shell.
virtual void del_deriv() const
Logical destructor of the derivatives.
Definition scalar.C:287
friend Scalar log(const Scalar &)
Neperian logarithm.
Scalar * p_lap
Pointer on the Laplacian of *this (0x0 if not up to date)
Definition scalar.h:448
friend Scalar cos(const Scalar &)
Cosine.
const Vector & derive_cov(const Metric &gam) const
Returns the gradient (1-form = covariant vector) of *this
void import_align(int nzet, const Scalar &ci)
Assignment to another Scalar defined on a different mapping, when the two mappings have aligned Carte...
Scalar * p_dsdy
Pointer on of *this , where (0x0 if not up to date)
Definition scalar.h:439
void sarra_filter_r(int lzmin, int lzmax, double p, double alpha=-1E-16)
Applies an exponential filter to the spectral coefficients in the radial direction.
Scalar sol_elliptic_pseudo_1d(Param_elliptic &) const
Solves a pseudo-1d elliptic equation with *this as a source.
Definition scalar_pde.C:430
void div_rsint_dzpuis(int ced_mult_r)
Division by but with the output flag dzpuis set to ced_mult_r .
const Scalar & deriv(int i) const
Returns of *this , where .
void div_r_dzpuis(int ced_mult_r)
Division by r everywhere but with the output flag dzpuis set to ced_mult_r .
const Scalar & dsdz() const
Returns of *this , where .
void import_asymy(const Scalar &ci)
Assignment to another Scalar defined on a different mapping.
Scalar scalar_out_bound(int n, bool leave_ylm=false)
Returns the Scalar containing the values of angular coefficients at the outer boundary.
void filtre(int n)
Sets the n lasts coefficients in r to 0 in the external domain.
const Scalar & dsdt() const
Returns of *this .
void div_sint()
Division by .
virtual ~Scalar()
Destructor.
Definition scalar.C:267
Scalar poisson() const
Solves the scalar Poisson equation with *this as a source.
Definition scalar_pde.C:136
friend Tbl norme(const Scalar &)
Sums of the absolute values of all the values of the Scalar in each domain.
int get_dzpuis() const
Returns dzpuis.
Definition scalar.h:557
const Scalar & laplacian(int ced_mult_r=4) const
Returns the Laplacian of *this.
const Scalar & dsdy() const
Returns of *this , where .
void smooth_decay(int k, int n)
Performs a matching of the last non-compactified shell with a decaying function where is the spher...
void import_anti_asymy(int nzet, const Scalar &ci)
Assignment to another Scalar defined on a different mapping, when the two mappings have anti-aligned ...
Scalar * p_dsdrho
Pointer on of *this
Definition scalar.h:458
void fixe_decroissance(int puis)
Substracts all the components behaving like in the external domain, with n strictly lower than puis ...
void operator=(const Scalar &a)
Assignment to another Scalar defined on the same mapping.
Definition scalar.C:446
Tbl tbl_in_bound(int n, bool leave_ylm=false)
Returns the Tbl containing the values of angular coefficients at the inner boundary.
virtual void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition scalar.C:353
double integrale() const
Computes the integral over all space of *this .
bool dz_nonzero() const
Returns true if the last domain is compactified and *this is not zero in this domain.
Definition scalar.C:814
void mult_cost()
Multiplication by .
friend Scalar operator%(const Scalar &, const Scalar &)
Scalar * Scalar with desaliasing.
Scalar poisson_angu(double lambda=0) const
Solves the (generalized) angular Poisson equation with *this as source.
Definition scalar_pde.C:200
friend Scalar asin(const Scalar &)
Arcsine.
const Scalar & srstdsdp() const
Returns of *this .
friend Scalar tan(const Scalar &)
Tangent.
friend double totalmax(const Scalar &)
Maximum values of a Scalar in each domain.
friend double totalmin(const Scalar &)
Minimum values of a Scalar in each domain.
void raccord_externe(int puis, int nbre, int lmax)
Matching of the external domain with the outermost shell.
Scalar * p_dsdr
Pointer on of *this (0x0 if not up to date)
Definition scalar.h:411
Scalar * p_stdsdp
Pointer on of *this (0x0 if not up to date)
Definition scalar.h:429
Scalar * p_dsdx
Pointer on of *this , where (0x0 if not up to date)
Definition scalar.h:434
virtual void std_spectral_base()
Sets the spectral bases of the Valeur va to the standard ones for a scalar field.
Definition scalar.C:784
double val_grid_point(int l, int k, int j, int i) const
Returns the value of the field at a specified grid point.
Definition scalar.h:637
const Scalar & dsdrho() const
Returns of *this .
Tbl tbl_out_bound(int l_dom, bool leave_ylm=false)
Returns the Tbl containing the values of angular coefficients at the outer boundary.
void import_gal(int nzet, const Scalar &ci)
Assignment to another Scalar defined on a different mapping, when the two mappings do not have a part...
friend Tbl min(const Scalar &)
Minimum values of a Scalar in each domain.
friend Scalar operator+(const Scalar &, const Scalar &)
Scalar + Scalar.
virtual void inc_dzpuis(int inc=1)
Increases by inc units the value of dzpuis and changes accordingly the values of the Scalar in the co...
friend Scalar Heaviside(const Scalar &)
Heaviside function.
Tbl test_poisson(const Scalar &uu, ostream &ostr, bool detail=false) const
Checks if a Poisson equation with *this as a source has been correctly solved.
bool is_nan(bool verbose=false) const
Looks for NaNs (not a number) in the scalar field.
Definition scalar.C:1008
void visu_box(double xmin, double xmax, double ymin, double ymax, double zmin, double zmax, const char *title0=0x0, const char *filename0=0x0, bool start_dx=true, int nx=40, int ny=40, int nz=40) const
3D visualization (volume rendering) via OpenDX.
const Scalar & stdsdp() const
Returns of *this .
friend Scalar operator*(const Scalar &, const Scalar &)
Scalar * Scalar.
void operator-=(const Scalar &)
-= Scalar
void match_tau_shell(Param &par_bc, Param *par_mat=0x0)
Method for matching accross domains and imposing boundary condition.
Scalar * p_dsdz
Pointer on of *this , where (0x0 if not up to date)
Definition scalar.h:444
friend ostream & operator<<(ostream &, const Scalar &)
Display.
Definition scalar.C:702
bool check_dzpuis(int dzi) const
Returns false if the last domain is compactified and *this is not zero in this domain and dzpuis is n...
Definition scalar.C:873
Scalar avance_dalembert(Param &par, const Scalar &fJm1, const Scalar &source) const
Performs one time-step integration (from ) of the scalar d'Alembert equation with *this being the val...
Definition scalar_pde.C:217
void div_r()
Division by r everywhere; dzpuis is not changed.
void exp_filter_r_all_domains(Scalar &ss, int p, double alpha=-16.)
Applies an exponential filter in radial direction in all domains.
void operator+=(const Scalar &)
+= Scalar
const Scalar & dsdx() const
Returns of *this , where .
virtual void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition scalar.C:324
void div_tant()
Division by .
const Scalar & dsdr() const
Returns of *this .
virtual void change_triad(const Base_vect &new_triad)
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.
Definition scalar.C:997
Tbl & set_domain(int l)
Read/write of the value in a given domain.
Definition scalar.h:615
const Scalar & dsdradial() const
Returns of *this if the mapping is affine (class Map_af) and of *this if the mapping is logarithmic...
int ind_lap
Power of r by which the last computed Laplacian has been multiplied in the compactified external doma...
Definition scalar.h:463
virtual void exponential_filter_ylm(int lzmin, int lzmax, int p, double alpha=-16.)
Applies an exponential filter to the spectral coefficients in the angular directions.
friend Scalar sin(const Scalar &)
Sine.
Definition scalar_math.C:72
virtual void annule(int l_min, int l_max)
Sets the Scalar to zero in several domains.
Definition scalar.C:391
Scalar * p_srstdsdp
Pointer on of *this (0x0 if not up to date)
Definition scalar.h:421
Valeur & set_spectral_va()
Returns va (read/write version)
Definition scalar.h:604
const Valeur & get_spectral_va() const
Returns va (read only version)
Definition scalar.h:601
void visu_section_anim(const char section_type, double aa, double umin, double umax, double vmin, double vmax, int jtime, double ttime, int jgraph=1, const char *title=0x0, const char *filename_root=0x0, bool start_dx=false, int nu=200, int nv=200) const
3D visualization via time evolving plane section (animation).
void annule_hard()
Sets the Scalar to zero in a hard way.
Definition scalar.C:380
void mult_sint()
Multiplication by .
void del_t()
Logical destructor.
Definition scalar.C:279
Scalar poisson_dir_neu(const Valeur &limite, int num, double fact_dir, double fact_neu) const
Is identicall to Scalar::poisson() .
void mult_r_ced()
Multiplication by r in the compactified external domain (CED), the dzpuis flag is not changed.
friend Scalar operator/(const Scalar &, const Scalar &)
Scalar / Scalar.
void filtre_tp(int nn, int nz1, int nz2)
Sets the n lasts coefficients in to 0 in the domain nz1 to nz2 when expressed in spherical harmonics...
int etat
The logical state ETATNONDEF (undefined), ETATZERO (null), ETATUN (one), or ETATQCQ (ordinary).
Definition scalar.h:396
int get_etat() const
Returns the logical state ETATNONDEF (undefined), ETATZERO (null) or ETATQCQ (ordinary).
Definition scalar.h:554
virtual void allocate_all()
Sets the logical state to ETATQCQ (ordinary state) and performs the memory allocation of all the elem...
Definition scalar.C:365
virtual void set_etat_nondef()
Sets the logical state to ETATNONDEF (undefined).
Definition scalar.C:344
void visu_section(const char section_type, double aa, double umin, double umax, double vmin, double vmax, const char *title=0x0, const char *filename=0x0, bool start_dx=true, int nu=200, int nv=200) const
3D visualization via a plane section.
Definition scalar_visu.C:78
const Base_val & get_spectral_base() const
Returns the spectral bases of the Valeur va
Definition scalar.h:1294
void mult_rsint()
Multiplication by everywhere; dzpuis is not changed.
void set_dzpuis(int)
Modifies the dzpuis flag.
Definition scalar.C:808
Scalar sol_elliptic_boundary(Param_elliptic &params, const Mtbl_cf &bound, double fact_dir, double fact_neu) const
Resolution of a general elliptic equation, putting zero at infinity and with inner boundary condition...
Definition scalar_pde.C:256
friend Tbl diffrelmax(const Scalar &a, const Scalar &b)
Relative difference between two Scalar (max version).
Scalar poisson_frontiere_double(const Valeur &, const Valeur &, int) const
Idem as Scalar::poisson_dirichlet , the boundary condition being on both the function and its radial ...
double & set_grid_point(int l, int k, int j, int i)
Setting the value of the field at a given grid point.
Definition scalar.h:684
void set_inner_boundary(int l, double x)
Sets the value of the Scalar at the inner boundary of a given domain.
void raccord_c1_zec(int puis, int nbre, int lmax)
Performs the matching of the external domain with respect to the last shell using function like wit...
void filtre_r(int *nn)
Sets the n lasts coefficients in r to 0 in all domains.
friend Scalar racine_cubique(const Scalar &)
Cube root.
void match_tau(Param &par_bc, Param *par_mat=0x0)
Method for matching accross domains and imposing boundary condition.
void import_align_symy(int nzet, const Scalar &ci)
Assignment to another Scalar defined on a different mapping, when the two mappings have aligned Carte...
friend Scalar abs(const Scalar &)
Absolute value.
Scalar * p_lapang
Pointer on the Laplacian of *this (0x0 if not up to date)
Definition scalar.h:452
void div_rsint()
Division by everywhere; dzpuis is not changed.
Scalar sol_elliptic(Param_elliptic &params) const
Resolution of a general elliptic equation, putting zero at infinity.
Definition scalar_pde.C:234
friend Scalar sqrt(const Scalar &)
Square root.
void mult_rsint_dzpuis(int ced_mult_r)
Multiplication by but with the output flag dzpuis set to ced_mult_r .
const Tbl & domain(int l) const
Read-only of the value in a given domain.
Definition scalar.h:625
friend Scalar log10(const Scalar &)
Basis 10 logarithm.
double val_point(double r, double theta, double phi) const
Computes the value of the field at an arbitrary point , by means of the spectral expansion.
Definition scalar.C:890
void mult_r_dzpuis(int ced_mult_r)
Multiplication by r everywhere but with the output flag dzpuis set to ced_mult_r .
Scalar sol_elliptic_only_zec(Param_elliptic &params, double val) const
Resolution of a general elliptic equation solving in the compactified domain and putting a given valu...
Definition scalar_pde.C:341
friend Scalar pow(const Scalar &, int)
Power .
void match_collocation(Param &par_bc, Param *par_mat=0x0)
Method for matching accross domains and imposing boundary condition.
void operator*=(const Scalar &)
*= Scalar
void filtre_phi(int n, int zone)
Sets the n lasts coefficients in to 0 in the domain zone .
Tbl multipole_spectrum() const
Gives the spectrum in terms of multipolar modes l .
Definition scalar.C:956
void sarra_filter_r_all_domains(double p, double alpha=1E-16)
Applies an exponential filter in radial direction in all domains for the case where p is a double (se...
virtual void std_spectral_base_odd()
Sets the spectral bases of the Valeur va to the standard odd ones for a scalar field.
Definition scalar.C:791
void set_outer_boundary(int l, double x)
Sets the value of the Scalar at the outer boundary of a given domain.
const Vector & derive_con(const Metric &gam) const
Returns the "contravariant" derivative of *this with respect to some metric , by raising the index of...
Scalar sol_elliptic_no_zec(Param_elliptic &params, double val=0) const
Resolution of a general elliptic equation, putting a given value at the outermost shell and not solvi...
Definition scalar_pde.C:314
Valeur va
The numerical value of the Scalar
Definition scalar.h:405
virtual void spectral_display(const char *comment=0x0, double threshold=1.e-7, int precision=4, ostream &ostr=cout) const
Displays the spectral coefficients and the associated basis functions.
Definition scalar.C:741
void poisson_regular(int k_div, int nzet, double unsgam1, Param &par, Scalar &uu, Scalar &uu_regu, Scalar &uu_div, Tensor &duu_div, Scalar &source_regu, Scalar &source_div) const
Solves the scalar Poisson equation with *this as a source (version with parameters to control the res...
void import_anti(int nzet, const Scalar &ci)
Assignment to another Scalar defined on a different mapping, when the two mappings have anti-aligned ...
Scalar * p_dsdradial
Pointer on of *this
Definition scalar.h:455
Tbl * p_integ
Pointer on the space integral of *this (values in each domain) (0x0 if not up to date)
Definition scalar.h:468
Scalar * p_dsdt
Pointer on of *this (0x0 if not up to date)
Definition scalar.h:424
Scalar poisson_dirichlet(const Valeur &limite, int num) const
Is identicall to Scalar::poisson() .
friend Tbl diffrel(const Scalar &a, const Scalar &b)
Relative difference between two Scalar (norme version).
const Tbl & integrale_domains() const
Computes the integral in each domain of *this .
friend Tbl max(const Scalar &)
Maximum values of a Scalar in each domain.
Scalar sol_elliptic_2d(Param_elliptic &) const
Solves the scalar 2-dimensional elliptic equation with *this as a source.
Definition scalar_pde.C:409
Scalar sol_elliptic_fixe_der_zero(double val, Param_elliptic &params) const
Resolution of a general elliptic equation fixing the dericative at the origin and relaxing one contin...
Definition scalar_pde.C:386
Scalar poisson_tau() const
Solves the scalar Poisson equation with *this as a source using a real Tau method The source of the ...
Definition scalar_pde.C:169
void set_der_0x0() const
Sets the pointers for derivatives to 0x0.
Definition scalar.C:306
virtual void dec_dzpuis(int dec=1)
Decreases by dec units the value of dzpuis and changes accordingly the values of the Scalar in the co...
void import_anti_symy(int nzet, const Scalar &ci)
Assignment to another Scalar defined on a different mapping, when the two mappings have anti-aligned ...
Scalar primr(bool null_infty=true) const
Computes the radial primitive which vanishes for .
Scalar sol_elliptic_sin_zec(Param_elliptic &params, double *coefs, double *phases) const
General elliptic solver.
Definition scalar_pde.C:363
friend Scalar acos(const Scalar &)
Arccosine.
void import_align_asymy(int nzet, const Scalar &ci)
Assignment to another Scalar defined on a different mapping, when the two mappings have aligned Carte...
void mult_r()
Multiplication by r everywhere; dzpuis is not changed.
Valeur ** asymptot(int n, const int flag=0) const
Asymptotic expansion at r = infinity.
Scalar * p_srdsdt
Pointer on of *this (0x0 if not up to date)
Definition scalar.h:416
friend Scalar operator-(const Scalar &)
- Scalar
Scalar derive_lie(const Vector &v) const
Computes the derivative of this along a vector field v.
void div_r_ced()
Division by r in the compactified external domain (CED), the dzpuis flag is not changed.
void import_symy(const Scalar &ci)
Assignment to another Scalar defined on a different mapping.
void set_spectral_base(const Base_val &)
Sets the spectral bases of the Valeur va
Definition scalar.C:797
int dzpuis
Power of r by which the quantity represented by this must be divided in the compactified external d...
Definition scalar.h:403
friend Scalar operator|(const Scalar &, const Scalar &)
Scalar * Scalar with desaliasing only in r.
virtual void exponential_filter_r(int lzmin, int lzmax, int p, double alpha=-16.)
Applies an exponential filter to the spectral coefficients in the radial direction.
friend Scalar atan(const Scalar &)
Arctangent.
Basic array class.
Definition tbl.h:161
Tensor handling.
Definition tensor.h:288
Values and coefficients of a (real-value) function.
Definition valeur.h:287
Tbl & set(int l)
Read/write of the value in a given domain (configuration space).
Definition valeur.h:363
Base_val base
Bases on which the spectral expansion is performed.
Definition valeur.h:305
Tensor field of valence 1.
Definition vector.h:188
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 atan(const Cmp &)
Arctangent.
Definition cmp_math.C:195
Cmp sqrt(const Cmp &)
Square root.
Definition cmp_math.C:220
Cmp log10(const Cmp &)
Basis 10 logarithm.
Definition cmp_math.C:322
Cmp operator/(const Cmp &, const Cmp &)
Cmp / Cmp.
Definition cmp_arithm.C:457
Cmp operator%(const Cmp &, const Cmp &)
Cmp * Cmp with desaliasing.
Definition cmp_arithm.C:364
Cmp exp(const Cmp &)
Exponential.
Definition cmp_math.C:270
Cmp sin(const Cmp &)
Sine.
Definition cmp_math.C:69
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
Cmp acos(const Cmp &)
Arccosine.
Definition cmp_math.C:169
Cmp asin(const Cmp &)
Arcsine.
Definition cmp_math.C:144
Cmp racine_cubique(const Cmp &)
Cube root.
Definition cmp_math.C:245
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 pow(const Cmp &, int)
Power .
Definition cmp_math.C:348
Cmp cos(const Cmp &)
Cosine.
Definition cmp_math.C:94
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
Cmp tan(const Cmp &)
Tangent.
Definition cmp_math.C:120
Cmp log(const Cmp &)
Neperian logarithm.
Definition cmp_math.C:296
Mtbl Heaviside(const Mtbl &)
Heaviside function.
Definition mtbl_math.C:317
double totalmin(const Mtbl &)
Minimum value of the Mtbl elements in all domain.
Definition mtbl_math.C:522
double totalmax(const Mtbl &)
Maximum value of the Mtbl elements in all domains.
Definition mtbl_math.C:494
Scalar operator|(const Scalar &, const Scalar &)
Scalar * Scalar with desaliasing only in r.
void exp_filter_ylm_all_domains(Scalar &ss, int p, double alpha=-16.)
Applies an exponential filter in angular directions in all domains.
Lorene prototypes.
Definition app_hor.h:64