KADATH
adapted_polar.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 
21 #ifndef __ADAPTED_POLAR_HPP_
22 #define __ADAPTED_POLAR_HPP_
23 
24 #include "space.hpp"
25 #include "term_eq.hpp"
26 #include "polar.hpp"
27 #include "val_domain.hpp"
28 
29 namespace Kadath {
30 
51 
52  private:
56  const Space& sp ;
59  double outer_radius ;
63  mutable Term_eq* normal_spher ;
67  mutable Term_eq* normal_cart ;
71  mutable Term_eq* rad_term_eq ;
80 
82 
83  public:
94  Domain_polar_shell_inner_adapted (const Space& sp, int num, int ttype, double rin, double rout, const Point& cr, const Dim_array& nbr) ;
105  Domain_polar_shell_inner_adapted (const Space& sp, int num, int ttype, const Val_domain& rin, double rout, const Point& cr, const Dim_array& nbr) ;
113  Domain_polar_shell_inner_adapted (const Space& sp, int num, FILE* fd) ;
114 
116  virtual void save (FILE*) const ;
121  void del_deriv() override;
122 
123  private:
124  virtual void do_absol () const ;
125  virtual void do_radius () const ;
126  virtual void do_cart () const ;
127  virtual void do_cart_surr () const ;
128 
129 
130  private:
131  virtual void set_cheb_base(Base_spectral&) const ;
132  virtual void set_legendre_base(Base_spectral&) const ;
133  virtual void set_anti_cheb_base(Base_spectral&) const ;
134  virtual void set_anti_legendre_base(Base_spectral&) const ;
135  virtual void set_cheb_base_with_m(Base_spectral&, int m) const ;
136  virtual void set_legendre_base_with_m(Base_spectral&, int m) const ;
137  virtual void set_anti_cheb_base_with_m(Base_spectral&, int m) const ;
138  virtual void set_anti_legendre_base_with_m(Base_spectral&, int m) const ;
139  virtual void set_cheb_r_base(Base_spectral&) const ;
140  virtual void set_legendre_r_base(Base_spectral&) const ;
141 
142  virtual void do_coloc () ;
143  virtual int give_place_var (char*) const ;
144 
145  virtual int nbr_unknowns_from_adapted() const ;
146  virtual void vars_to_terms() const ;
147  virtual void affecte_coef (int&, int, bool&) const ;
148  virtual void xx_to_vars_from_adapted (Val_domain&, const Array<double>&, int&) const ;
149  virtual void xx_to_ders_from_adapted (const Array<double>&, int&) const ;
150  virtual void update_term_eq (Term_eq*) const ;
151  virtual void update_variable (const Val_domain&, const Scalar&, Scalar&) const ;
152  virtual void update_constante (const Val_domain&, const Scalar&, Scalar&) const ;
153  virtual void update_mapping(const Val_domain&) ;
154 
155  public:
159  void set_mapping(const Val_domain& so) const ;
163  void update() const ;
164 
165  public:
166  virtual Point get_center () const {return center ;} ;
167  virtual bool is_in(const Point& xx, double prec=1e-13) const ;
168  virtual const Point absol_to_num(const Point&) const;
169  virtual void do_der_abs_from_der_var(const Val_domain *const *const der_var, Val_domain **const der_abs) const ;
170 
171  virtual Base_spectral mult (const Base_spectral&, const Base_spectral&) const ;
172 
173 
174  public:
175  virtual Val_domain mult_cos_theta (const Val_domain&) const ;
176  virtual Val_domain mult_sin_theta (const Val_domain&) const ;
177  virtual Val_domain div_sin_theta (const Val_domain&) const ;
178  virtual Val_domain mult_r (const Val_domain&) const ;
179  virtual Val_domain div_r (const Val_domain&) const ;
180  virtual Val_domain laplacian (const Val_domain&, int) const ;
181  virtual Val_domain laplacian2 (const Val_domain&, int) const ;
182  virtual Val_domain der_r (const Val_domain&) const ;
183 
184 
185  virtual double val_boundary (int, const Val_domain&, const Index&) const ;
186  virtual void find_other_dom (int, int, int&, int&) const ;
187  virtual Val_domain der_normal (const Val_domain&, int) const ;
188 
194  Term_eq flat_grad_spher (const Term_eq& so) const ;
195 
196  virtual Term_eq partial_spher (const Term_eq&) const ;
197  virtual Term_eq partial_cart (const Term_eq&) const ;
198  virtual const Term_eq* give_normal(int, int) const ;
199  virtual Term_eq grad_term_eq (const Term_eq&) const ;
200 
206  Term_eq derive_r (const Term_eq& so) const ;
212  Term_eq derive_t (const Term_eq& so) const ;
216  void do_normal_spher () const ;
220  void do_normal_cart () const ;
221 
222  virtual Term_eq der_normal_term_eq (const Term_eq&, int) const ;
223  virtual Term_eq dr_term_eq (const Term_eq&) const ;
224  virtual Term_eq lap_term_eq (const Term_eq&, int) const ;
225  virtual Term_eq lap2_term_eq (const Term_eq&, int) const ;
226  virtual Term_eq mult_r_term_eq (const Term_eq&) const ;
227  virtual Term_eq div_r_term_eq (const Term_eq&) const ;
228  virtual Term_eq integ_term_eq (const Term_eq&, int) const ;
229 
230  virtual int nbr_unknowns (const Tensor&, int) const ;
238  int nbr_unknowns_val_domain (const Val_domain& so, int mquant) const ;
239  virtual Array<int> nbr_conditions (const Tensor&, int, int, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
248  int nbr_conditions_val_domain (const Val_domain& so, int mquant, int order) const ;
249  virtual Array<int> nbr_conditions_boundary (const Tensor&, int, int, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
258  int nbr_conditions_val_domain_boundary (const Val_domain& eq, int mquant) const ;
259  virtual void export_tau (const Tensor&, int, int, Array<double>&, int&, const Array<int>&, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
270  void export_tau_val_domain (const Val_domain& eq, int mquant, int order, Array<double>& res, int& pos_res, int ncond) const ;
271  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 ;
282  void export_tau_val_domain_boundary (const Val_domain& eq, int mquant, int bound, Array<double>& res, int& pos_res, int ncond) const ;
283  virtual void affecte_tau (Tensor&, int, const Array<double>&, int&) const ;
292  void affecte_tau_val_domain (Val_domain& so, int mquant, const Array<double>& cf, int& pos_cf) const ;
293  virtual void affecte_tau_one_coef (Tensor&, int, int, int&) const ;
302  void affecte_tau_one_coef_val_domain (Val_domain& so, int mquant, int cc, int& pos_cf) const ;
303 
304 public:
305  virtual ostream& print (ostream& o) const ;
306 
307  friend class Space_polar_adapted ;
308 } ;
309 
329 
330  private:
334  const Space& sp ;
337  double inner_radius ;
341  mutable Term_eq* normal_spher ;
345  mutable Term_eq* normal_cart ;
349  mutable Term_eq* rad_term_eq ;
358 
360 
361  public:
372  Domain_polar_shell_outer_adapted (const Space& sp, int num, int ttype, double rin, double rout, const Point& cr, const Dim_array& nbr) ;
383  Domain_polar_shell_outer_adapted (const Space& sp, int num, int ttype, double rin, const Val_domain& rout, const Point& cr, const Dim_array& nbr) ;
391  Domain_polar_shell_outer_adapted (const Space& sp, int num, FILE* fd) ;
392 
393  virtual ~Domain_polar_shell_outer_adapted () ;
394 
395  void del_deriv() override;
396 
397  virtual void save (FILE*) const ;
402  private:
403  virtual void do_absol () const ;
404  virtual void do_radius () const ;
405  virtual void do_cart () const ;
406  virtual void do_cart_surr () const ;
407 
408 
409  private:
410  virtual void set_cheb_base(Base_spectral&) const ;
411  virtual void set_legendre_base(Base_spectral&) const ;
412  virtual void set_anti_cheb_base(Base_spectral&) const ;
413  virtual void set_anti_legendre_base(Base_spectral&) const ;
414  virtual void set_cheb_base_with_m(Base_spectral&, int m) const ;
415  virtual void set_legendre_base_with_m(Base_spectral&, int m) const ;
416  virtual void set_anti_cheb_base_with_m(Base_spectral&, int m) const ;
417  virtual void set_anti_legendre_base_with_m(Base_spectral&, int m) const ;
418  virtual void set_cheb_r_base(Base_spectral&) const ;
419  virtual void set_legendre_r_base(Base_spectral&) const ;
420 
421  virtual void do_coloc() ;
422 
423 
424  virtual int give_place_var (char*) const ;
425 
426  virtual int nbr_unknowns_from_adapted() const ;
427  virtual void vars_to_terms() const ;
428  virtual void affecte_coef (int&, int, bool&) const ;
429  virtual void xx_to_vars_from_adapted (Val_domain&, const Array<double>&, int&) const ;
430  virtual void xx_to_ders_from_adapted (const Array<double>&, int&) const ;
431  virtual void update_term_eq (Term_eq*) const ;
432  virtual void update_variable (const Val_domain&, const Scalar&, Scalar&) const ;
433  virtual void update_constante (const Val_domain&, const Scalar&, Scalar&) const ;
434  virtual void update_mapping(const Val_domain&) ;
435 
436  public:
440  void set_mapping(const Val_domain& so) const ;
444  void update() const ;
445 
446  public:
447 
448  virtual Point get_center () const {return center ;} ;
449 
450  virtual bool is_in(const Point& xx, double prec=1e-13) const ;
451  virtual const Point absol_to_num(const Point&) const;
452  virtual void do_der_abs_from_der_var(const Val_domain *const *const der_var, Val_domain **const der_abs) const ;
453 
454  virtual Base_spectral mult (const Base_spectral&, const Base_spectral&) const ;
455 
456  public:
457 
458  virtual Val_domain mult_cos_theta (const Val_domain&) const ;
459  virtual Val_domain mult_sin_theta (const Val_domain&) const ;
460  virtual Val_domain div_sin_theta (const Val_domain&) const ;
461 
462  virtual Val_domain mult_r (const Val_domain&) const ;
463  virtual Val_domain div_r (const Val_domain&) const ;
464  virtual Val_domain laplacian (const Val_domain&, int) const ;
465  virtual Val_domain laplacian2 (const Val_domain&, int) const ;
466  virtual Val_domain der_r (const Val_domain&) const ;
467 
473  Term_eq flat_grad_spher (const Term_eq&) const ;
474 
475  virtual Term_eq partial_spher (const Term_eq&) const ;
476  virtual Term_eq partial_cart (const Term_eq&) const ;
477  virtual const Term_eq* give_normal(int, int) const ;
478  virtual Term_eq grad_term_eq (const Term_eq&) const ;
479  virtual Term_eq integ_term_eq (const Term_eq&, int) const ;
485  Term_eq derive_r (const Term_eq& so) const ;
491  Term_eq derive_t (const Term_eq& so) const ;
495  void do_normal_spher () const ;
499  void do_normal_cart () const ;
500 
501  virtual Term_eq der_normal_term_eq (const Term_eq&, int) const ;
502  virtual Term_eq dr_term_eq (const Term_eq&) const ;
503  virtual Term_eq lap_term_eq (const Term_eq&, int) const ;
504  virtual Term_eq lap2_term_eq (const Term_eq&, int) const ;
505  virtual Term_eq mult_r_term_eq (const Term_eq&) const ;
506  virtual Term_eq div_r_term_eq (const Term_eq&) const ;
507 
508  virtual double val_boundary (int, const Val_domain&, const Index&) const ;
509  virtual void find_other_dom (int, int, int&, int&) const ;
510  virtual Val_domain der_normal (const Val_domain&, int) const ;
511 
512  virtual int nbr_unknowns (const Tensor&, int) const ;
520  int nbr_unknowns_val_domain (const Val_domain& so, int mquant) const ;
521  virtual Array<int> nbr_conditions (const Tensor&, int, int, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
530  int nbr_conditions_val_domain (const Val_domain& so, int mquant, int order) const ;
531  virtual Array<int> nbr_conditions_boundary (const Tensor&, int, int, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
540  int nbr_conditions_val_domain_boundary (const Val_domain& eq, int mquant) const ;
541  virtual void export_tau (const Tensor&, int, int, Array<double>&, int&, const Array<int>&, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
552  void export_tau_val_domain (const Val_domain& eq, int mquant, int order, Array<double>& res, int& pos_res, int ncond) const ;
553  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 ;
564  void export_tau_val_domain_boundary (const Val_domain& eq, int mquant, int bound, Array<double>& res, int& pos_res, int ncond) const ;
565  virtual void affecte_tau (Tensor&, int, const Array<double>&, int&) const ;
574  void affecte_tau_val_domain (Val_domain& so, int mquant, const Array<double>& cf, int& pos_cf) const ;
575  virtual void affecte_tau_one_coef (Tensor&, int, int, int&) const ;
584  void affecte_tau_one_coef_val_domain (Val_domain& so, int mquant, int cf, int& pos_cf) const ;
585 
586 
587 public:
588  virtual ostream& print (ostream& o) const ;
589 
590  friend class Space_polar_adapted ;
591 } ;
592 
598 class Space_polar_adapted : public Space {
599  public:
607  Space_polar_adapted (int ttyp, const Point& cr, const Dim_array& nbr, const Array<double>& bounds) ;
608  Space_polar_adapted (FILE*) ;
609  virtual ~Space_polar_adapted() ;
610  virtual void save(FILE*) const ;
611 
612  virtual int nbr_unknowns_from_variable_domains() const ;
613  virtual void affecte_coef_to_variable_domains(int& , int, Array<int>&) const ;
614  virtual void xx_to_ders_variable_domains(const Array<double>&, int&) const ;
615  virtual void xx_to_vars_variable_domains(System_of_eqs*, const Array<double>&, int&) const ;
616 
622  void add_eq_ori (System_of_eqs& syst, const char* eq) ;
632  void add_eq (System_of_eqs& syst, const char* eq, const char* rac, const char* rac_der, int nused=-1, Array<int>** pused=0x0) const ;
633 
634 } ;
635 }
636 #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-like domain, having a symmetry with respect to the plane .
void export_tau_val_domain(const Val_domain &eq, int mquant, int order, Array< double > &res, int &pos_res, int ncond) const
Exports a residual equation in the bulk.
int nbr_conditions_val_domain_boundary(const Val_domain &eq, int mquant) const
Computes number of discretized equations associated with a given equation on a boundary.
virtual void set_anti_cheb_base(Base_spectral &) const
Gives the base using Chebyshev polynomials, for functions antisymetric with respect to .
virtual Term_eq lap2_term_eq(const Term_eq &, int) const
Returns the flat 2d-Laplacian of Term_eq, for a given harmonic.
virtual Val_domain der_normal(const Val_domain &, int) const
Normal derivative with respect to a given surface.
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.
void set_mapping(const Val_domain &so) const
Affects the inner radius.
int nbr_conditions_val_domain(const Val_domain &so, int mquant, int order) const
Computes number of discretized equations associated with a given tensorial 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 der_normal_term_eq(const Term_eq &, int) const
Returns the normal derivative of a Term_eq.
virtual void do_cart_surr() const
Computes the Cartesian coordinates over the radius.
virtual const Point absol_to_num(const Point &) const
Computes the numerical coordinates from the physical ones.
Term_eq * der_rad_term_eq
Pointer on the Term_eq containing the .
virtual Val_domain laplacian(const Val_domain &, int) const
Computes the ordinary flat Laplacian for a scalar field with an harmonic index m.
virtual void xx_to_ders_from_adapted(const Array< double > &, int &) const
Affects the derivative part of variable a Domain from a set of values.
virtual void update_constante(const Val_domain &, const Scalar &, Scalar &) const
Update the value of a scalar, after the shape of the Domain has been changed by the system.
virtual void affecte_tau_one_coef(Tensor &, int, int, int &) const
Sets at most one coefficient of a Tensor to 1.
virtual int nbr_unknowns(const Tensor &, int) const
Computes the number of true unknowns of a Tensor, in a given domain.
virtual void set_cheb_base_with_m(Base_spectral &, int m) const
Gives the standard base using Chebyshev polynomials.
virtual void update_variable(const Val_domain &, const Scalar &, Scalar &) const
Update the value of a scalar, after the shape of the Domain has been changed by the system.
Val_domain get_inner_radius() const
Returns the inner variable boundary.
virtual Val_domain div_sin_theta(const Val_domain &) const
Division by .
Point center
Absolute coordinates of the center.
virtual void update_mapping(const Val_domain &)
Updates the variables parts of the Domain.
void del_deriv() override
Destroys the derivated members (like coloc, cart and radius), when changing the type of colocation po...
void do_normal_cart() const
Computes the normal wrt the inner boundary, in Cartesian coordinates.
Val_domain * inner_radius
Pointer on the inner boundary , as a Val_domain.
virtual void update_term_eq(Term_eq *) const
Update the value of a field, after the shape of the Domain has been changed by the system.
virtual double val_boundary(int, const Val_domain &, const Index &) const
Computes the value of a field at a boundary.
Term_eq * dt_rad_term_eq
Pointer on the Term_eq containing the .
Term_eq * rad_term_eq
Pointer on the Term_eq containing the radius.
virtual void set_legendre_base(Base_spectral &) const
Gives the standard base for Legendre polynomials.
Term_eq * inner_radius_term_eq
Pointer on the inner boundary , as a Term_eq.
virtual void vars_to_terms() const
The Term_eq describing the variable shape of the Domain are updated.
virtual void set_legendre_r_base(Base_spectral &) const
Gives the base using odd Legendre polynomials$ for the radius.
virtual void set_legendre_base_with_m(Base_spectral &, int m) const
Gives the stnadard base using Legendre polynomials.
virtual Val_domain der_r(const Val_domain &) const
Compute the radial derivative of a scalar field.
void update() const
Updates all the quantities that depend on the inner radius (like the normal vectors).
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 ...
Term_eq derive_t(const Term_eq &so) const
Computes .
const Space & sp
The corresponding Space ; required for updating fields whene the mapping changes.
virtual const Term_eq * give_normal(int, int) const
Returns the vector normal to a surface.
Term_eq * normal_cart
Pointer on the Term_eq containing the normal vector to the inner boundary, in Cartesian coordinates.
virtual void save(FILE *) const
Saving function.
virtual void xx_to_vars_from_adapted(Val_domain &, const Array< double > &, int &) const
Computes the new boundary of a Domain from a set of values.
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_anti_cheb_base_with_m(Base_spectral &, int m) const
Gives the base using Chebyshev polynomials, for functions antisymetric with respect to .
virtual Val_domain mult_sin_theta(const Val_domain &) const
Multiplication by .
virtual ostream & print(ostream &o) const
Delegate function to virtualize the << operator.
virtual Val_domain laplacian2(const Val_domain &, int) const
Computes the ordinary flat 2dè- Laplacian for a scalar field with an harmonic index m.
virtual void do_radius() const
Computes the generalized radius.
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 void set_cheb_r_base(Base_spectral &) const
Gives the base using odd Chebyshev polynomials$ for the radius.
void affecte_tau_val_domain(Val_domain &so, int mquant, const Array< double > &cf, int &pos_cf) const
Affects some coefficients to a Val_domain.
void export_tau_val_domain_boundary(const Val_domain &eq, int mquant, 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 Term_eq dr_term_eq(const Term_eq &) const
Radial derivative of a Term_eq.
Term_eq derive_r(const Term_eq &so) const
Computes .
virtual Val_domain mult_r(const Val_domain &) const
Multiplication by .
virtual Term_eq lap_term_eq(const Term_eq &, int) const
Returns the flat Laplacian of Term_eq, for a given harmonic.
void affecte_tau_one_coef_val_domain(Val_domain &so, int mquant, int cc, int &pos_cf) const
Sets at most one coefficient of a Val_domain to 1.
virtual void do_absol() const
Computes the absolute coordinates.
virtual void affecte_tau(Tensor &, int, const Array< double > &, int &) const
Affects some coefficients to a Tensor.
virtual void set_anti_legendre_base_with_m(Base_spectral &, int m) const
Gives the base using Legendre polynomials, for functions antisymetric with respect to .
virtual int give_place_var(char *) const
Translates a name of a coordinate into its corresponding numerical name.
virtual void affecte_coef(int &, int, bool &) const
The variation of the functions describing the shape of the Domain are affected from the unknowns of t...
virtual Term_eq div_r_term_eq(const Term_eq &) const
Division by of a Term_eq.
int nbr_unknowns_val_domain(const Val_domain &so, int mquant) const
Computes the number of true unknowns of a Val_domain.
virtual void set_anti_legendre_base(Base_spectral &) const
Gives the base using Legendre polynomials, for functions antisymetric with respect to .
virtual int nbr_unknowns_from_adapted() const
Gives the number of unknowns coming from the variable shape of the domain.
Domain_polar_shell_inner_adapted(const Space &sp, int num, int ttype, double rin, double rout, const Point &cr, const Dim_array &nbr)
Constructor :
Term_eq flat_grad_spher(const Term_eq &so) const
Computes the flat gradient of a field, in orthonormal spherical coordinates.
virtual Term_eq partial_spher(const Term_eq &) const
Computes the part of the gradient containing the partial derivative of the field, in spherical orthon...
virtual Term_eq partial_cart(const Term_eq &) const
Computes the part of the gradient containing the partial derivative of the field, in Cartesian coordi...
virtual void do_coloc()
Computes the colocation points.
virtual Term_eq mult_r_term_eq(const Term_eq &) const
Multiplication by of a Term_eq.
virtual void find_other_dom(int, int, int &, int &) const
Gives the informations corresponding the a touching neighboring domain.
Term_eq * normal_spher
Pointer on the Term_eq containing the normal vector to the inner boundary, in orthonormal spherical c...
virtual Base_spectral mult(const Base_spectral &, const Base_spectral &) const
Method for the multiplication of two Base_spectral.
virtual Term_eq integ_term_eq(const Term_eq &, int) const
Surface integral of a Term_eq.
virtual void set_cheb_base(Base_spectral &) const
Gives the standard base for Chebyshev polynomials.
virtual void do_cart() const
Computes the Cartesian coordinates.
virtual Val_domain div_r(const Val_domain &) const
Division by .
virtual Term_eq grad_term_eq(const Term_eq &) const
Gradient of Term_eq.
virtual Val_domain mult_cos_theta(const Val_domain &) const
Multiplication by .
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 do_normal_spher() const
Computes the normal wrt the inner boundary, in orthonormal spherical coordinates.
Class for a spherical-like domain, having a symmetry with respect to the plane .
virtual void vars_to_terms() const
The Term_eq describing the variable shape of the Domain are updated.
virtual void affecte_tau_one_coef(Tensor &, int, int, int &) const
Sets at most one coefficient of a Tensor to 1.
virtual ostream & print(ostream &o) const
Delegate function to virtualize the << operator.
virtual Point get_center() const
Returns the center.
virtual Term_eq lap_term_eq(const Term_eq &, int) const
Returns the flat Laplacian of Term_eq, for a given harmonic.
virtual Term_eq dr_term_eq(const Term_eq &) const
Radial derivative of a Term_eq.
Term_eq * normal_cart
Pointer on the Term_eq containing the normal vector to the outer boundary, in Cartesian coordinates.
virtual void set_legendre_base(Base_spectral &) const
Gives the standard base for Legendre polynomials.
Point center
Absolute coordinates of the center.
virtual Term_eq der_normal_term_eq(const Term_eq &, int) const
Returns the normal derivative of a Term_eq.
virtual Val_domain mult_r(const Val_domain &) const
Multiplication by .
virtual void update_variable(const Val_domain &, const Scalar &, Scalar &) const
Update the value of a scalar, after the shape of the Domain has been changed by the system.
virtual Term_eq integ_term_eq(const Term_eq &, int) const
Surface integral of a Term_eq.
void do_normal_cart() const
Computes the normal wrt the inner boundary, in Cartesian coordinates.
virtual Val_domain mult_cos_theta(const Val_domain &) const
Multiplication by .
int nbr_unknowns_val_domain(const Val_domain &so, int mquant) const
Computes the number of true unknowns of a Val_domain.
const Space & sp
The corresponding Space ; required for updating fields whene the mapping changes.
virtual Term_eq lap2_term_eq(const Term_eq &, int) const
Returns the flat 2d-Laplacian of Term_eq, for a given harmonic.
virtual void do_coloc()
Computes the colocation points.
virtual const Point absol_to_num(const Point &) const
Computes the numerical coordinates from the physical ones.
virtual Term_eq partial_spher(const Term_eq &) const
Computes the part of the gradient containing the partial derivative of the field, in spherical orthon...
virtual void set_anti_cheb_base(Base_spectral &) const
Gives the base using Chebyshev polynomials, for functions antisymetric with respect to .
virtual void xx_to_ders_from_adapted(const Array< double > &, int &) const
Affects the derivative part of variable a Domain from a set of values.
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 ...
Term_eq derive_r(const Term_eq &so) const
Computes .
virtual Val_domain laplacian2(const Val_domain &, int) const
Computes the ordinary flat 2dè- Laplacian for a scalar field with an harmonic index m.
void export_tau_val_domain(const Val_domain &eq, int mquant, int order, Array< double > &res, int &pos_res, int ncond) const
Exports a residual equation in the bulk.
virtual void save(FILE *) const
Saving function.
void do_normal_spher() const
Computes the normal wrt the inner boundary, in orthonormal spherical coordinates.
virtual Term_eq partial_cart(const Term_eq &) const
Computes the part of the gradient containing the partial derivative of the field, in Cartesian coordi...
Term_eq flat_grad_spher(const Term_eq &) const
Computes the flat gradient of a field, in orthonormal spherical coordinates.
Domain_polar_shell_outer_adapted(const Space &sp, int num, int ttype, double rin, double rout, const Point &cr, const Dim_array &nbr)
Constructor :
virtual Term_eq mult_r_term_eq(const Term_eq &) const
Multiplication by of a Term_eq.
virtual void update_mapping(const Val_domain &)
Updates the variables parts of the Domain.
Term_eq * dt_rad_term_eq
Pointer on the Term_eq containing the .
void del_deriv() override
Destroys the derivated members (like coloc, cart and radius), when changing the type of colocation po...
virtual void set_anti_legendre_base(Base_spectral &) const
Gives the base using Legendre polynomials, for functions antisymetric with respect to .
virtual int give_place_var(char *) const
Translates a name of a coordinate into its corresponding numerical name.
virtual void affecte_tau(Tensor &, int, const Array< double > &, int &) const
Affects some coefficients to a Tensor.
virtual void set_anti_legendre_base_with_m(Base_spectral &, int m) const
Gives the base using Legendre polynomials, for functions antisymetric with respect to .
Term_eq * normal_spher
Pointer on the Term_eq containing the normal vector to the outer boundary, in orthonormal spherical c...
virtual Val_domain der_normal(const Val_domain &, int) const
Normal derivative with respect to a given surface.
void affecte_tau_one_coef_val_domain(Val_domain &so, int mquant, int cf, int &pos_cf) const
Affects some coefficients to a Val_domain.
virtual Base_spectral mult(const Base_spectral &, const Base_spectral &) const
Method for the multiplication of two Base_spectral.
virtual void set_anti_cheb_base_with_m(Base_spectral &, int m) const
Gives the base using Chebyshev polynomials, for functions antisymetric with respect to .
virtual void set_legendre_base_with_m(Base_spectral &, int m) const
Gives the stnadard base using Legendre polynomials.
int nbr_conditions_val_domain(const Val_domain &so, int mquant, int order) const
Computes number of discretized equations associated with a given tensorial equation in the bulk.
void export_tau_val_domain_boundary(const Val_domain &eq, int mquant, 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 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 xx_to_vars_from_adapted(Val_domain &, const Array< double > &, int &) const
Computes the new boundary of a Domain from a set of values.
virtual void set_cheb_base(Base_spectral &) const
Gives the standard base for Chebyshev polynomials.
virtual Val_domain der_r(const Val_domain &) const
Compute the radial derivative of a scalar field.
virtual void do_radius() const
Computes the generalized radius.
virtual const Term_eq * give_normal(int, int) const
Returns the vector normal to a surface.
virtual void set_cheb_r_base(Base_spectral &) const
Gives the base using odd Chebyshev polynomials$ for the radius.
Term_eq derive_t(const Term_eq &so) const
Computes .
virtual void update_term_eq(Term_eq *) const
Update the value of a field, after the shape of the Domain has been changed by the system.
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 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...
int nbr_conditions_val_domain_boundary(const Val_domain &eq, int mquant) const
Computes number of discretized equations associated with a given equation on a boundary.
virtual Term_eq grad_term_eq(const Term_eq &) const
Gradient of Term_eq.
virtual int nbr_unknowns_from_adapted() const
Gives the number of unknowns coming from the variable shape of the domain.
virtual void do_absol() const
Computes the absolute coordinates.
Val_domain * outer_radius
Pointer on the outer boundary , as a Val_domain.
virtual void do_cart_surr() const
Computes the Cartesian coordinates over the radius.
virtual Term_eq div_r_term_eq(const Term_eq &) const
Division by of a Term_eq.
virtual Val_domain div_r(const Val_domain &) const
Division by .
void update() const
Updates all the quantities that depend on the inner radius (like the normal vectors).
virtual void do_cart() const
Computes the Cartesian coordinates.
Term_eq * rad_term_eq
Pointer on the Term_eq containing the radius.
Val_domain get_outer_radius() const
Returns the outer variable boundary.
virtual void set_cheb_base_with_m(Base_spectral &, int m) const
Gives the standard base using Chebyshev polynomials.
void affecte_tau_val_domain(Val_domain &so, int mquant, const Array< double > &cf, int &pos_cf) const
Affects some coefficients to a Val_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 double val_boundary(int, const Val_domain &, const Index &) const
Computes the value of a field at a boundary.
virtual Val_domain div_sin_theta(const Val_domain &) const
Division by .
virtual void find_other_dom(int, int, int &, int &) const
Gives the informations corresponding the a touching neighboring domain.
virtual Val_domain laplacian(const Val_domain &, int) const
Computes the ordinary flat Laplacian for a scalar field with an harmonic index m.
virtual Val_domain mult_sin_theta(const Val_domain &) const
Multiplication by .
void set_mapping(const Val_domain &so) const
Affects the outer radius.
virtual bool is_in(const Point &xx, double prec=1e-13) const
Check whether a point lies inside Domain.
Term_eq * der_rad_term_eq
Pointer on the Term_eq containing the .
virtual int nbr_unknowns(const Tensor &, int) const
Computes the number of true unknowns of a Tensor, in a given domain.
virtual void update_constante(const Val_domain &, const Scalar &, Scalar &) const
Update the value of a scalar, after the shape of the Domain has been changed by the system.
virtual void affecte_coef(int &, int, bool &) const
The variation of the functions describing the shape of the Domain are affected from the unknowns of t...
Term_eq * outer_radius_term_eq
Pointer on the outer boundary , as a Term_eq.
virtual void set_legendre_r_base(Base_spectral &) const
Gives the base using odd Legendre polynomials$ for the radius.
Abstract class that implements the fonctionnalities common to all the type of domains.
Definition: space.hpp:60
Class that gives the position inside a multi-dimensional Array.
Definition: index.hpp:38
The class Point is used to store the coordinates of a point.
Definition: point.hpp:30
The class Scalar does not really implements scalars in the mathematical sense but rather tensorial co...
Definition: scalar.hpp:67
The Space_polar_adapted class fills the space with one polar nucleus, one polar shell adapted on the ...
virtual void affecte_coef_to_variable_domains(int &, int, Array< int > &) const
The variation of the functions describing the shape of the Domain are affected from the unknowns of t...
virtual void save(FILE *) const
Saving function.
virtual int nbr_unknowns_from_variable_domains() const
Gives the number of unknowns coming from the variable shape of the domain.
virtual void xx_to_vars_variable_domains(System_of_eqs *, const Array< double > &, int &) const
Update the variables of a system, from the variation of the shape of the domains.
void add_eq(System_of_eqs &syst, const char *eq, const char *rac, const char *rac_der, int nused=-1, Array< int > **pused=0x0) const
Adds a bulk equation and two matching conditions.
Space_polar_adapted(int ttyp, const Point &cr, const Dim_array &nbr, const Array< double > &bounds)
virtual void xx_to_ders_variable_domains(const Array< double > &, int &) const
Update the vairable domains from a set of values.
virtual ~Space_polar_adapted()
Destructor
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