KADATH
adapted.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 __ADAPTED_HPP_
21 #define __ADAPTED_HPP_
22 
23 #include "space.hpp"
24 #include "term_eq.hpp"
25 #include "spheric.hpp"
26 #include "metric.hpp"
27 
28 namespace Kadath {
29 
52 
53  protected:
57  const Space& sp ;
60  double outer_radius ;
64  mutable Term_eq* normal_spher ;
68  mutable Term_eq* normal_cart ;
72  mutable Term_eq* rad_term_eq ;
80  mutable Term_eq* dt_rad_term_eq ;
85 
87 
88  public:
99  Domain_shell_inner_adapted (const Space& sp, int num, int ttype, double rin, double rout, const Point& cr, const Dim_array& nbr) ;
110  Domain_shell_inner_adapted (const Space& sp, int num, int ttype, const Val_domain& rin, double rout, const Point& cr, const Dim_array& nbr) ;
118  Domain_shell_inner_adapted (const Space& sp, int num, FILE* fd) ;
119 
120  virtual ~Domain_shell_inner_adapted () ;
121  void del_deriv() override;
122  virtual void save (FILE*) const ;
127 
128  private:
129  virtual void do_absol () const ;
130  virtual void do_radius () const ;
131  virtual void do_cart () const ;
132  virtual void do_cart_surr () const ;
133 
134 
135  protected:
136 
137  virtual void set_cheb_base(Base_spectral&) const ;
138  virtual void set_legendre_base(Base_spectral&) const ;
139 
140  virtual void set_anti_cheb_base(Base_spectral&) const ;
141  virtual void set_anti_legendre_base(Base_spectral&) const ;
142 
143  virtual Tensor change_basis_cart_to_spher (int, const Tensor&) const ;
144  virtual Tensor change_basis_spher_to_cart (int, const Tensor&) const ;
145 
146  virtual void set_cheb_base_r_spher(Base_spectral&) const ;
147  virtual void set_cheb_base_t_spher(Base_spectral&) const ;
148  virtual void set_cheb_base_p_spher(Base_spectral&) const ;
149  virtual void set_legendre_base_r_spher(Base_spectral&) const ;
150  virtual void set_legendre_base_t_spher(Base_spectral&) const ;
151  virtual void set_legendre_base_p_spher(Base_spectral&) const ;
152 
153  virtual void set_cheb_r_base(Base_spectral&) const ;
154  virtual void set_legendre_r_base(Base_spectral&) const ;
155 
156  virtual void do_coloc () ;
157  virtual int give_place_var (char*) const ;
158 
159  virtual int nbr_unknowns_from_adapted() const ;
160  virtual void vars_to_terms() const ;
161  virtual void affecte_coef (int&, int, bool&) const ;
162  virtual void xx_to_vars_from_adapted (Val_domain&, const Array<double>&, int&) const ;
163  virtual void xx_to_ders_from_adapted (const Array<double>&, int&) const ;
164  virtual void update_term_eq (Term_eq*) const ;
165  virtual void update_variable (const Val_domain&, const Scalar&, Scalar&) const ;
166  virtual void update_constante (const Val_domain&, const Scalar&, Scalar&) const ;
167  virtual void update_mapping(const Val_domain&) ;
168 
169  public:
173  void set_mapping(const Val_domain& so) const ;
177  void update() const ;
178 
179  public:
180  virtual Point get_center () const {return center ;} ;
181  virtual bool is_in(const Point&xx, double prec=1e-13) const ;
182  virtual const Point absol_to_num(const Point&) const;
183 
184  virtual const Point absol_to_num_bound(const Point&, int) const;
185 
186 
187  virtual void do_der_abs_from_der_var(const Val_domain *const *const der_var, Val_domain **const der_abs) const ;
188  virtual Base_spectral mult (const Base_spectral&, const Base_spectral&) const ;
189 
190  public:
191 
192  virtual Val_domain mult_cos_phi (const Val_domain&) const ;
193  virtual Val_domain mult_sin_phi (const Val_domain&) const ;
194  virtual Val_domain mult_cos_theta (const Val_domain&) const ;
195  virtual Val_domain mult_sin_theta (const Val_domain&) const ;
196  virtual Val_domain div_sin_theta (const Val_domain&) const ;
197  virtual Val_domain div_cos_theta (const Val_domain&) const ;
198  virtual Val_domain ddp (const Val_domain&) const ;
199  virtual Val_domain der_r (const Val_domain&) const ;
200  virtual Val_domain div_r (const Val_domain&) const ;
201  virtual Val_domain laplacian (const Val_domain&, int) const ;
202  virtual Val_domain laplacian2 (const Val_domain&, int) const ;
203 
204  virtual double val_boundary (int, const Val_domain&, const Index&) const ;
205  virtual void find_other_dom (int, int, int&, int&) const ;
206  virtual Val_domain der_normal (const Val_domain&, int) const ;
207 
208  virtual int nbr_unknowns (const Tensor&, int) const ;
216  int nbr_unknowns_val_domain (const Val_domain& so, int mlim) const ;
217  virtual Array<int> nbr_conditions (const Tensor&, int, int, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
226  int nbr_conditions_val_domain (const Val_domain& so, int mlim, int order) const ;
227  virtual Array<int> nbr_conditions_boundary (const Tensor&, int, int, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
236  int nbr_conditions_val_domain_boundary (const Val_domain& eq, int mlim) const ;
237  virtual void export_tau (const Tensor&, int, int, Array<double>&, int&, const Array<int>&, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
248  void export_tau_val_domain (const Val_domain& eq, int mlim, int order, Array<double>& res, int& pos_res, int ncond) const ;
249  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 ;
260  void export_tau_val_domain_boundary (const Val_domain& eq, int mlim, int bound, Array<double>& res, int& pos_res, int ncond) const ;
261  virtual void affecte_tau (Tensor&, int, const Array<double>&, int&) const ;
270  void affecte_tau_val_domain (Val_domain& so, int mlim, const Array<double>& cf, int& pos_cf) const ;
271  virtual void affecte_tau_one_coef (Tensor&, int, int, int&) const ;
280  void affecte_tau_one_coef_val_domain (Val_domain& so, int mlim, int cc, int& pos_cf) const ;
281 
282  virtual int nbr_points_boundary (int, const Base_spectral&) const ;
283  virtual void do_which_points_boundary (int, const Base_spectral&, Index**, int) const ;
284 
290  Term_eq flat_grad_spher (const Term_eq&) const ;
291 
292  virtual Term_eq partial_spher (const Term_eq&) const ;
293  virtual Term_eq partial_cart (const Term_eq&) const ;
294  virtual Term_eq connection_spher (const Term_eq&) const ;
295  virtual const Term_eq* give_normal(int, int) const ;
296 
302  Term_eq derive_r (const Term_eq& so) const ;
308  Term_eq derive_t (const Term_eq& so) const ;
314  Term_eq derive_p (const Term_eq& so) const ;
318  void do_normal_spher () const ;
322  void do_normal_cart () const ;
323 
324  virtual Term_eq der_normal_term_eq (const Term_eq&, int) const ;
325  virtual Term_eq dr_term_eq (const Term_eq&) const ;
326  virtual Term_eq lap_term_eq (const Term_eq&, int) const ;
327  virtual Term_eq mult_r_term_eq (const Term_eq&) const ;
328  virtual Term_eq integ_volume_term_eq (const Term_eq&) const ;
329  virtual Term_eq integ_term_eq (const Term_eq&, int) const ;
330 
331  virtual Term_eq derive_flat_spher (int, char, const Term_eq&, const Metric*) const ;
332  virtual Term_eq derive_flat_cart (int, char, const Term_eq&, const Metric*) const ;
333 
334  virtual Tensor import (int, int, int, const Array<int>&, Tensor**) const ;
335 
336  virtual double integ_volume (const Val_domain&) const ;
337 
338 public:
339  virtual ostream& print (ostream& o) const ;
340 
341  friend class Space_spheric_adapted ;
342  friend class Space_bin_ns ;
343  friend class Space_bin_bh ;
344  friend class Space_adapted_bh ;
345  friend class Space_bbh ;
346 } ;
347 
368 
369  protected:
373  const Space& sp ;
376  double inner_radius ;
384  mutable Term_eq* normal_cart ;
388  mutable Term_eq* rad_term_eq ;
401 
403 
404  public:
415  Domain_shell_outer_adapted (const Space& sp, int num, int ttype, double rin, double rout, const Point& cr, const Dim_array& nbr) ;
426  Domain_shell_outer_adapted (const Space& sp, int num, int ttype, double rin, const Val_domain& rout, const Point& cr, const Dim_array& nbr) ;
434  Domain_shell_outer_adapted (const Space& sp, int num, FILE* fd) ;
435 
436  virtual ~Domain_shell_outer_adapted () ;
437  void del_deriv() override;
438  virtual void save (FILE*) const ;
443 
444  private:
445  virtual void do_absol () const ;
446  virtual void do_radius () const ;
447  virtual void do_cart () const ;
448  virtual void do_cart_surr () const ;
449 
450 
451  protected:
452 
453  virtual void set_cheb_base(Base_spectral&) const ;
454  virtual void set_legendre_base(Base_spectral&) const ;
455 
456  virtual void set_anti_cheb_base(Base_spectral&) const ;
457  virtual void set_anti_legendre_base(Base_spectral&) const ;
458 
459  virtual void set_cheb_base_r_spher(Base_spectral&) const ;
460  virtual void set_cheb_base_t_spher(Base_spectral&) const ;
461  virtual void set_cheb_base_p_spher(Base_spectral&) const ;
462  virtual void set_legendre_base_r_spher(Base_spectral&) const ;
463  virtual void set_legendre_base_t_spher(Base_spectral&) const ;
464  virtual void set_legendre_base_p_spher(Base_spectral&) const ;
465 
466  virtual void set_cheb_r_base(Base_spectral&) const ;
467  virtual void set_legendre_r_base(Base_spectral&) const ;
468 
469  virtual void do_coloc () ;
470  virtual int give_place_var (char*) const ;
471 
472  virtual int nbr_unknowns_from_adapted() const ;
473  virtual void vars_to_terms() const ;
474  virtual void affecte_coef (int&, int, bool&) const ;
475  virtual void xx_to_vars_from_adapted (Val_domain&, const Array<double>&, int&) const ;
476  virtual void xx_to_ders_from_adapted (const Array<double>&, int&) const ;
477  virtual void update_term_eq (Term_eq*) const ;
478  virtual void update_variable (const Val_domain&, const Scalar&, Scalar&) const ;
479  virtual void update_constante (const Val_domain&, const Scalar&, Scalar&) const ;
480  virtual void update_mapping(const Val_domain&) ;
481  public:
485  void set_mapping(const Val_domain& so) const ;
489  void update() const ;
490 
491  public:
492  virtual Point get_center () const {return center ;} ;
493 
494  virtual bool is_in(const Point&xx, double prec=1e-13) const ;
495  virtual const Point absol_to_num(const Point&) const;
496  virtual void do_der_abs_from_der_var(const Val_domain *const *const der_var, Val_domain **const der_abs) const ;
497  virtual Base_spectral mult (const Base_spectral&, const Base_spectral&) const ;
498 
499  public:
500  virtual Val_domain mult_cos_phi (const Val_domain&) const ;
501  virtual Val_domain mult_sin_phi (const Val_domain&) const ;
502  virtual Val_domain mult_cos_theta (const Val_domain&) const ;
503  virtual Val_domain mult_sin_theta (const Val_domain&) const ;
504  virtual Val_domain div_sin_theta (const Val_domain&) const ;
505  virtual Val_domain div_cos_theta (const Val_domain&) const ;
506  virtual Val_domain laplacian (const Val_domain&, int) const ;
507  virtual Val_domain laplacian2 (const Val_domain&, int) const ;
508 
509  virtual Tensor change_basis_cart_to_spher (int dd, const Tensor&) const ;
510  virtual Tensor change_basis_spher_to_cart (int dd, const Tensor&) const ;
511 
512  virtual Val_domain ddp (const Val_domain&) const ;
513  virtual Val_domain der_r (const Val_domain&) const ;
514  virtual Val_domain div_r (const Val_domain&) const ;
515 
516  virtual double val_boundary (int, const Val_domain&, const Index&) const ;
517  virtual void find_other_dom (int, int, int&, int&) const ;
518  virtual Val_domain der_normal (const Val_domain&, int) const ;
519 
520  virtual int nbr_unknowns (const Tensor&, int) const ;
528  int nbr_unknowns_val_domain (const Val_domain& so, int mlim) const ;
529  virtual Array<int> nbr_conditions (const Tensor&, int, int, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
538  int nbr_conditions_val_domain (const Val_domain& so, int mlim, int order) const ;
539  virtual Array<int> nbr_conditions_boundary (const Tensor&, int, int, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
548  int nbr_conditions_val_domain_boundary (const Val_domain& eq, int mlim) const ;
549  virtual void export_tau (const Tensor&, int, int, Array<double>&, int&, const Array<int>&, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
560  void export_tau_val_domain (const Val_domain& eq, int mlim, int order, Array<double>& res, int& pos_res, int ncond) const ;
561  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 ;
572  void export_tau_val_domain_boundary (const Val_domain& eq, int mlim, int bound, Array<double>& res, int& pos_res, int ncond) const ;
573  virtual void affecte_tau (Tensor&, int, const Array<double>&, int&) const ;
582  void affecte_tau_val_domain (Val_domain& so, int mlim, const Array<double>& cf, int& pos_cf) const ;
583  virtual void affecte_tau_one_coef (Tensor&, int, int, int&) const ;
592  void affecte_tau_one_coef_val_domain (Val_domain& so, int mlim, int cf, int& pos_cf) const ;
593 
594  virtual int nbr_points_boundary (int, const Base_spectral&) const ;
595  virtual void do_which_points_boundary (int, const Base_spectral&, Index**, int) const ;
596 
602  Term_eq flat_grad_spher (const Term_eq&) const ;
603  virtual Term_eq partial_spher (const Term_eq&) const ;
604  virtual Term_eq partial_cart (const Term_eq&) const ;
605  virtual Term_eq connection_spher (const Term_eq&) const ;
606  virtual const Term_eq* give_normal(int, int) const ;
607 
613  Term_eq derive_r (const Term_eq& so) const ;
619  Term_eq derive_t (const Term_eq& so) const ;
625  Term_eq derive_p (const Term_eq& so) const ;
629  void do_normal_spher () const ;
633  void do_normal_cart () const ;
634 
635  virtual Term_eq der_normal_term_eq (const Term_eq&, int) const ;
636  virtual Term_eq dr_term_eq (const Term_eq&) const ;
637  virtual Term_eq lap_term_eq (const Term_eq&, int) const ;
638  virtual Term_eq mult_r_term_eq (const Term_eq&) const ;
639  virtual Term_eq integ_volume_term_eq (const Term_eq&) const ;
640  virtual Term_eq integ_term_eq (const Term_eq&, int) const ;
641  virtual Term_eq derive_flat_spher (int, char, const Term_eq&, const Metric*) const ;
642  virtual Term_eq derive_flat_cart (int, char, const Term_eq&, const Metric*) const ;
643 
644  virtual double integ_volume (const Val_domain&) const ;
645  virtual Tensor import (int, int, int, const Array<int>&, Tensor**) const ;
646 
647 public:
648  virtual ostream& print (ostream& o) const ;
649 
650  friend class Space_spheric_adapted ;
651  friend class Space_bin_ns ;
652  friend class Space_bin_bh ;
653  friend class Space_adapted_bh ;
654  friend class Space_bbh ;
655 } ;
656 
661 class Space_spheric_adapted : public Space {
662  public:
670  Space_spheric_adapted (int ttyp, const Point& cr, const Dim_array& nbr, const Array<double>& bounds) ;
671 
672  Space_spheric_adapted (FILE*) ;
673  virtual ~Space_spheric_adapted() ;
674  virtual void save(FILE*) const ;
675 
676  virtual int nbr_unknowns_from_variable_domains() const ;
677  virtual void affecte_coef_to_variable_domains(int& , int, Array<int>&) const ;
678  virtual void xx_to_ders_variable_domains(const Array<double>&, int&) const ;
679  virtual void xx_to_vars_variable_domains(System_of_eqs*, const Array<double>&, int&) const ;
680 
686  void add_eq_ori (System_of_eqs& syst, const char* eq) ;
696  void add_eq (System_of_eqs& syst, const char* eq, const char* rac, const char* rac_der, int nused=-1, Array<int>** pused=0x0) ;
702  void add_eq_int_inf (System_of_eqs& syst, const char* eq) ;
709  void add_eq_int_volume (System_of_eqs& syst, int nz, const char* eq) ;
710 
711  virtual Array<int> get_indices_matching_non_std(int, int) const ;
712 } ;
713 }
714 #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 .
Definition: adapted.hpp:51
virtual Val_domain div_sin_theta(const Val_domain &) const
Division by .
Val_domain get_inner_radius() const
Returns the inner variable boundary.
Definition: adapted.hpp:126
Term_eq * der_rad_term_eq
Pointer on the Term_eq containing the .
Definition: adapted.hpp:76
virtual void set_legendre_r_base(Base_spectral &) const
Gives the base using odd Legendre polynomials$ for the radius.
Term_eq flat_grad_spher(const Term_eq &) const
Computes the flat gradient of a field, in orthonormal spherical coordinates.
virtual int nbr_unknowns(const Tensor &, int) const
Computes the number of true unknowns of a Tensor, in a given domain.
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...
int nbr_unknowns_val_domain(const Val_domain &so, int mlim) const
Computes the number of true unknowns of a Val_domain.
virtual Base_spectral mult(const Base_spectral &, const Base_spectral &) const
Method for the multiplication of two Base_spectral.
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 void set_legendre_base_t_spher(Base_spectral &) const
Gives the base using Legendre polynomials, for the component of a vector.
Term_eq derive_t(const Term_eq &so) const
Computes .
virtual double val_boundary(int, const Val_domain &, const Index &) const
Computes the value of a field at a boundary.
virtual Tensor change_basis_cart_to_spher(int, const Tensor &) const
Changes the tensorial basis from Cartsian to spherical in a given domain.
void update() const
Updates all the quantities that depend on the inner radius (like the normal vectors).
Term_eq derive_p(const Term_eq &so) const
Computes .
virtual const Point absol_to_num(const Point &) const
Computes the numerical coordinates from the physical ones.
int nbr_conditions_val_domain_boundary(const Val_domain &eq, int mlim) const
Computes number of discretized equations associated with a given equation on a boundary.
virtual int give_place_var(char *) const
Translates a name of a coordinate into its corresponding numerical name.
virtual void update_mapping(const Val_domain &)
Updates the variables parts of the Domain.
virtual Term_eq connection_spher(const Term_eq &) const
Computes the part of the gradient involving the connections, in spherical orthonormal coordinates.
virtual void set_cheb_base_r_spher(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the radial component of a vector.
virtual Term_eq integ_volume_term_eq(const Term_eq &) const
Volume integral of a Term_eq.
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 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_legendre_base(Base_spectral &) const
Gives the standard base for Legendre polynomials.
virtual Val_domain div_cos_theta(const Val_domain &) const
Division by .
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 ddp(const Val_domain &) const
Compute the second derivative with respect to of a scalar field.
virtual int nbr_unknowns_from_adapted() const
Gives the number of unknowns coming from the variable shape of the domain.
Term_eq * rad_term_eq
Pointer on the Term_eq containing the radius.
Definition: adapted.hpp:72
virtual void do_coloc()
Computes the colocation points.
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.
double outer_radius
The outer radius .
Definition: adapted.hpp:60
virtual Term_eq dr_term_eq(const Term_eq &) const
Radial derivative of a Term_eq.
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.
void set_mapping(const Val_domain &so) const
Affects the inner radius.
virtual void affecte_tau(Tensor &, int, const Array< double > &, int &) const
Affects some coefficients to a Tensor.
virtual int nbr_points_boundary(int, const Base_spectral &) const
Computes the number of relevant collocation points on a boundary.
virtual Term_eq der_normal_term_eq(const Term_eq &, int) const
Returns the normal derivative of a Term_eq.
virtual Term_eq lap_term_eq(const Term_eq &, int) const
Returns the flat Laplacian of Term_eq, for a given harmonic.
virtual void set_cheb_r_base(Base_spectral &) const
Gives the base using odd Chebyshev polynomials$ for the radius.
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 void set_cheb_base_t_spher(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
virtual bool is_in(const Point &xx, double prec=1e-13) const
Check whether a point lies inside Domain.
virtual void set_cheb_base(Base_spectral &) const
Gives the standard base for Chebyshev polynomials.
Term_eq * normal_spher
Pointer on the Term_eq containing the normal vector to the inner boundary, in orthonormal spherical c...
Definition: adapted.hpp:64
void do_normal_cart() const
Computes the normal wrt the inner boundary, in Cartesian coordinates.
virtual void affecte_tau_one_coef(Tensor &, int, int, int &) const
Sets at most one coefficient of a Tensor to 1.
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.
void export_tau_val_domain_boundary(const Val_domain &eq, int mlim, 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 do_which_points_boundary(int, const Base_spectral &, Index **, int) const
Lists all the indices corresponding to true collocation points on a boundary.
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 set_legendre_base_r_spher(Base_spectral &) const
Gives the base using Legendre polynomials, for the radial component of a vector.
virtual Val_domain div_r(const Val_domain &) const
Division by .
virtual Point get_center() const
Returns the center.
Definition: adapted.hpp:180
virtual Val_domain der_normal(const Val_domain &, int) const
Normal derivative with respect to a given surface.
virtual Val_domain der_r(const Val_domain &) const
Compute the radial derivative of a scalar field.
virtual void set_anti_legendre_base(Base_spectral &) const
Gives the base using Legendre polynomials, for functions antisymetric with respect to .
Val_domain * inner_radius
Pointer on the inner boundary , as a Val_domain.
Definition: adapted.hpp:58
const Space & sp
The corresponding Space ; required for updating fields whene the mapping changes.
Definition: adapted.hpp:57
Term_eq * dt_rad_term_eq
Pointer on the Term_eq containing the .
Definition: adapted.hpp:80
virtual Term_eq integ_term_eq(const Term_eq &, int) const
Surface integral of a Term_eq.
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 Val_domain mult_sin_theta(const Val_domain &) const
Multiplication by .
virtual Val_domain laplacian(const Val_domain &, int) const
Computes the ordinary flat Laplacian for a scalar field with an harmonic index m.
Term_eq * normal_cart
Pointer on the Term_eq containing the normal vector to the inner boundary, in Cartesian coordinates.
Definition: adapted.hpp:68
virtual Val_domain mult_cos_phi(const Val_domain &) const
Multiplication by .
void del_deriv() override
Destroys the derivated members (like coloc, cart and radius), when changing the type of colocation po...
void do_normal_spher() const
Computes the normal wrt the inner boundary, in orthonormal spherical coordinates.
Term_eq * inner_radius_term_eq
Pointer on the inner boundary , as a Term_eq.
Definition: adapted.hpp:59
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 const Term_eq * give_normal(int, int) const
Returns the vector normal to a surface.
virtual void set_anti_cheb_base(Base_spectral &) const
Gives the base using Chebyshev polynomials, for functions antisymetric with respect to .
int nbr_conditions_val_domain(const Val_domain &so, int mlim, int order) const
Computes number of discretized equations associated with a given tensorial equation in the bulk.
virtual double integ_volume(const Val_domain &) const
Volume integral.
virtual void vars_to_terms() const
The Term_eq describing the variable shape of the Domain are updated.
virtual Val_domain mult_sin_phi(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.
virtual void set_legendre_base_p_spher(Base_spectral &) const
Gives the base using Legendre polynomials, for the component of a vector.
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...
Domain_shell_inner_adapted(const Space &sp, int num, int ttype, double rin, double rout, const Point &cr, const Dim_array &nbr)
Constructor :
void affecte_tau_val_domain(Val_domain &so, int mlim, const Array< double > &cf, int &pos_cf) const
Affects some coefficients to a Val_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 ostream & print(ostream &o) const
Delegate function to virtualize the << operator.
virtual Tensor change_basis_spher_to_cart(int, const Tensor &) const
Changes the tensorial basis from spherical to Cartesian in a given domain.
virtual void do_cart_surr() const
Computes the Cartesian coordinates over the radius.
virtual void do_cart() const
Computes the Cartesian coordinates.
virtual Term_eq mult_r_term_eq(const Term_eq &) const
Multiplication by of a Term_eq.
Point center
Absolute coordinates of the center.
Definition: adapted.hpp:86
virtual Val_domain mult_cos_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...
void affecte_tau_one_coef_val_domain(Val_domain &so, int mlim, int cc, int &pos_cf) const
Sets at most one coefficient of a Val_domain to 1.
virtual void save(FILE *) const
Saving function.
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 do_radius() const
Computes the generalized radius.
void export_tau_val_domain(const Val_domain &eq, int mlim, int order, Array< double > &res, int &pos_res, int ncond) const
Exports a residual equation in the bulk.
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 .
Term_eq * dp_rad_term_eq
Pointer on the Term_eq containing the .
Definition: adapted.hpp:84
virtual void find_other_dom(int, int, int &, int &) const
Gives the informations corresponding the a touching neighboring domain.
virtual void do_absol() const
Computes the absolute coordinates.
Class for a spherical-like domain, having a symmetry with respect to the plane .
Definition: adapted.hpp:367
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 der_normal_term_eq(const Term_eq &, int) const
Returns the normal derivative of a Term_eq.
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 mult_sin_phi(const Val_domain &) const
Multiplication by .
void set_mapping(const Val_domain &so) const
Affects the outer radius.
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.
Domain_shell_outer_adapted(const Space &sp, int num, int ttype, double rin, double rout, const Point &cr, const Dim_array &nbr)
Constructor :
virtual Val_domain div_cos_theta(const Val_domain &) const
Division by .
virtual Term_eq integ_volume_term_eq(const Term_eq &) const
Volume integral of a Term_eq.
Val_domain get_outer_radius() const
Returns the outer variable boundary.
Definition: adapted.hpp:442
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 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_absol() const
Computes the absolute coordinates.
virtual int nbr_unknowns(const Tensor &, int) const
Computes the number of true unknowns of a Tensor, in a given domain.
virtual void do_cart_surr() const
Computes the Cartesian coordinates over the radius.
virtual Term_eq dr_term_eq(const Term_eq &) const
Radial derivative of a Term_eq.
virtual void set_cheb_base(Base_spectral &) const
Gives the standard base for Chebyshev polynomials.
Term_eq * outer_radius_term_eq
Pointer on the outer boundary , as a Term_eq.
Definition: adapted.hpp:375
virtual void do_coloc()
Computes the colocation points.
virtual Base_spectral mult(const Base_spectral &, const Base_spectral &) const
Method for the multiplication of two Base_spectral.
virtual void affecte_tau(Tensor &, int, const Array< double > &, int &) const
Affects some coefficients to a Tensor.
virtual void set_legendre_base_r_spher(Base_spectral &) const
Gives the base using Legendre polynomials, for the radial component of a vector.
virtual double val_boundary(int, const Val_domain &, const Index &) const
Computes the value of a field at a boundary.
Term_eq * dp_rad_term_eq
Pointer on the Term_eq containing the .
Definition: adapted.hpp:400
void affecte_tau_one_coef_val_domain(Val_domain &so, int mlim, int cf, int &pos_cf) const
Affects some coefficients to a Val_domain.
virtual Val_domain mult_cos_phi(const Val_domain &) const
Multiplication by .
Point center
Absolute coordinates of the center.
Definition: adapted.hpp:402
virtual void do_radius() const
Computes the generalized radius.
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 int give_place_var(char *) const
Translates a name of a coordinate into its corresponding numerical name.
virtual void save(FILE *) const
Saving function.
Term_eq flat_grad_spher(const Term_eq &) const
Computes the flat gradient of a field, in orthonormal spherical coordinates.
void do_normal_cart() const
Computes the normal wrt the inner boundary, in Cartesian 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...
virtual void set_cheb_base_p_spher(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
Term_eq derive_t(const Term_eq &so) const
Computes .
virtual void set_cheb_r_base(Base_spectral &) const
Gives the base using odd Chebyshev polynomials$ for the radius.
void update() const
Updates all the quantities that depend on the inner radius (like the normal vectors).
virtual void set_legendre_base_t_spher(Base_spectral &) const
Gives the base using Legendre polynomials, for the component of a vector.
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.
Term_eq derive_r(const Term_eq &so) const
Computes .
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.
Term_eq * normal_spher
Pointer on the Term_eq containing the normal vector to the outer boundary, in orthonormal spherical c...
Definition: adapted.hpp:380
void del_deriv() override
Destroys the derivated members (like coloc, cart and radius), when changing the type of colocation po...
virtual void set_legendre_base(Base_spectral &) const
Gives the standard base for Legendre polynomials.
void export_tau_val_domain(const Val_domain &eq, int mlim, int order, Array< double > &res, int &pos_res, int ncond) const
Exports a residual equation in the bulk.
virtual Tensor change_basis_cart_to_spher(int dd, const Tensor &) const
Changes the tensorial basis from Cartsian to spherical in a given domain.
Term_eq * normal_cart
Pointer on the Term_eq containing the normal vector to the outer boundary, in Cartesian coordinates.
Definition: adapted.hpp:384
virtual void find_other_dom(int, int, int &, int &) const
Gives the informations corresponding the a touching neighboring domain.
virtual bool is_in(const Point &xx, double prec=1e-13) const
Check whether a point lies inside Domain.
virtual int nbr_unknowns_from_adapted() const
Gives the number of unknowns coming from the variable shape of the domain.
virtual void set_cheb_base_t_spher(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
virtual Term_eq connection_spher(const Term_eq &) const
Computes the part of the gradient involving the connections, in spherical orthonormal coordinates.
virtual Val_domain div_sin_theta(const Val_domain &) const
Division by .
Term_eq * der_rad_term_eq
Pointer on the Term_eq containing the .
Definition: adapted.hpp:392
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 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 double integ_volume(const Val_domain &) const
Volume integral.
Val_domain * outer_radius
Pointer on the outer boundary , as a Val_domain.
Definition: adapted.hpp:374
virtual Val_domain mult_sin_theta(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 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 Val_domain der_r(const Val_domain &) const
Compute the radial derivative of a scalar field.
Term_eq * rad_term_eq
Pointer on the Term_eq containing the radius.
Definition: adapted.hpp:388
virtual void affecte_tau_one_coef(Tensor &, int, int, int &) const
Sets at most one coefficient of a Tensor to 1.
virtual Term_eq integ_term_eq(const Term_eq &, int) const
Surface integral of a Term_eq.
int nbr_conditions_val_domain_boundary(const Val_domain &eq, int mlim) const
Computes number of discretized equations associated with a given equation on a boundary.
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.
int nbr_conditions_val_domain(const Val_domain &so, int mlim, int order) const
Computes number of discretized equations associated with a given tensorial equation 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...
virtual Val_domain laplacian(const Val_domain &, int) const
Computes the ordinary flat Laplacian for a scalar field with an harmonic index m.
virtual const Point absol_to_num(const Point &) const
Computes the numerical coordinates from the physical ones.
double inner_radius
The inner radius .
Definition: adapted.hpp:376
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 xx_to_ders_from_adapted(const Array< double > &, int &) const
Affects the derivative part of variable a Domain from a set of values.
virtual Tensor change_basis_spher_to_cart(int dd, const Tensor &) const
Changes the tensorial basis from spherical to Cartesian in a given domain.
void affecte_tau_val_domain(Val_domain &so, int mlim, const Array< double > &cf, int &pos_cf) const
Affects some coefficients to a Val_domain.
void do_normal_spher() const
Computes the normal wrt the inner boundary, in orthonormal spherical coordinates.
virtual ostream & print(ostream &o) const
Delegate function to virtualize the << operator.
virtual void set_legendre_r_base(Base_spectral &) const
Gives the base using odd Legendre polynomials$ for the radius.
virtual void update_mapping(const Val_domain &)
Updates the variables parts of the Domain.
virtual void set_anti_legendre_base(Base_spectral &) const
Gives the base using Legendre polynomials, for functions antisymetric with respect to .
Term_eq * dt_rad_term_eq
Pointer on the Term_eq containing the .
Definition: adapted.hpp:396
virtual const Term_eq * give_normal(int, int) const
Returns the vector normal to a surface.
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 vars_to_terms() const
The Term_eq describing the variable shape of the Domain are updated.
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 Val_domain ddp(const Val_domain &) const
Compute the second derivative with respect to of a scalar field.
virtual void set_anti_cheb_base(Base_spectral &) const
Gives the base using Chebyshev polynomials, for functions antisymetric with respect to .
virtual Val_domain der_normal(const Val_domain &, int) const
Normal derivative with respect to a given surface.
int nbr_unknowns_val_domain(const Val_domain &so, int mlim) const
Computes the number of true unknowns of a Val_domain.
virtual Point get_center() const
Returns the center.
Definition: adapted.hpp:492
void export_tau_val_domain_boundary(const Val_domain &eq, int mlim, 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 do_which_points_boundary(int, const Base_spectral &, Index **, int) const
Lists all the indices corresponding to true collocation points on a boundary.
virtual Term_eq mult_r_term_eq(const Term_eq &) const
Multiplication by of a Term_eq.
virtual Val_domain div_r(const Val_domain &) const
Division by .
virtual Val_domain mult_cos_theta(const Val_domain &) const
Multiplication by .
const Space & sp
The corresponding Space ; required for updating fields whene the mapping changes.
Definition: adapted.hpp:373
virtual int nbr_points_boundary(int, const Base_spectral &) const
Computes the number of relevant collocation points on a boundary.
virtual void do_cart() const
Computes the Cartesian coordinates.
Term_eq derive_p(const Term_eq &so) const
Computes .
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
Purely abstract class for metric handling.
Definition: metric.hpp:39
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_spheric_adapted class fills the space with one shell adapted on the inside,...
Definition: adapted_bh.hpp:36
Spacetime intended for binary black hole configurations in full general relativity (see constructor f...
Definition: bbh.hpp:33
Spacetime intended for binary black hole configurations (see constructor for details about the domain...
Definition: bin_bh.hpp:35
Spacetime intended for binary neutron stars configurations (see constructor for details about the dom...
Definition: bin_ns.hpp:35
The Space_spheric_adapted class fills the space with one nucleus, one shell adapted on the outside,...
Definition: adapted.hpp:661
void add_eq_ori(System_of_eqs &syst, const char *eq)
Adds an equation being the value of some field at the origin.
void add_eq(System_of_eqs &syst, const char *eq, const char *rac, const char *rac_der, int nused=-1, Array< int > **pused=0x0)
Adds a bulk equation and two matching conditions.
virtual void xx_to_ders_variable_domains(const Array< double > &, int &) const
Update the vairable domains from a set of values.
virtual ~Space_spheric_adapted()
Destructor
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...
void add_eq_int_inf(System_of_eqs &syst, const char *eq)
Adds an equation being a surface integral at infinity.
virtual Array< int > get_indices_matching_non_std(int, int) const
Gives the number of the other domains, touching a given boundary.
virtual void save(FILE *) const
Saving function.
Space_spheric_adapted(int ttyp, const Point &cr, const Dim_array &nbr, const Array< double > &bounds)
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_int_volume(System_of_eqs &syst, int nz, const char *eq)
Adds an equation being a volume integral in the domains below a given number.
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