KADATH
spheric_symphi.hpp
1 /*
2  Copyright 2017 Philippe Grandclement
3 
4  This file is part of Kadath.
5 
6  Kadath 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 3 of the License, or
9  (at your option) any later version.
10 
11  Kadath 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 Kadath. If not, see <http://www.gnu.org/licenses/>.
18 */
19 
20 #ifndef __SPHERIC_SYMPHI_HPP_
21 #define __SPHERIC_SYMPHI_HPP_
22 
23 #include "space.hpp"
24 #include "term_eq.hpp"
25 namespace Kadath {
55 class Domain_nucleus_symphi : public Domain {
56 
57  private:
58  double alpha ;
60 
61  public:
70  Domain_nucleus_symphi (int num, int ttype, double radius, const Point& cr, const Dim_array& nbr) ;
80  Domain_nucleus_symphi (int num, FILE* fd) ;
81 
82  virtual ~Domain_nucleus_symphi() ;
83  virtual void save (FILE*) const ;
84 
85  private:
86  virtual void do_absol () const ;
87  virtual void do_radius () const ;
88  virtual void do_cart () const ;
89  virtual void do_cart_surr () const ;
90 
91  virtual void set_cheb_base(Base_spectral&) const;
92  virtual void set_legendre_base(Base_spectral&) const ;
93 
94  virtual void set_cheb_base_r_spher(Base_spectral&) const ;
95  virtual void set_cheb_base_t_spher(Base_spectral&) const ;
96  virtual void set_cheb_base_p_spher(Base_spectral&) const ;
97  virtual void set_cheb_base_rt_spher(Base_spectral&) const ;
98  virtual void set_cheb_base_rp_spher(Base_spectral&) const ;
99  virtual void set_cheb_base_tp_spher(Base_spectral&) const ;
100  virtual void set_legendre_base_r_spher(Base_spectral&) const ;
101  virtual void set_legendre_base_t_spher(Base_spectral&) const ;
102  virtual void set_legendre_base_p_spher(Base_spectral&) const ;
103 
104  virtual void set_cheb_base_x_cart(Base_spectral&) const ;
105  virtual void set_cheb_base_y_cart(Base_spectral&) const ;
106  virtual void set_cheb_base_z_cart(Base_spectral&) const ;
107  virtual void set_cheb_base_xy_cart(Base_spectral&) const ;
108  virtual void set_cheb_base_xz_cart(Base_spectral&) const ;
109  virtual void set_cheb_base_yz_cart(Base_spectral&) const ;
110  virtual void set_legendre_base_x_cart(Base_spectral&) const ;
111  virtual void set_legendre_base_y_cart(Base_spectral&) const ;
112  virtual void set_legendre_base_z_cart(Base_spectral&) const ;
113 
119  void set_cheb_base_forr(Base_spectral& so) const ;
125  void set_cheb_base_forx_cart(Base_spectral& so) const ;
131  void set_cheb_base_fory_cart(Base_spectral& so) const ;
137  void set_cheb_base_forz_cart(Base_spectral& so) const ;
143  void set_legendre_base_forr(Base_spectral& so) const ;
149  void set_legendre_base_forx_cart(Base_spectral& so ) const ;
155  void set_legendre_base_fory_cart(Base_spectral& so ) const ;
162 
163 
164  virtual void do_coloc () ;
165  virtual int give_place_var (char*) const ;
166 
167  public:
168  virtual double get_rmax() const {return alpha ; } ;
169  virtual Point get_center () const {return center ;} ;
170 
177  virtual bool is_in(const Point&xx, double prec=1e-13) const ;
186  virtual const Point absol_to_num(const Point& xxx) const;
187 
188  virtual const Point absol_to_num_bound(const Point&, int) const;
189 
207  virtual void do_der_abs_from_der_var(const Val_domain *const *const der_var, Val_domain **const der_abs) const ;
213  virtual Base_spectral mult (const Base_spectral&, const Base_spectral&) const ;
214 
221  int mult_base_angle_int(int basea, int baseb) const;
228  int mult_base_r_int(int basea, int baseb) const;
229 
230  public:
231 
232  virtual Val_domain mult_cos_phi (const Val_domain&) const ;
233  virtual Val_domain mult_sin_phi (const Val_domain&) const ;
234  virtual Val_domain mult_cos_theta (const Val_domain&) const ;
235  virtual Val_domain mult_sin_theta (const Val_domain&) const ;
236  virtual Val_domain div_sin_theta (const Val_domain&) const ;
237  virtual Val_domain div_cos_theta (const Val_domain&) const ;
238 
239  virtual Tensor change_basis_cart_to_spher (int dd, const Tensor&) const ;
240  virtual Tensor change_basis_spher_to_cart (int dd, const Tensor&) const ;
241 
242  virtual Val_domain div_x (const Val_domain&) const ;
243  virtual Val_domain mult_r (const Val_domain&) const ;
244  virtual Val_domain div_r (const Val_domain&) const ;
245  virtual Val_domain der_r (const Val_domain&) const ;
246  virtual Val_domain ddp (const Val_domain&) const ;
247  virtual Val_domain srdr (const Val_domain&) const ;
248  virtual Val_domain div_1mx2 (const Val_domain&) const ;
249 
250  virtual double val_boundary (int, const Val_domain&, const Index&) const ;
251  virtual void find_other_dom (int, int, int&, int&) const ;
252  virtual Val_domain der_normal (const Val_domain&, int) const ;
253 
254  virtual int nbr_unknowns (const Tensor&, int) const ;
261  int nbr_unknowns_val_domain (const Val_domain& so) const ;
262 
263  virtual Array<int> nbr_conditions (const Tensor&, int, int, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
271  int nbr_conditions_val_domain (const Val_domain& eq, int order) const ;
272 
273  virtual Array<int> nbr_conditions_boundary (const Tensor&, int, int, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
281  int nbr_conditions_val_domain_boundary (const Val_domain& eq) const ;
282 
283  virtual void export_tau (const Tensor&, int, int, Array<double>&, int&, const Array<int>&, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
293  void export_tau_val_domain (const Val_domain& eq, int order, Array<double>& res, int& pos_res, int ncond) const ;
294 
295  virtual void export_tau_boundary (const Tensor&, int, int, Array<double>&, int&, const Array<int>&, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
305  void export_tau_val_domain_boundary (const Val_domain& eq, int bound, Array<double>& res, int& pos_res, int ncond) const ;
306 
307  virtual void affecte_tau (Tensor&, int, const Array<double>&, int&) const ;
315  void affecte_tau_val_domain (Val_domain& so, const Array<double>& cf, int& pos_cf) const ;
316 
317  virtual void affecte_tau_one_coef (Tensor&, int, int, int&) const ;
325  void affecte_tau_one_coef_val_domain (Val_domain& so, int cc, int& pos_cf) const ;
326 
327  virtual double integ_volume (const Val_domain& so) const ;
328 
329  virtual Term_eq derive_flat_spher (int, char, const Term_eq&, const Metric*) const ;
330  virtual Term_eq derive_flat_cart (int, char, const Term_eq&, const Metric*) const ;
331  virtual double integ(const Val_domain& so, int bound) const ;
332 
339  void filter_val_domain (Val_domain& so, double threshold) const ;
340  virtual void filter (Tensor& tt, int dom, double threshold) const ;
341 
342 public:
343  virtual ostream& print (ostream& o) const ;
344 } ;
345 
375 class Domain_shell_symphi : public Domain {
376 
377  private:
378  double alpha ;
379  double beta ;
381 
382  public:
392  Domain_shell_symphi (int num, int ttype, double r_int, double r_ext, const Point& cr, const Dim_array& nbr) ;
399  Domain_shell_symphi (int num, FILE* fd) ;
400 
401  virtual ~Domain_shell_symphi() ;
402  virtual void save (FILE*) const ;
403 
404  private:
405  virtual void do_absol () const ;
406  virtual void do_radius () const ;
407  virtual void do_cart () const ;
408  virtual void do_cart_surr () const ;
409 
410 
411  private:
412  virtual void set_cheb_base(Base_spectral&) const ;
413  virtual void set_legendre_base(Base_spectral&) const ;
414 
415  virtual void set_cheb_base_r_spher(Base_spectral&) const ;
416  virtual void set_cheb_base_t_spher(Base_spectral&) const ;
417  virtual void set_cheb_base_p_spher(Base_spectral&) const ;
418  virtual void set_cheb_base_rt_spher(Base_spectral&) const ;
419  virtual void set_cheb_base_rp_spher(Base_spectral&) const ;
420  virtual void set_cheb_base_tp_spher(Base_spectral&) const ;
421  virtual void set_legendre_base_r_spher(Base_spectral&) const ;
422  virtual void set_legendre_base_t_spher(Base_spectral&) const ;
423  virtual void set_legendre_base_p_spher(Base_spectral&) const ;
424 
425  virtual void set_cheb_base_x_cart(Base_spectral&) const ;
426  virtual void set_cheb_base_y_cart(Base_spectral&) const ;
427  virtual void set_cheb_base_z_cart(Base_spectral&) const ;
428  virtual void set_cheb_base_xy_cart(Base_spectral&) const ;
429  virtual void set_cheb_base_xz_cart(Base_spectral&) const ;
430  virtual void set_cheb_base_yz_cart(Base_spectral&) const ;
431  virtual void set_legendre_base_x_cart(Base_spectral&) const ;
432  virtual void set_legendre_base_y_cart(Base_spectral&) const ;
433  virtual void set_legendre_base_z_cart(Base_spectral&) const ;
439  void set_cheb_base_forr(Base_spectral& so) const ;
445  void set_cheb_base_forx_cart(Base_spectral& so) const ;
451  void set_cheb_base_fory_cart(Base_spectral& so) const ;
457  void set_cheb_base_forz_cart(Base_spectral& so) const ;
463  void set_legendre_base_forr(Base_spectral& so) const ;
482 
483 
484  virtual void do_coloc () ;
485  virtual int give_place_var (char*) const ;
486  public:
487  virtual double get_rmax() const {return alpha+beta ; } ;
488  virtual Point get_center () const {return center ;} ;
489 
496  virtual bool is_in(const Point&xx, double prec=1e-13) const ;
505  virtual const Point absol_to_num(const Point& xxx) const;
506 
507  virtual const Point absol_to_num_bound(const Point&, int) const;
508 
526  virtual void do_der_abs_from_der_var(const Val_domain *const *const der_var, Val_domain **const der_abs) const ;
532  virtual Base_spectral mult (const Base_spectral&, const Base_spectral&) const ;
539  int mult_base_angle_int(int basea, int baseb) const;
540 
541  public:
542  virtual Val_domain mult_cos_phi (const Val_domain&) const ;
543  virtual Val_domain mult_sin_phi (const Val_domain&) const ;
544  virtual Val_domain mult_cos_theta (const Val_domain&) const ;
545  virtual Val_domain mult_sin_theta (const Val_domain&) const ;
546  virtual Val_domain div_sin_theta (const Val_domain&) const ;
547  virtual Val_domain div_cos_theta (const Val_domain&) const ;
548 
549  virtual Tensor change_basis_cart_to_spher (int dd, const Tensor&) const ;
550  virtual Tensor change_basis_spher_to_cart (int dd, const Tensor&) const ;
551 
552 
553  virtual Val_domain mult_r (const Val_domain&) const ;
554  virtual Val_domain div_r (const Val_domain&) const ;
555  virtual Val_domain der_r (const Val_domain&) const ;
556  virtual Val_domain der_partial_var (const Val_domain&, int) const ;
557  virtual Val_domain ddp (const Val_domain&) const ;
558  virtual Val_domain div_xm1 (const Val_domain&) const ;
559  virtual Val_domain mult_xm1 (const Val_domain&) const ;
560  virtual Val_domain div_1mx2 (const Val_domain&) const ;
561  virtual Val_domain div_1mrsL (const Val_domain&) const ;
562  virtual Val_domain mult_1mrsL (const Val_domain&) const ;
563 
564 
565  virtual double val_boundary (int, const Val_domain&, const Index&) const ;
566  virtual void find_other_dom (int, int, int&, int&) const ;
567  virtual Val_domain der_normal (const Val_domain&, int) const ;
568 
569  virtual int nbr_unknowns (const Tensor&, int) const ;
576  int nbr_unknowns_val_domain (const Val_domain& so) const ;
577 
578  virtual Array<int> nbr_conditions (const Tensor&, int, int, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
586  int nbr_conditions_val_domain (const Val_domain& eq, int order) const ;
587 
588  virtual Array<int> nbr_conditions_boundary (const Tensor&, int, int, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
596  int nbr_conditions_val_domain_boundary (const Val_domain& eq) const ;
597 
598  virtual void export_tau (const Tensor&, int, int, Array<double>&, int&, const Array<int>&, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
608  void export_tau_val_domain (const Val_domain& eq, int order, Array<double>& res, int& pos_res, int ncond) const ;
609 
610  virtual void export_tau_boundary (const Tensor&, int, int, Array<double>&, int&, const Array<int>&, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
620  void export_tau_val_domain_boundary (const Val_domain& eq, int bound, Array<double>& res, int& pos_res, int ncond) const ;
621 
622  virtual void affecte_tau (Tensor&, int, const Array<double>&, int&) const ;
630  void affecte_tau_val_domain (Val_domain& so, const Array<double>& cf, int& pos_cf) const ;
631 
632  virtual void affecte_tau_one_coef (Tensor&, int, int, int&) const ;
640  void affecte_tau_one_coef_val_domain (Val_domain& so, int cc, int& pos_cf) const ;
641 
642  virtual void export_tau_boundary_exception (const Tensor&, int, int, Array<double>&, int&, const Array<int>&, const Param&, int,
643  const Tensor&, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
644 
658  void export_tau_val_domain_boundary_exception (const Val_domain& eq, int mlim, int bound, Array<double>& res, int& pos_res, int ncond, const Param& param, int type_exception, const Val_domain& exception) const ;
659 
660  virtual double integ_volume (const Val_domain& so) const ;
661 
662  virtual Term_eq derive_flat_spher (int, char, const Term_eq&, const Metric*) const ;
663  virtual Term_eq derive_flat_cart (int, char, const Term_eq&, const Metric*) const ;
664  virtual double integ(const Val_domain& so, int bound) const ;
665 
666 public:
667  virtual ostream& print (ostream& o) const ;
668 } ;
669 
670 
701 
702  private:
703  double alpha ;
705 
706  public:
715  Domain_compact_symphi (int num, int ttype, double r_int, const Point& cr, const Dim_array& nbr) ;
722  Domain_compact_symphi (int num, FILE* fd) ;
723 
724  virtual ~Domain_compact_symphi () ;
725  virtual void save(FILE*) const ;
726 
727  private:
728  virtual void do_absol () const ;
729  virtual void do_radius () const ;
730  virtual void do_cart () const ;
731  virtual void do_cart_surr () const ;
732 
733  private:
734 
735  virtual void set_cheb_base(Base_spectral&) const ;
736  virtual void set_legendre_base(Base_spectral&) const ;
737 
738  virtual void set_cheb_base_r_spher(Base_spectral&) const ;
739  virtual void set_cheb_base_t_spher(Base_spectral&) const ;
740  virtual void set_cheb_base_p_spher(Base_spectral&) const ;
741  virtual void set_legendre_base_r_spher(Base_spectral&) const ;
742  virtual void set_legendre_base_t_spher(Base_spectral&) const ;
743  virtual void set_legendre_base_p_spher(Base_spectral&) const ;
744 
745  virtual void set_cheb_base_x_cart(Base_spectral&) const ;
746  virtual void set_cheb_base_y_cart(Base_spectral&) const ;
747  virtual void set_cheb_base_z_cart(Base_spectral&) const ;
748  virtual void set_legendre_base_x_cart(Base_spectral&) const ;
749  virtual void set_legendre_base_y_cart(Base_spectral&) const ;
750  virtual void set_legendre_base_z_cart(Base_spectral&) const ;
756  void set_cheb_base_forx_cart(Base_spectral& so) const ;
762  void set_cheb_base_fory_cart(Base_spectral& so) const ;
768  void set_cheb_base_forz_cart(Base_spectral& so) const ;
787 
788 
789  virtual void do_coloc () ;
790  virtual int give_place_var (char*) const ;
791  public:
792  virtual double get_rmax() const {return alpha ; } ;
793  virtual Point get_center () const {return center ;} ;
794 
801  virtual bool is_in(const Point&xx, double prec=1e-13) const ;
810  virtual const Point absol_to_num(const Point& xxx) const;
811 
812  virtual const Point absol_to_num_bound(const Point&, int) const;
813 
831  virtual void do_der_abs_from_der_var(const Val_domain *const *const der_var, Val_domain **const der_abs) const ;
837  virtual Base_spectral mult (const Base_spectral&, const Base_spectral&) const ;
838 
839  public:
840 
841  virtual Tensor change_basis_cart_to_spher (int dd, const Tensor&) const ;
842  virtual Tensor change_basis_spher_to_cart (int dd, const Tensor&) const ;
843 
844 
845 
846  virtual Val_domain mult_cos_phi (const Val_domain&) const ;
847  virtual Val_domain mult_sin_phi (const Val_domain&) const ;
848  virtual Val_domain mult_cos_theta (const Val_domain&) const ;
849  virtual Val_domain mult_sin_theta (const Val_domain&) const ;
850  virtual Val_domain div_sin_theta (const Val_domain&) const ;
851  virtual Val_domain div_cos_theta (const Val_domain&) const ;
852  virtual Val_domain div_xm1 (const Val_domain&) const ;
853  virtual Val_domain mult_xm1 (const Val_domain&) const ;
854  virtual Val_domain mult_r (const Val_domain&) const ;
855  virtual Val_domain div_r (const Val_domain&) const ;
856  virtual Val_domain der_r (const Val_domain&) const ;
857  virtual Val_domain der_r_rtwo (const Val_domain&) const ;
858  virtual Val_domain der_partial_var (const Val_domain&, int) const ;
859 
860  virtual void set_val_inf (Val_domain& so, double xx) const ;
861 
862  virtual double val_boundary (int, const Val_domain&, const Index&) const ;
863  virtual void find_other_dom (int, int, int&, int&) const ;
864  virtual Val_domain der_normal (const Val_domain&, int) const ;
865 
866  virtual int nbr_unknowns (const Tensor&, int) const ;
873  int nbr_unknowns_val_domain (const Val_domain& so) const ;
874 
875  virtual Array<int> nbr_conditions (const Tensor&, int, int, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
883  int nbr_conditions_val_domain (const Val_domain& eq, int order) const ;
884 
885  virtual Array<int> nbr_conditions_boundary (const Tensor&, int, int, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
893  int nbr_conditions_val_domain_boundary (const Val_domain& eq) const ;
894 
895  virtual void export_tau (const Tensor&, int, int, Array<double>&, int&, const Array<int>&, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
905  void export_tau_val_domain (const Val_domain& eq, int order, Array<double>& res, int& pos_res, int ncond) const ;
906 
907  virtual void export_tau_boundary (const Tensor&, int, int, Array<double>&, int&, const Array<int>&, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
917  void export_tau_val_domain_boundary (const Val_domain& eq, int bound, Array<double>& res, int& pos_res, int ncond) const ;
918 
919  virtual void affecte_tau (Tensor&, int, const Array<double>&, int&) const ;
927  void affecte_tau_val_domain (Val_domain& so, const Array<double>& cf, int& pos_cf) const ;
928 
929  virtual void affecte_tau_one_coef (Tensor&, int, int, int&) const ;
937  void affecte_tau_one_coef_val_domain (Val_domain& so, int cc, int& pos_cf) const ;
938 
939  virtual double integ_volume (const Val_domain& so) const ;
940 
941  virtual Term_eq derive_flat_spher (int, char, const Term_eq&, const Metric*) const ;
942  virtual Term_eq derive_flat_cart (int, char, const Term_eq&, const Metric*) const ;
943  virtual double integ(const Val_domain& so, int bound) const ;
944 
945 public:
946  virtual ostream& print (ostream& o) const ;
947 } ;
948 
949 
954 class Space_spheric_symphi : public Space {
955  public:
963  Space_spheric_symphi (int ttype, const Point& cr, const Dim_array& nbr, double bound) ;
964 
973  Space_spheric_symphi (int ttype, const Point& cr, const Dim_array& nbr, const Array<double>& bounds, bool withzec=false) ;
974 
975 
976  Space_spheric_symphi (FILE* fd, bool withzec=false) ;
977  virtual ~Space_spheric_symphi() ;
978  virtual void save(FILE*) const ;
979 
985  void add_eq_ori (System_of_eqs& syst, const char* eq) ;
986 } ;
987 
988 }
989 #endif
Class for storing the basis of decompositions of a field.
Class for storing the dimensions of an array.
Definition: dim_array.hpp:34
Class for a spherical compactified domain and a symmetry with respect to the plane and a quadrant sy...
virtual Val_domain mult_cos_phi(const Val_domain &) const
Multiplication by .
virtual void set_cheb_base(Base_spectral &) const
Gives the standard base for Chebyshev polynomials.
double alpha
Relates the numerical to the physical radii.
virtual Val_domain div_xm1(const Val_domain &) const
Division by .
virtual void set_legendre_base(Base_spectral &) const
Gives the standard base for Legendre polynomials.
virtual double get_rmax() const
Returns the maximum radius.
virtual Val_domain mult_sin_theta(const Val_domain &) const
Multiplication by .
virtual double integ(const Val_domain &so, int bound) const
Surface integral on a given boundary.
void set_legendre_base_forx_cart(Base_spectral &so) const
Sets the base to the standard one for Legendre polynomials for a field like the component of a vecto...
virtual Tensor change_basis_spher_to_cart(int dd, const Tensor &) const
Changes the tensorial basis from spherical to Cartesian in a given domain.
virtual const Point absol_to_num(const Point &xxx) const
Computes the numerical coordinates from the physical ones.
virtual Base_spectral mult(const Base_spectral &, const Base_spectral &) const
Method for the multiplication of two Base_spectral.
virtual ostream & print(ostream &o) const
Delegate function to virtualize the << operator.
virtual void affecte_tau_one_coef(Tensor &, int, int, int &) const
Sets at most one coefficient of a Tensor to 1.
virtual Val_domain div_cos_theta(const Val_domain &) const
Division by .
virtual Val_domain der_partial_var(const Val_domain &, int) const
Partial derivative with respect to a coordinate.
virtual Val_domain mult_cos_theta(const Val_domain &) const
Multiplication by .
virtual void do_radius() const
Computes the generalized radius.
virtual void do_cart_surr() const
Computes the Cartesian coordinates over the radius.
virtual void set_legendre_base_p_spher(Base_spectral &) const
Gives the base using Legendre polynomials, for the component of a vector.
virtual Val_domain der_normal(const Val_domain &, int) const
Normal derivative with respect to a given surface.
void set_legendre_base_forz_cart(Base_spectral &so) const
Sets the base to the standard one for Legendre polynomials for a field like the component of a vecto...
virtual void set_val_inf(Val_domain &so, double xx) const
Sets the value at infinity of a Val_domain : not implemented for this type of Domain.
virtual void save(FILE *) const
Saving function.
virtual void export_tau(const Tensor &, int, int, Array< double > &, int &, const Array< int > &, int n_cmp=-1, Array< int > **p_cmp=0x0) const
Exports all the residual equations corresponding to a tensorial one in the bulk.
virtual Val_domain div_sin_theta(const Val_domain &) const
Division by .
virtual void set_legendre_base_y_cart(Base_spectral &) const
Gives the base using Legendre polynomials, for the component of a vector.
virtual void set_legendre_base_t_spher(Base_spectral &) const
Gives the base using Legendre polynomials, for the component of a vector.
virtual void set_cheb_base_r_spher(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the radial component of a vector.
virtual void do_cart() const
Computes the Cartesian coordinates.
virtual void find_other_dom(int, int, int &, int &) const
Gives the informations corresponding the a touching neighboring domain.
virtual Array< int > nbr_conditions_boundary(const Tensor &, int, int, int n_cmp=-1, Array< int > **p_cmp=0x0) const
Computes number of discretized equations associated with a given tensorial equation on a boundary.
virtual void set_cheb_base_p_spher(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
virtual Val_domain der_r(const Val_domain &) const
Compute the radial derivative of a scalar field.
virtual double val_boundary(int, const Val_domain &, const Index &) const
Computes the value of a field at a boundary.
virtual void set_cheb_base_y_cart(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
virtual void set_legendre_base_x_cart(Base_spectral &) const
Gives the base using Legendre polynomials, for the component of a vector.
virtual void set_legendre_base_r_spher(Base_spectral &) const
Gives the base using Legendre polynomials, for the radial component of a vector.
virtual void set_cheb_base_z_cart(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
void set_cheb_base_forx_cart(Base_spectral &so) const
Sets the base to the standard one for Chebyshev polynomials for a field like the component of a vect...
virtual void set_legendre_base_z_cart(Base_spectral &) const
Gives the base using Legendre polynomials, for the component of a vector.
virtual Tensor change_basis_cart_to_spher(int dd, const Tensor &) const
Changes the tensorial basis from Cartsian to spherical in a given domain.
void affecte_tau_one_coef_val_domain(Val_domain &so, int cc, int &pos_cf) const
Sets at most one coefficient of a Val_domain to 1.
virtual double integ_volume(const Val_domain &so) const
Volume integral.
virtual int nbr_unknowns(const Tensor &, int) const
Computes the number of true unknowns of a Tensor, in a given domain.
virtual int give_place_var(char *) const
Translates a name of a coordinate into its corresponding numerical name.
Point center
Absolute coordinates of the center.
virtual const Point absol_to_num_bound(const Point &, int) const
Computes the numerical coordinates from the physical ones for a point lying on a boundary.
int nbr_conditions_val_domain_boundary(const Val_domain &eq) const
Computes number of discretized equations associated with a given equation on a boundary.
void set_legendre_base_fory_cart(Base_spectral &so) const
Sets the base to the standard one for Legendre polynomials for a field like the component of a vecto...
void export_tau_val_domain_boundary(const Val_domain &eq, int bound, Array< double > &res, int &pos_res, int ncond) const
Exports all the residual equations corresponding to a tensorial one on a given boundary It makes use ...
virtual Val_domain mult_xm1(const Val_domain &) const
Multiplication by .
virtual Term_eq derive_flat_spher(int, char, const Term_eq &, const Metric *) const
Computes the flat derivative of a Term_eq, in spherical orthonormal coordinates.
void set_cheb_base_forz_cart(Base_spectral &so) const
Sets the base to the standard one for Chebyshev polynomials for a field like the component of a vect...
virtual Val_domain mult_r(const Val_domain &) const
Multiplication by .
virtual void export_tau_boundary(const Tensor &, int, int, Array< double > &, int &, const Array< int > &, int n_cmp=-1, Array< int > **p_cmp=0x0) const
Exports all the residual equations corresponding to a tensorial one on a given boundary It makes use ...
virtual void do_coloc()
Computes the colocation points.
int nbr_conditions_val_domain(const Val_domain &eq, int order) const
Computes number of discretized equations associated with a given tensorial equation in the bulk.
virtual bool is_in(const Point &xx, double prec=1e-13) const
Check whether a point lies inside Domain.
virtual void do_absol() const
Computes the absolute coordinates.
virtual void set_cheb_base_t_spher(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
virtual void affecte_tau(Tensor &, int, const Array< double > &, int &) const
Affects some coefficients to a Tensor.
virtual Val_domain mult_sin_phi(const Val_domain &) const
Multiplication by .
void set_cheb_base_fory_cart(Base_spectral &so) const
Sets the base to the standard one for Chebyshev polynomials for a field like the component of a vect...
virtual Val_domain div_r(const Val_domain &) const
Division by .
int nbr_unknowns_val_domain(const Val_domain &so) const
Computes the number of true unknowns of a Val_domain.
Domain_compact_symphi(int num, int ttype, double r_int, const Point &cr, const Dim_array &nbr)
Standard constructor :
void affecte_tau_val_domain(Val_domain &so, const Array< double > &cf, int &pos_cf) const
Affects some coefficients to a Val_domain.
virtual void do_der_abs_from_der_var(const Val_domain *const *const der_var, Val_domain **const der_abs) const
Computes the derivative with respect to the absolute Cartesian coordinates from the derivative with r...
virtual void set_cheb_base_x_cart(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
virtual Point get_center() const
Returns the center.
void export_tau_val_domain(const Val_domain &eq, int order, Array< double > &res, int &pos_res, int ncond) const
Exports a residual equation in the bulk.
virtual Array< int > nbr_conditions(const Tensor &, int, int, int n_cmp=-1, Array< int > **p_cmp=0x0) const
Computes number of discretized equations associated with a given tensorial equation in the bulk.
virtual Term_eq derive_flat_cart(int, char, const Term_eq &, const Metric *) const
Computes the flat derivative of a Term_eq, in Cartesian coordinates.
virtual Val_domain der_r_rtwo(const Val_domain &) const
Compute the radial derivative multiplied by of a scalar field.
Class for a spherical domain containing the origin a symmetry with respect to the plane and an quadr...
virtual Val_domain div_sin_theta(const Val_domain &) const
Division by .
void set_cheb_base_fory_cart(Base_spectral &so) const
Sets the base to the standard one for Chebyshev polynomials for a field like the component of a vect...
virtual Val_domain div_cos_theta(const Val_domain &) const
Division by .
virtual void set_cheb_base_x_cart(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
int nbr_unknowns_val_domain(const Val_domain &so) const
Computes the number of true unknowns of a Val_domain.
virtual void do_coloc()
Computes the colocation points.
virtual void do_cart() const
Computes the Cartesian coordinates.
virtual Point get_center() const
Returns the center.
void export_tau_val_domain_boundary(const Val_domain &eq, int bound, Array< double > &res, int &pos_res, int ncond) const
Exports all the residual equations corresponding to a tensorial one on a given boundary It makes use ...
virtual void set_cheb_base_yz_cart(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
virtual Val_domain ddp(const Val_domain &) const
Compute the second derivative with respect to of a scalar field.
virtual void do_cart_surr() const
Computes the Cartesian coordinates over the radius.
void affecte_tau_val_domain(Val_domain &so, const Array< double > &cf, int &pos_cf) const
Affects some coefficients to a Val_domain.
void filter_val_domain(Val_domain &so, double threshold) const
Sets to zero all the coefficients smaller than a given treshold.
virtual void set_legendre_base_y_cart(Base_spectral &) const
Gives the base using Legendre polynomials, for the component of a vector.
virtual Term_eq derive_flat_cart(int, char, const Term_eq &, const Metric *) const
Computes the flat derivative of a Term_eq, in Cartesian coordinates.
virtual Val_domain mult_sin_theta(const Val_domain &) const
Multiplication by .
virtual void do_der_abs_from_der_var(const Val_domain *const *const der_var, Val_domain **const der_abs) const
Computes the derivative with respect to the absolute Cartesian coordinates from the derivative with r...
virtual double integ(const Val_domain &so, int bound) const
Surface integral on a given boundary.
virtual int give_place_var(char *) const
Translates a name of a coordinate into its corresponding numerical name.
virtual Val_domain div_x(const Val_domain &) const
Division by .
virtual void set_cheb_base_xz_cart(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
int mult_base_r_int(int basea, int baseb) const
Multiply two radial basis.
void export_tau_val_domain(const Val_domain &eq, int order, Array< double > &res, int &pos_res, int ncond) const
Exports a residual equation in the bulk.
virtual Term_eq derive_flat_spher(int, char, const Term_eq &, const Metric *) const
Computes the flat derivative of a Term_eq, in spherical orthonormal coordinates.
virtual Val_domain mult_cos_theta(const Val_domain &) const
Multiplication by .
virtual Val_domain mult_sin_phi(const Val_domain &) const
Multiplication by .
virtual ostream & print(ostream &o) const
Delegate function to virtualize the << operator.
virtual Val_domain der_r(const Val_domain &) const
Compute the radial derivative of a scalar field.
int nbr_conditions_val_domain(const Val_domain &eq, int order) const
Computes number of discretized equations associated with a given tensorial equation in the bulk.
virtual Val_domain srdr(const Val_domain &) const
Compute the of a scalar field .
virtual void save(FILE *) const
Saving function.
virtual const Point absol_to_num_bound(const Point &, int) const
Computes the numerical coordinates from the physical ones for a point lying on a boundary.
virtual double get_rmax() const
Returns the maximum radius.
virtual Val_domain mult_r(const Val_domain &) const
Multiplication by .
virtual void find_other_dom(int, int, int &, int &) const
Gives the informations corresponding the a touching neighboring domain.
virtual int nbr_unknowns(const Tensor &, int) const
Computes the number of true unknowns of a Tensor, in a given domain.
virtual void set_legendre_base_p_spher(Base_spectral &) const
Gives the base using Legendre polynomials, for the component of a vector.
virtual void filter(Tensor &tt, int dom, double threshold) const
Puts to zero all the coefficients below a given treshold.
virtual Array< int > nbr_conditions(const Tensor &, int, int, int n_cmp=-1, Array< int > **p_cmp=0x0) const
Computes number of discretized equations associated with a given tensorial equation in the bulk.
virtual const Point absol_to_num(const Point &xxx) const
Computes the numerical coordinates from the physical ones.
void set_cheb_base_forr(Base_spectral &so) const
Sets the base to the standard one for Chebyshev polynomials for a field like the radius .
void set_cheb_base_forx_cart(Base_spectral &so) const
Sets the base to the standard one for Chebyshev polynomials for a field like the component of a vect...
int mult_base_angle_int(int basea, int baseb) const
Multiply two angular basis.
virtual void set_cheb_base_y_cart(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
virtual Val_domain mult_cos_phi(const Val_domain &) const
Multiplication by .
double alpha
Relates the numerical to the physical radii.
virtual void set_legendre_base_x_cart(Base_spectral &) const
Gives the base using Legendre polynomials, for the component of a vector.
virtual void do_absol() const
Computes the absolute coordinates.
virtual bool is_in(const Point &xx, double prec=1e-13) const
Check whether a point lies inside Domain.
virtual void set_legendre_base_t_spher(Base_spectral &) const
Gives the base using Legendre polynomials, for the component of a vector.
virtual void set_cheb_base_p_spher(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
int nbr_conditions_val_domain_boundary(const Val_domain &eq) const
Computes number of discretized equations associated with a given equation on a boundary.
virtual void set_legendre_base(Base_spectral &) const
Gives the standard base for Legendre polynomials.
virtual void set_cheb_base_r_spher(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the radial component of a vector.
void set_legendre_base_forr(Base_spectral &so) const
Sets the base to the standard one for Legendre polynomials for a field like the radius .
Domain_nucleus_symphi(int num, int ttype, double radius, const Point &cr, const Dim_array &nbr)
Standard constructor :
void set_legendre_base_forx_cart(Base_spectral &so) const
Sets the base to the standard one for Legendre polynomials for a field like the component of a vecto...
virtual void set_cheb_base_rp_spher(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a 2-tensor.
void set_legendre_base_fory_cart(Base_spectral &so) const
Sets the base to the standard one for Legendre polynomials for a field like the component of a vecto...
virtual void set_cheb_base_xy_cart(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
virtual void affecte_tau_one_coef(Tensor &, int, int, int &) const
Sets at most one coefficient of a Tensor to 1.
virtual void set_cheb_base_z_cart(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
virtual void export_tau_boundary(const Tensor &, int, int, Array< double > &, int &, const Array< int > &, int n_cmp=-1, Array< int > **p_cmp=0x0) const
Exports all the residual equations corresponding to a tensorial one on a given boundary It makes use ...
virtual Tensor change_basis_spher_to_cart(int dd, const Tensor &) const
Changes the tensorial basis from spherical to Cartesian in a given domain.
virtual Array< int > nbr_conditions_boundary(const Tensor &, int, int, int n_cmp=-1, Array< int > **p_cmp=0x0) const
Computes number of discretized equations associated with a given tensorial equation on a boundary.
virtual void set_legendre_base_z_cart(Base_spectral &) const
Gives the base using Legendre polynomials, for the component of a vector.
virtual void set_cheb_base(Base_spectral &) const
Gives the standard base for Chebyshev polynomials.
virtual void set_cheb_base_t_spher(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
Point center
Absolute coordinates of the center.
virtual void export_tau(const Tensor &, int, int, Array< double > &, int &, const Array< int > &, int n_cmp=-1, Array< int > **p_cmp=0x0) const
Exports all the residual equations corresponding to a tensorial one in the bulk.
virtual double val_boundary(int, const Val_domain &, const Index &) const
Computes the value of a field at a boundary.
virtual Val_domain der_normal(const Val_domain &, int) const
Normal derivative with respect to a given surface.
virtual void set_cheb_base_tp_spher(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a 2-tensor.
virtual Base_spectral mult(const Base_spectral &, const Base_spectral &) const
Method for the multiplication of two Base_spectral.
virtual Val_domain div_r(const Val_domain &) const
Division by .
virtual void set_legendre_base_r_spher(Base_spectral &) const
Gives the base using Legendre polynomials, for the radial component of a vector.
virtual void set_cheb_base_rt_spher(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a 2-tensor.
virtual void do_radius() const
Computes the generalized radius.
void set_cheb_base_forz_cart(Base_spectral &so) const
Sets the base to the standard one for Chebyshev polynomials for a field like the component of a vect...
virtual Val_domain div_1mx2(const Val_domain &) const
Division by .
virtual void affecte_tau(Tensor &, int, const Array< double > &, int &) const
Affects some coefficients to a Tensor.
virtual double integ_volume(const Val_domain &so) const
Volume integral.
void set_legendre_base_forz_cart(Base_spectral &so) const
Sets the base to the standard one for Legendre polynomials for a field like the component of a vecto...
virtual Tensor change_basis_cart_to_spher(int dd, const Tensor &) const
Changes the tensorial basis from Cartsian to spherical in a given domain.
void affecte_tau_one_coef_val_domain(Val_domain &so, int cc, int &pos_cf) const
Sets at most one coefficient of a Val_domain to 1.
Class for a spherical shell and a symmetry with respect to the plane and an quadrant symmetry wrt .
void export_tau_val_domain_boundary_exception(const Val_domain &eq, int mlim, int bound, Array< double > &res, int &pos_res, int ncond, const Param &param, int type_exception, const Val_domain &exception) const
Exports all the residual equations corresponding to one tensorial one on a given boundary,...
virtual Array< int > nbr_conditions_boundary(const Tensor &, int, int, int n_cmp=-1, Array< int > **p_cmp=0x0) const
Computes number of discretized equations associated with a given tensorial equation on a boundary.
double alpha
Relates the numerical to the physical radii.
virtual void find_other_dom(int, int, int &, int &) const
Gives the informations corresponding the a touching neighboring domain.
virtual const Point absol_to_num(const Point &xxx) const
Computes the numerical coordinates from the physical ones.
virtual double get_rmax() const
Returns the maximum radius.
virtual Val_domain ddp(const Val_domain &) const
Compute the second derivative with respect to of a scalar field.
double beta
Relates the numerical to the physical radii.
virtual void set_legendre_base_r_spher(Base_spectral &) const
Gives the base using Legendre polynomials, for the radial component of a vector.
virtual void set_legendre_base_z_cart(Base_spectral &) const
Gives the base using Legendre polynomials, for the component of a vector.
virtual void set_cheb_base_z_cart(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
virtual ~Domain_shell_symphi()
Destructor.
int mult_base_angle_int(int basea, int baseb) const
Multiply two angular basis.
virtual Val_domain der_normal(const Val_domain &, int) const
Normal derivative with respect to a given surface.
virtual void set_cheb_base_p_spher(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
virtual void export_tau_boundary(const Tensor &, int, int, Array< double > &, int &, const Array< int > &, int n_cmp=-1, Array< int > **p_cmp=0x0) const
Exports all the residual equations corresponding to a tensorial one on a given boundary It makes use ...
virtual void do_der_abs_from_der_var(const Val_domain *const *const der_var, Val_domain **const der_abs) const
Computes the derivative with respect to the absolute Cartesian coordinates from the derivative with r...
virtual Val_domain mult_xm1(const Val_domain &) const
Multiplication by .
virtual int nbr_unknowns(const Tensor &, int) const
Computes the number of true unknowns of a Tensor, in a given domain.
void set_cheb_base_forx_cart(Base_spectral &so) const
Sets the base to the standard one for Chebyshev polynomials for a field like the component of a vect...
void set_legendre_base_forz_cart(Base_spectral &so) const
Sets the base to the standard one for Legendre polynomials for a field like the component of a vecto...
virtual double integ(const Val_domain &so, int bound) const
Surface integral on a given boundary.
virtual Val_domain div_r(const Val_domain &) const
Division by .
int nbr_conditions_val_domain(const Val_domain &eq, int order) const
Computes number of discretized equations associated with a given tensorial equation in the bulk.
virtual void set_cheb_base(Base_spectral &) const
Gives the standard base for Chebyshev polynomials.
virtual Val_domain mult_r(const Val_domain &) const
Multiplication by .
void set_legendre_base_fory_cart(Base_spectral &so) const
Sets the base to the standard one for Legendre polynomials for a field like the component of a vecto...
virtual Tensor change_basis_spher_to_cart(int dd, const Tensor &) const
Changes the tensorial basis from spherical to Cartesian in a given domain.
virtual Val_domain der_r(const Val_domain &) const
Compute the radial derivative of a scalar field.
virtual void set_cheb_base_tp_spher(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a 2-tensor.
virtual void set_legendre_base_y_cart(Base_spectral &) const
Gives the base using Legendre polynomials, for the component of a vector.
virtual Val_domain mult_sin_phi(const Val_domain &) const
Multiplication by .
void set_legendre_base_forr(Base_spectral &so) const
Sets the base to the standard one for Legendre polynomials for a field like the radius .
virtual Val_domain div_sin_theta(const Val_domain &) const
Division by .
virtual void set_cheb_base_xz_cart(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
virtual void set_cheb_base_r_spher(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the radial component of a vector.
void set_cheb_base_fory_cart(Base_spectral &so) const
Sets the base to the standard one for Chebyshev polynomials for a field like the component of a vect...
virtual Term_eq derive_flat_cart(int, char, const Term_eq &, const Metric *) const
Computes the flat derivative of a Term_eq, in Cartesian coordinates.
void affecte_tau_one_coef_val_domain(Val_domain &so, int cc, int &pos_cf) const
Sets at most one coefficient of a Val_domain to 1.
virtual Val_domain div_cos_theta(const Val_domain &) const
Division by .
virtual int give_place_var(char *) const
Translates a name of a coordinate into its corresponding numerical name.
virtual void do_coloc()
Computes the colocation points.
virtual void set_legendre_base(Base_spectral &) const
Gives the standard base for Legendre polynomials.
virtual void set_legendre_base_t_spher(Base_spectral &) const
Gives the base using Legendre polynomials, for the component of a vector.
void affecte_tau_val_domain(Val_domain &so, const Array< double > &cf, int &pos_cf) const
Affects some coefficients to a Val_domain.
virtual Val_domain mult_1mrsL(const Val_domain &) const
Multiplication by .
virtual double val_boundary(int, const Val_domain &, const Index &) const
Computes the value of a field at a boundary.
virtual void export_tau_boundary_exception(const Tensor &, int, int, Array< double > &, int &, const Array< int > &, const Param &, int, const Tensor &, int n_cmp=-1, Array< int > **p_cmp=0x0) const
Exports all the residual equations corresponding to one tensorial one on a given boundary,...
virtual Val_domain der_partial_var(const Val_domain &, int) const
Partial derivative with respect to a coordinate.
virtual void set_legendre_base_x_cart(Base_spectral &) const
Gives the base using Legendre polynomials, for the component of a vector.
Point center
Absolute coordinates of the center.
virtual Array< int > nbr_conditions(const Tensor &, int, int, int n_cmp=-1, Array< int > **p_cmp=0x0) const
Computes number of discretized equations associated with a given tensorial equation in the bulk.
virtual void set_cheb_base_yz_cart(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
virtual void affecte_tau_one_coef(Tensor &, int, int, int &) const
Sets at most one coefficient of a Tensor to 1.
virtual void set_cheb_base_rp_spher(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a 2-tensor.
virtual void do_cart() const
Computes the Cartesian coordinates.
virtual void set_cheb_base_t_spher(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
virtual void set_cheb_base_x_cart(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
virtual void export_tau(const Tensor &, int, int, Array< double > &, int &, const Array< int > &, int n_cmp=-1, Array< int > **p_cmp=0x0) const
Exports all the residual equations corresponding to a tensorial one in the bulk.
void set_legendre_base_forx_cart(Base_spectral &so) const
Sets the base to the standard one for Legendre polynomials for a field like the component of a vecto...
virtual void set_cheb_base_xy_cart(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
virtual Base_spectral mult(const Base_spectral &, const Base_spectral &) const
Method for the multiplication of two Base_spectral.
void export_tau_val_domain(const Val_domain &eq, int order, Array< double > &res, int &pos_res, int ncond) const
Exports a residual equation in the bulk.
int nbr_unknowns_val_domain(const Val_domain &so) const
Computes the number of true unknowns of a Val_domain.
Domain_shell_symphi(int num, int ttype, double r_int, double r_ext, const Point &cr, const Dim_array &nbr)
Standard constructor :
virtual void set_cheb_base_y_cart(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
virtual Tensor change_basis_cart_to_spher(int dd, const Tensor &) const
Changes the tensorial basis from Cartsian to spherical in a given domain.
virtual Term_eq derive_flat_spher(int, char, const Term_eq &, const Metric *) const
Computes the flat derivative of a Term_eq, in spherical orthonormal coordinates.
void export_tau_val_domain_boundary(const Val_domain &eq, int bound, Array< double > &res, int &pos_res, int ncond) const
Exports all the residual equations corresponding to a tensorial one on a given boundary It makes use ...
int nbr_conditions_val_domain_boundary(const Val_domain &eq) const
Computes number of discretized equations associated with a given equation on a boundary.
virtual void save(FILE *) const
Saving function.
virtual void set_legendre_base_p_spher(Base_spectral &) const
Gives the base using Legendre polynomials, for the component of a vector.
virtual void set_cheb_base_rt_spher(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a 2-tensor.
void set_cheb_base_forr(Base_spectral &so) const
Sets the base to the standard one for Chebyshev polynomials for a field like the radius .
virtual ostream & print(ostream &o) const
Delegate function to virtualize the << operator.
virtual Val_domain div_xm1(const Val_domain &) const
Division by .
virtual Val_domain div_1mx2(const Val_domain &) const
Division by .
virtual Point get_center() const
Returns the center.
virtual bool is_in(const Point &xx, double prec=1e-13) const
Check whether a point lies inside Domain.
virtual Val_domain mult_cos_phi(const Val_domain &) const
Multiplication by .
virtual Val_domain div_1mrsL(const Val_domain &) const
Division by .
virtual void do_cart_surr() const
Computes the Cartesian coordinates over the radius.
virtual void do_absol() const
Computes the absolute coordinates.
virtual void do_radius() const
Computes the generalized radius.
virtual double integ_volume(const Val_domain &so) const
Volume integral.
void set_cheb_base_forz_cart(Base_spectral &so) const
Sets the base to the standard one for Chebyshev polynomials for a field like the component of a vect...
virtual void affecte_tau(Tensor &, int, const Array< double > &, int &) const
Affects some coefficients to a Tensor.
virtual Val_domain mult_sin_theta(const Val_domain &) const
Multiplication by .
virtual Val_domain mult_cos_theta(const Val_domain &) const
Multiplication by .
virtual const Point absol_to_num_bound(const Point &, int) const
Computes the numerical coordinates from the physical ones for a point lying on a boundary.
Abstract class that implements the fonctionnalities common to all the type of domains.
Definition: space.hpp:60
Val_domain * radius
The generalized radius.
Definition: space.hpp:78
Class that gives the position inside a multi-dimensional Array.
Definition: index.hpp:38
Purely abstract class for metric handling.
Definition: metric.hpp:39
Parameter storage.
Definition: param.hpp:30
The class Point is used to store the coordinates of a point.
Definition: point.hpp:30
The Space_spheric_symphi_symphi class fills the space with spherical domain(s) with symmetry in phi.
virtual void save(FILE *) const
Saving function.
Space_spheric_symphi(int ttype, const Point &cr, const Dim_array &nbr, double bound)
Standard constructor with only one nucleus.
void add_eq_ori(System_of_eqs &syst, const char *eq)
Adds an equation being the value of some field at the origin.
The Space class is an ensemble of domains describing the whole space of the computation.
Definition: space.hpp:1362
Class used to describe and solve a system of equations.
Tensor handling.
Definition: tensor.hpp:149
This class is intended to describe the manage objects appearing in the equations.
Definition: term_eq.hpp:62
Class for storing the basis of decompositions of a field and its values on both the configuration and...
Definition: val_domain.hpp:69