KADATH
spheric_periodic.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_PERIODIC_HPP_
21 #define __SPHERIC_PERIODIC_HPP_
22 
23 #include "space.hpp"
24 
25 #define TO_PI 0
26 #define TO_PI_OVER_2 1
27 namespace Kadath {
28 
29 
51 
52  private:
53  double alpha ;
54  double ome ;
60  int type_time ;
61  double maxt ;
62 
63  public:
73  Domain_spheric_periodic_nucleus (int num, int ttype, int typet, double radius, double ome, const Dim_array& nbr) ;
80  Domain_spheric_periodic_nucleus (int num, FILE* ff) ;
81 
83  virtual void save (FILE*) const ;
84 
85  private:
86  virtual void do_absol () const ;
87 
88  private:
89  virtual void set_cheb_base(Base_spectral&) const ;
90  virtual void set_legendre_base(Base_spectral&) const ;
91  virtual void set_anti_cheb_base(Base_spectral&) const ;
92  virtual void set_anti_legendre_base(Base_spectral&) const ;
93  virtual void set_cheb_base_odd(Base_spectral&) const ;
94  virtual void set_legendre_base_odd(Base_spectral&) const ;
95  virtual void do_coloc () ;
96  virtual void do_radius () const ;
97  public:
98  virtual bool is_in(const Point&xx, double prec=1e-13) const ;
99  virtual const Point absol_to_num(const Point&) const;
100  virtual void do_der_abs_from_der_var(const Val_domain *const *const der_var, Val_domain **const der_abs) const ;
101  virtual Base_spectral mult (const Base_spectral&, const Base_spectral&) const ;
102 
103  public:
104  virtual Val_domain div_x (const Val_domain&) const ;
105  virtual Val_domain mult_r (const Val_domain&) const ;
106  virtual Val_domain div_r (const Val_domain&) const ;
107  virtual Val_domain laplacian (const Val_domain&, int) const ;
108  virtual Val_domain der_r (const Val_domain&) const ;
109  virtual Val_domain srdr (const Val_domain&) const ;
110  virtual Val_domain ddtime (const Val_domain&) const ;
111  virtual Val_domain dtime (const Val_domain&) const ;
112 
113  virtual double val_boundary (int, const Val_domain&, const Index&) const ;
114  virtual void find_other_dom (int, int, int&, int&) const ;
115  virtual Val_domain der_normal (const Val_domain&, int) const ;
116 
117  virtual int nbr_unknowns (const Tensor&, int) const ;
124  int nbr_unknowns_val_domain (const Val_domain& so) const ;
125  virtual Array<int> nbr_conditions (const Tensor&, int, int, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
133  int nbr_conditions_val_domain (const Val_domain& so, int order) const ;
134  virtual Array<int> nbr_conditions_boundary (const Tensor&, int, int, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
135 
143  int nbr_conditions_val_domain_boundary (const Val_domain& eq) const ;
144  virtual void export_tau (const Tensor&, int, int, Array<double>&, int&, const Array<int>&,int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
154  void export_tau_val_domain (const Val_domain& eq, int order, Array<double>& res, int& pos_res, int ncond) const ;
155  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 ;
165  void export_tau_val_domain_boundary (const Val_domain& eq, int bound, Array<double>& res, int& pos_res, int ncond) const ;
166  virtual void affecte_tau (Tensor&, int, const Array<double>&, int&) const ;
174  void affecte_tau_val_domain (Val_domain& so, const Array<double>& cf, int& pos_cf) const ;
175  virtual void affecte_tau_one_coef (Tensor&, int, int, int&) const ;
183  void affecte_tau_one_coef_val_domain (Val_domain& so, int cc, int& pos_cf) const ;
184 
185 public:
186  virtual ostream& print (ostream& o) const ;
187 } ;
188 
210 
211  private:
212  double alpha ;
213  double beta ;
214  double ome ;
220  int type_time ;
221  double maxt ;
222 
223  public:
234  Domain_spheric_periodic_shell (int num, int ttype, int typet, double r_int, double r_ext, double ome, const Dim_array& nbr) ;
241  Domain_spheric_periodic_shell (int num, FILE* ff) ;
242 
243  virtual ~Domain_spheric_periodic_shell() ;
244  virtual void save (FILE*) const ;
245 
246  private:
247  virtual void do_absol () const ;
248 
249  private :
250  virtual void set_cheb_base(Base_spectral&) const ;
251  virtual void set_legendre_base(Base_spectral&) const ;
252  virtual void set_anti_cheb_base(Base_spectral&) const ;
253  virtual void set_anti_legendre_base(Base_spectral&) const ;
254  virtual void set_cheb_base_odd(Base_spectral&) const ;
255  virtual void set_legendre_base_odd(Base_spectral&) const ;
256  virtual void do_coloc () ;
257  virtual void do_radius () const ;
258  public:
259  virtual bool is_in(const Point& xx, double prec=1e-13) const ;
260  virtual const Point absol_to_num(const Point&) const;
261  virtual void do_der_abs_from_der_var(const Val_domain *const *const der_var, Val_domain **const der_abs) const ;
262 
263  virtual Base_spectral mult (const Base_spectral&, const Base_spectral&) const ;
264 
265  public:
266  virtual Val_domain mult_r (const Val_domain&) const ;
267  virtual Val_domain div_r (const Val_domain&) const ;
268  virtual Val_domain laplacian (const Val_domain&, int) const ;
269  virtual Val_domain der_r (const Val_domain&) const ;
270  virtual Val_domain ddtime (const Val_domain&) const ;
271  virtual Val_domain dtime (const Val_domain&) const ;
272  virtual Val_domain div_xm1 (const Val_domain&) const ;
273 
274  virtual void find_other_dom (int, int, int&, int&) const ;
275  virtual Val_domain der_normal (const Val_domain&, int) const ;
276  virtual double val_boundary (int, const Val_domain&, const Index&) const ;
277 
278  virtual int nbr_unknowns (const Tensor&, int) const ;
285  int nbr_unknowns_val_domain (const Val_domain& so) const ;
286  virtual Array<int> nbr_conditions (const Tensor&, int, int, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
294  int nbr_conditions_val_domain (const Val_domain& so, int order) const ;
295  virtual Array<int> nbr_conditions_boundary (const Tensor&, int, int, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
303  int nbr_conditions_val_domain_boundary (const Val_domain& eq) const ;
304  virtual void export_tau (const Tensor&, int, int, Array<double>&, int&, const Array<int>&,int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
314  void export_tau_val_domain (const Val_domain& eq, int order, Array<double>& res, int& pos_res, int ncond) const ;
315  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 ;
325  void export_tau_val_domain_boundary (const Val_domain& eq, int bound, Array<double>& res, int& pos_res, int ncond) const ;
326  virtual void affecte_tau (Tensor&, int, const Array<double>&, int&) const ;
334  void affecte_tau_val_domain (Val_domain& so, const Array<double>& cf, int& pos_cf) const ;
335  virtual void affecte_tau_one_coef (Tensor&, int, int, int&) const ;
343  void affecte_tau_one_coef_val_domain (Val_domain& so, int cc, int& pos_cf) const ;
344 
354  double give_harmonique (const Val_domain& so, int nn, double phase, int dim, double fact) const ;
363  Val_domain fitwaves (const Val_domain& so, const Array<double>& phases, int dim, double fact) const ;
372  Term_eq fitwaves (const Term_eq& so, const Array<double>& phases, int dim, double fact) const ;
373 
380  Val_domain fitschwarz (const Val_domain&, int) const ;
387  Term_eq fitschwarz (const Term_eq&, int) const ;
388 
400  double give_harmonique_nonflat (const Val_domain& so, double mass, int nn, double phase, int dim, double fact) const ;
411  Term_eq fitwaves_nonflat (const Term_eq& so, const Term_eq& field, const Array<double>& phases, int dim, double factor) const;
412 
413 
414 public:
415  virtual ostream& print (ostream& o) const ;
416 } ;
417 
439 
440  private:
441  double alpha ;
442  double ome ;
448  int type_time ;
449  double maxt ;
450 
451 
452  public:
462  Domain_spheric_periodic_compact (int num, int ttype, int typet, double r_int, double ome, const Dim_array& nbr) ;
469  Domain_spheric_periodic_compact (int num, FILE* ff) ;
470 
472  virtual void save(FILE*) const ;
473 
474  private:
475  virtual void set_cheb_base(Base_spectral&) const ;
476  virtual void set_legendre_base(Base_spectral&) const ;
477  virtual void set_anti_cheb_base(Base_spectral&) const ;
478  virtual void set_anti_legendre_base(Base_spectral&) const ;
479  virtual void set_cheb_base_odd(Base_spectral&) const ;
480  virtual void set_legendre_base_odd(Base_spectral&) const ;
481 
482  virtual void do_coloc() ;
483  virtual void do_absol () const ;
484  virtual void do_radius () const ;
485 
486  public:
490  double get_alpha() const {return alpha ;} ;
491 
492 
493  virtual bool is_in(const Point& xx, double prec=1e-13) const ;
494  virtual const Point absol_to_num(const Point&) const;
495  virtual void do_der_abs_from_der_var(const Val_domain *const *const der_var, Val_domain **const der_abs) const ;
496 
497  virtual Base_spectral mult (const Base_spectral&, const Base_spectral&) const ;
498 
499  public:
500 
501  virtual Val_domain div_xm1 (const Val_domain&) const ;
502  virtual Val_domain mult_xm1 (const Val_domain&) const ;
503  virtual Val_domain mult_r (const Val_domain&) const ;
504  virtual Val_domain div_r (const Val_domain&) const ;
505  virtual Val_domain laplacian (const Val_domain&, int) const ;
506  virtual Val_domain der_r (const Val_domain&) const ;
507  virtual Val_domain ddtime (const Val_domain&) const ;
508  virtual Val_domain dtime (const Val_domain&) const ;
509 
510  virtual void set_val_inf (Val_domain& so, double xx) const ;
511 
512  virtual void find_other_dom (int, int, int&, int&) const ;
513  virtual Val_domain der_normal (const Val_domain&, int) const ;
514  virtual double val_boundary (int, const Val_domain&, const Index&) const ;
515 
516  virtual int nbr_unknowns (const Tensor&, int) const ;
523  int nbr_unknowns_val_domain (const Val_domain& so) const ;
524  virtual Array<int> nbr_conditions (const Tensor&, int, int, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
532  int nbr_conditions_val_domain (const Val_domain& so, int order) const ;
533  virtual Array<int> nbr_conditions_boundary (const Tensor&, int, int, int n_cmp=-1, Array<int>** p_cmp=0x0) const ;
541  int nbr_conditions_val_domain_boundary (const Val_domain& eq) const ;
542  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 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 ;
563  void export_tau_val_domain_boundary (const Val_domain& eq, int bound, Array<double>& res, int& pos_res, int ncond) const ;
564  virtual void affecte_tau (Tensor&, int, const Array<double>&, int&) const ;
572  void affecte_tau_val_domain (Val_domain& so, const Array<double>& cf, int& pos_cf) const ;
573  virtual void affecte_tau_one_coef (Tensor&, int, int, int&) const ;
581  void affecte_tau_one_coef_val_domain (Val_domain& so, int cc, int& pos_cf) const ;
582 
583 public:
584  virtual ostream& print (ostream& o) const ;
585 } ;
586 
593 
594  protected:
595  double omega ;
596 
597  public:
606  Space_spheric_periodic (int ttype, int typet, const Dim_array& nbr, const Array<double>& bounds, double ome) ;
607  Space_spheric_periodic (FILE*) ;
608  virtual ~Space_spheric_periodic() ;
609  virtual void save(FILE*) const ;
610 
618  void add_outer_bc (System_of_eqs& syst, const char* eq, int nused=-1, Array<int>** pused=0x0) ;
626  void add_eq (System_of_eqs& syst, const char* eq, int nused=-1, Array<int>** pused=0x0) ;
634  void add_eq_full (System_of_eqs& syst, const char* eq, int nused=-1, Array<int>** pused=0x0) ;
642  void add_matching (System_of_eqs& syst, const char* eq, int nused=-1, Array<int>** pused=0x0) ;
652  void add_eq (System_of_eqs& syst, const char* eq, const char* rac, const char* rac_der, int nused=-1, Array<int>** pused=0x0) ;
662  void add_eq_inverted (System_of_eqs& syst, const char* eq, const char* rac, const char* rac_der, int nused=-1, Array<int>** pused=0x0) ;
668  void add_eq_ori (System_of_eqs& syst, const char* eq) ;
669 } ;
670 }
671 #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 2-dimensional compactified spherical domain and a symetry with respect to the plane .
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 affecte_tau_one_coef(Tensor &, int, int, int &) const
Sets at most one coefficient of a Tensor to 1.
virtual void set_anti_cheb_base(Base_spectral &) const
Gives the base using Chebyshev polynomials, for functions antisymetric with respect to .
virtual void set_anti_legendre_base(Base_spectral &) const
Gives the base using Legendre polynomials, for functions antisymetric with respect to .
virtual void affecte_tau(Tensor &, int, const Array< double > &, int &) const
Affects some coefficients to a Tensor.
virtual void do_radius() const
Computes the generalized radius.
virtual void do_absol() const
Computes the absolute coordinates.
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 save(FILE *) const
Saving function.
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 Val_domain div_r(const Val_domain &) const
Division by .
double ome
Relates the numerical time to the physical one.
virtual void find_other_dom(int, int, int &, int &) const
Gives the informations corresponding the a touching neighboring domain.
virtual Val_domain dtime(const Val_domain &) const
Computes the time derivative of a field.
virtual Base_spectral mult(const Base_spectral &, const Base_spectral &) const
Method for the multiplication of two Base_spectral.
Domain_spheric_periodic_compact(int num, int ttype, int typet, double r_int, double ome, const Dim_array &nbr)
Standard constructor :
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_conditions_val_domain(const Val_domain &so, int order) const
Computes number of discretized equations associated with a given tensorial equation in the bulk.
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(Base_spectral &) const
Gives the standard base for Legendre polynomials.
virtual void do_coloc()
Computes the colocation points.
int type_time
Gives the type of time periodicity.
virtual Val_domain ddtime(const Val_domain &) const
Computes the second time derivative of a field.
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 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 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 const Point absol_to_num(const Point &) const
Computes the numerical coordinates from the physical ones.
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 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 set_cheb_base_odd(Base_spectral &) const
Gives the base using odd Chebyshev polynomials$.
virtual Val_domain mult_xm1(const Val_domain &) const
Multiplication by .
virtual Val_domain der_normal(const Val_domain &, int) const
Normal derivative with respect to a given surface.
virtual Val_domain div_xm1(const Val_domain &) const
Division by .
double alpha
Relates the numerical radius to the physical one.
virtual ostream & print(ostream &o) const
Delegate function to virtualize the << operator.
virtual void set_legendre_base_odd(Base_spectral &) const
Gives the base using odd Legendre polynomials$.
virtual void set_cheb_base(Base_spectral &) const
Gives the standard base for Chebyshev polynomials.
int nbr_unknowns_val_domain(const Val_domain &so) const
Computes the number of true unknowns of a Val_domain.
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.
double maxt
Upper bound of which is or .
virtual bool is_in(const Point &xx, double prec=1e-13) const
Check whether a point lies inside Domain.
virtual Val_domain der_r(const Val_domain &) const
Compute the radial derivative of a scalar field.
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 double val_boundary(int, const Val_domain &, const Index &) const
Computes the value of a field at a boundary.
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_r(const Val_domain &) const
Multiplication by .
Class for a 2-dimensional spherical domain containing the origin and a symetry with respect to the pl...
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 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.
double maxt
Upper bound of which is or .
virtual void affecte_tau_one_coef(Tensor &, int, int, int &) const
Sets at most one coefficient of a Tensor to 1.
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 void do_coloc()
Computes the colocation points.
virtual void set_cheb_base(Base_spectral &) const
Gives the standard base for Chebyshev polynomials.
virtual Val_domain srdr(const Val_domain &) const
Compute the of a scalar field .
virtual Val_domain div_x(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.
double ome
Relates the numerical time to the physical one.
virtual ostream & print(ostream &o) const
Delegate function to virtualize the << operator.
double alpha
Relates the numerical radius to the physical one.
virtual void set_legendre_base_odd(Base_spectral &) const
Gives the base using odd Legendre polynomials$.
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_r(const Val_domain &) const
Compute the radial derivative of a scalar field.
virtual Val_domain mult_r(const Val_domain &) const
Multiplication by .
Domain_spheric_periodic_nucleus(int num, int ttype, int typet, double radius, double ome, const Dim_array &nbr)
Standard constructor :
virtual void set_anti_legendre_base(Base_spectral &) const
Gives the base using Legendre polynomials, for functions antisymetric with respect to .
virtual Val_domain ddtime(const Val_domain &) const
Computes the second time derivative of a field.
virtual Val_domain laplacian(const Val_domain &, int) const
Computes the ordinary flat Laplacian for a scalar field with an harmonic index m.
virtual Base_spectral mult(const Base_spectral &, const Base_spectral &) const
Method for the multiplication of two Base_spectral.
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_anti_cheb_base(Base_spectral &) const
Gives the base using Chebyshev polynomials, for functions antisymetric with respect to .
virtual const Point absol_to_num(const Point &) const
Computes the numerical coordinates from the physical ones.
virtual Val_domain div_r(const Val_domain &) const
Division by .
int nbr_conditions_val_domain(const Val_domain &so, 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...
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 void find_other_dom(int, int, int &, int &) const
Gives the informations corresponding the a touching neighboring domain.
int type_time
Gives the type of time periodicity.
virtual Val_domain dtime(const Val_domain &) const
Computes the time derivative of a field.
virtual Val_domain der_normal(const Val_domain &, int) const
Normal derivative with respect to a given surface.
virtual void set_legendre_base(Base_spectral &) const
Gives the standard base for Legendre polynomials.
virtual void set_cheb_base_odd(Base_spectral &) const
Gives the base using odd Chebyshev polynomials$.
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 int nbr_unknowns(const Tensor &, int) const
Computes the number of true unknowns of a Tensor, in a given domain.
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_absol() const
Computes the absolute coordinates.
virtual void do_radius() const
Computes the generalized radius.
virtual void save(FILE *) const
Saving function.
virtual void affecte_tau(Tensor &, int, const Array< double > &, int &) const
Affects some coefficients to a Tensor.
virtual bool is_in(const Point &xx, double prec=1e-13) const
Check whether a point lies inside Domain.
int nbr_conditions_val_domain_boundary(const Val_domain &eq) const
Computes number of discretized equations associated with a given equation on a boundary.
Class for a 2-dimensional spherical domain bounded between two fixed radii and a symetry with respect...
double maxt
Upper bound of which is or .
virtual void set_anti_cheb_base(Base_spectral &) const
Gives the base using Chebyshev polynomials, for functions antisymetric with respect to .
virtual void do_absol() const
Computes the absolute coordinates.
virtual void find_other_dom(int, int, int &, int &) const
Gives the informations corresponding the a touching neighboring domain.
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 Base_spectral mult(const Base_spectral &, const Base_spectral &) const
Method for the multiplication of two Base_spectral.
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_legendre_base(Base_spectral &) const
Gives the base using Legendre polynomials, for functions antisymetric with respect to .
virtual void save(FILE *) const
Saving function.
virtual void set_cheb_base(Base_spectral &) const
Gives the standard base for Chebyshev polynomials.
virtual Val_domain div_xm1(const Val_domain &) const
Division by .
virtual void do_radius() const
Computes the generalized radius.
virtual void affecte_tau(Tensor &, int, const Array< double > &, int &) const
Affects some coefficients to a Tensor.
Val_domain fitschwarz(const Val_domain &, int) const
Fit some field with a decay (Val_domain version).
Definition: fitschwarz.cpp:29
virtual void set_legendre_base_odd(Base_spectral &) const
Gives the base using odd Legendre polynomials$.
double give_harmonique_nonflat(const Val_domain &so, double mass, int nn, double phase, int dim, double fact) const
Computes the coefficient of a given outgoing wave harmonique inside some field, at the outer boundary...
virtual Val_domain ddtime(const Val_domain &) const
Computes the second time derivative of a field.
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 dtime(const Val_domain &) const
Computes the time derivative of a field.
double alpha
Relates the numerical to the physical radii.
virtual void do_coloc()
Computes the colocation points.
Val_domain fitwaves(const Val_domain &so, const Array< double > &phases, int dim, double fact) const
Fit some field with homogenous solutions of the wave equation (Val_domain version).
Definition: fitwaves.cpp:66
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.
void affecte_tau_val_domain(Val_domain &so, const Array< double > &cf, int &pos_cf) const
Affects some coefficients to a Val_domain.
int nbr_conditions_val_domain_boundary(const Val_domain &eq) const
Computes number of discretized equations associated with a given equation on a boundary.
Term_eq fitwaves_nonflat(const Term_eq &so, const Term_eq &field, const Array< double > &phases, int dim, double factor) const
Fit some field with homogenous solutions of the wave equation (Term_eq version).
virtual void affecte_tau_one_coef(Tensor &, int, int, int &) const
Sets at most one coefficient of a Tensor to 1.
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 ...
double ome
Relates the numerical time to the physical one.
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 der_r(const Val_domain &) const
Compute the radial derivative of a scalar field.
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 void set_legendre_base(Base_spectral &) const
Gives the standard base for Legendre polynomials.
int nbr_unknowns_val_domain(const Val_domain &so) const
Computes the number of true unknowns of a Val_domain.
double give_harmonique(const Val_domain &so, int nn, double phase, int dim, double fact) const
Computes the coefficient of a given outgoing wave harmonique inside some field, at the outer boundary...
Definition: fitwaves.cpp:30
virtual int nbr_unknowns(const Tensor &, int) const
Computes the number of true unknowns of a Tensor, in a given domain.
double beta
Relates the numerical to the physical radii.
virtual bool is_in(const Point &xx, double prec=1e-13) const
Check whether a point lies inside Domain.
int type_time
Gives the type of time periodicity.
virtual double val_boundary(int, const Val_domain &, const Index &) const
Computes the value of a field at a boundary.
virtual const Point absol_to_num(const Point &) const
Computes the numerical coordinates from the physical ones.
virtual Val_domain div_r(const Val_domain &) const
Division by .
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_r(const Val_domain &) const
Multiplication by .
virtual void set_cheb_base_odd(Base_spectral &) const
Gives the base using odd Chebyshev polynomials$.
virtual ostream & print(ostream &o) const
Delegate function to virtualize the << operator.
Domain_spheric_periodic_shell(int num, int ttype, int typet, double r_int, double r_ext, double ome, const Dim_array &nbr)
Standard constructor :
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.
int nbr_conditions_val_domain(const Val_domain &so, int order) const
Computes number of discretized equations associated with a given tensorial equation in the bulk.
virtual Val_domain der_normal(const Val_domain &, int) const
Normal derivative with respect to a given surface.
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
The class Point is used to store the coordinates of a point.
Definition: point.hpp:30
The Space_oned class fills the space with 2-dimensioanl domains : spherically symmetric and time peri...
void add_eq_ori(System_of_eqs &syst, const char *eq)
Adds an equation being the value of some field at the origin.
double omega
The pulsation of the space.
void add_matching(System_of_eqs &syst, const char *eq, int nused=-1, Array< int > **pused=0x0)
Adds a matching condition, at all the interface present in a given system.
void add_eq_inverted(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 (from outer domains to inner ones)
void add_eq_full(System_of_eqs &syst, const char *eq, int nused=-1, Array< int > **pused=0x0)
Adds a bulk equation in all the domains of a given system (equation is assumed to be 0th order)
void add_eq(System_of_eqs &syst, const char *eq, int nused=-1, Array< int > **pused=0x0)
Adds a bulk equation in all the domains of a given system (equation is assumed to be second order)
virtual void save(FILE *) const
Saving function.
void add_outer_bc(System_of_eqs &syst, const char *eq, int nused=-1, Array< int > **pused=0x0)
Sets a boundary condition at the outer radius of the compactified domain.
Space_spheric_periodic(int ttype, int typet, const Dim_array &nbr, const Array< double > &bounds, double ome)
Standard constructor.
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