20 #include "headcpp.hpp"
21 #include "utilities.hpp"
22 #include "spheric_symphi.hpp"
24 #include "array_math.hpp"
25 #include "val_domain.hpp"
28 void coef_1d (
int, Array<double>&) ;
29 void coef_i_1d (
int, Array<double>&) ;
30 int der_1d (
int, Array<double>&) ;
34 Domain(num, ttype, nbr), alpha(r),center(cr) {
46 fread_be (&
alpha,
sizeof(
double), 1, fd) ;
51 Domain_nucleus_symphi::~Domain_nucleus_symphi() {}
56 fwrite_be (&
ndim,
sizeof(
int), 1, fd) ;
57 fwrite_be (&
type_base,
sizeof(
int), 1, fd) ;
59 fwrite_be (&
alpha,
sizeof(
double), 1, fd) ;
63 o <<
"Nucleus Sym Phi" << endl ;
64 o <<
"Rmax = " <<
alpha << endl ;
65 o <<
"Center = " <<
center << endl ;
80 cerr <<
"Unknown boundary case in Domain_nucleus_symphi::der_normal" << endl ;
88 for (
int i=0 ; i<3 ; i++)
89 assert (
coloc[i] != 0x0) ;
90 for (
int i=0 ; i<3 ; i++)
91 assert (
absol[i] == 0x0) ;
92 for (
int i=0 ; i<3 ; i++) {
94 absol[i]->allocate_conf() ;
104 while (index.
inc()) ;
111 for (
int i=0 ; i<3 ; i++)
112 assert (
coloc[i] != 0x0) ;
119 while (index.
inc()) ;
129 cerr <<
"Unknown type of basis in Domain_nucleus_symphi::do_radius" << endl ;
137 for (
int i=0 ; i<3 ; i++)
138 assert (
coloc[i] != 0x0) ;
139 for (
int i=0 ; i<3 ; i++)
140 assert (
cart[i] == 0x0) ;
141 for (
int i=0 ; i<3 ; i++) {
143 cart[i]->allocate_conf() ;
153 while (index.
inc()) ;
168 cerr <<
"Unknown type of basis in Domain_nucleus_symphi::do_cart" << endl ;
175 for (
int i=0 ; i<3 ; i++)
176 assert (
coloc[i] != 0x0) ;
177 for (
int i=0 ; i<3 ; i++)
179 for (
int i=0 ; i<3 ; i++) {
189 while (index.
inc()) ;
197 double x_loc = xx(1) -
center(1) ;
198 double y_loc = xx(2) -
center(2) ;
199 double z_loc = xx(3) -
center(3) ;
200 double air_loc = sqrt (x_loc*x_loc + y_loc*y_loc + z_loc*z_loc) ;
202 bool res = (air_loc <=
alpha+prec) ?
true :
false ;
210 assert (
is_in(abs)) ;
213 double x_loc = abs(1) -
center(1) ;
214 double y_loc = abs(2) -
center(2) ;
215 double z_loc = abs(3) -
center(3) ;
216 double air = sqrt(x_loc*x_loc+y_loc*y_loc+z_loc*z_loc) ;
218 double rho = sqrt(x_loc*x_loc+y_loc*y_loc) ;
222 num.
set(2) = (z_loc>=0) ? 0 : M_PI ;
226 num.
set(2) = atan(rho/z_loc) ;
227 num.
set(3) = atan2 (y_loc, x_loc) ;
231 num.
set(2) = M_PI + num(2) ;
240 assert (bound==OUTER_BC) ;
241 assert (
is_in(abs, 1e-3)) ;
244 double x_loc = abs(1) -
center(1) ;
245 double y_loc = abs(2) -
center(2) ;
246 double z_loc = abs(3) -
center(3) ;
248 double rho = sqrt(x_loc*x_loc+y_loc*y_loc) ;
252 num.
set(2) = (z_loc>=0) ? 0 : M_PI ;
256 num.
set(2) = atan(rho/z_loc) ;
257 num.
set(3) = atan2 (y_loc, x_loc) ;
261 num.
set(2) = M_PI + num(2) ;
266 double coloc_leg_parity(
int,
int) ;
273 for (
int i=0 ; i<
ndim ; i++)
285 for (
int i=0 ; i<
ndim ; i++)
295 cerr <<
"Unknown type of basis in Domain_nucleus_symphi::do_coloc" << endl ;
306 base.
bases_1d[2]->set(0) = COS_EVEN ;
308 base.
bases_1d[1]->set(k) = COS_EVEN ;
310 index.
set(0) = j ; index.
set(1) = k ;
311 base.
bases_1d[0]->set(index) = CHEB_EVEN ;
322 base.
bases_1d[2]->set(0) = SIN_EVEN ;
324 base.
bases_1d[1]->set(k) = COS_EVEN ;
326 index.
set(0) = j ; index.
set(1) = k ;
327 base.
bases_1d[0]->set(index) = CHEB_ODD ;
338 base.
bases_1d[2]->set(0) = SIN_EVEN ;
340 base.
bases_1d[1]->set(k) = SIN_EVEN ;
342 index.
set(0) = j ; index.
set(1) = k ;
343 base.
bases_1d[0]->set(index) = CHEB_ODD ;
354 base.
bases_1d[2]->set(0) = COS_EVEN ;
356 base.
bases_1d[1]->set(k) = SIN_ODD ;
358 index.
set(0) = j ; index.
set(1) = k ;
359 base.
bases_1d[0]->set(index) = CHEB_ODD ;
370 base.
bases_1d[2]->set(0) = SIN_ODD ;
372 base.
bases_1d[1]->set(k) = SIN_ODD ;
374 index.
set(0) = j ; index.
set(1) = k ;
375 base.
bases_1d[0]->set(index) = CHEB_ODD ;
386 base.
bases_1d[2]->set(0) = COS_ODD ;
388 base.
bases_1d[1]->set(k) = SIN_ODD ;
390 index.
set(0) = j ; index.
set(1) = k ;
391 base.
bases_1d[0]->set(index) = CHEB_ODD ;
402 base.
bases_1d[2]->set(0) = SIN_EVEN ;
404 base.
bases_1d[1]->set(k) = COS_ODD ;
406 index.
set(0) = j ; index.
set(1) = k ;
407 base.
bases_1d[0]->set(index) = CHEB_ODD ;
418 base.
bases_1d[2]->set(0) = COS_EVEN ;
420 base.
bases_1d[1]->set(k) = SIN_EVEN ;
422 index.
set(0) = j ; index.
set(1) = k ;
423 base.
bases_1d[0]->set(index) = CHEB_EVEN ;
434 base.
bases_1d[2]->set(0) = SIN_EVEN ;
436 base.
bases_1d[1]->set(k) = SIN_ODD ;
438 index.
set(0) = j ; index.
set(1) = k ;
439 base.
bases_1d[0]->set(index) = CHEB_EVEN ;
450 base.
bases_1d[2]->set(0) = SIN_EVEN ;
452 base.
bases_1d[1]->set(k) = COS_ODD ;
454 index.
set(0) = j ; index.
set(1) = k ;
455 base.
bases_1d[0]->set(index) = CHEB_EVEN ;
466 base.
bases_1d[2]->set(0) = SIN_EVEN ;
468 base.
bases_1d[1]->set(k) = COS_EVEN ;
470 index.
set(0) = j ; index.
set(1) = k ;
471 base.
bases_1d[0]->set(index) = CHEB_EVEN ;
482 base.
bases_1d[2]->set(0) = COS_ODD ;
484 base.
bases_1d[1]->set(k) = SIN_EVEN ;
486 index.
set(0) = j ; index.
set(1) = k ;
487 base.
bases_1d[0]->set(index) = CHEB_EVEN ;
498 base.
bases_1d[2]->set(0) = SIN_ODD ;
500 base.
bases_1d[1]->set(k) = SIN_EVEN ;
502 index.
set(0) = j ; index.
set(1) = k ;
503 base.
bases_1d[0]->set(index) = CHEB_EVEN ;
514 base.
bases_1d[2]->set(0) = COS_EVEN ;
516 base.
bases_1d[1]->set(k) = COS_EVEN ;
518 index.
set(0) = j ; index.
set(1) = k ;
519 base.
bases_1d[0]->set(index) = CHEB_ODD ;
532 base.
bases_1d[2]->set(0) = COS_ODD ;
534 base.
bases_1d[1]->set(k) = SIN_ODD ;
536 index.
set(0) = j ; index.
set(1) = k ;
537 base.
bases_1d[0]->set(index) = CHEB_ODD ;
550 base.
bases_1d[2]->set(0) = SIN_ODD ;
552 base.
bases_1d[1]->set(k) = SIN_ODD ;
554 index.
set(0) = j ; index.
set(1) = k ;
555 base.
bases_1d[0]->set(index) = CHEB_ODD ;
569 base.
bases_1d[2]->set(0) = COS_EVEN ;
571 base.
bases_1d[1]->set(k) = COS_ODD ;
573 index.
set(0) = j ; index.
set(1) = k ;
574 base.
bases_1d[0]->set(index) = CHEB_ODD ;
588 base.
bases_1d[2]->set(0) = COS_EVEN ;
590 base.
bases_1d[1]->set(k) = COS_EVEN ;
592 index.
set(0) = j ; index.
set(1) = k ;
593 base.
bases_1d[0]->set(index) = LEG_EVEN ;
607 base.
bases_1d[2]->set(0) = SIN_EVEN ;
609 base.
bases_1d[1]->set(k) = COS_EVEN ;
611 index.
set(0) = j ; index.
set(1) = k ;
612 base.
bases_1d[0]->set(index) = LEG_ODD ;
625 base.
bases_1d[2]->set(0) = SIN_EVEN ;
627 base.
bases_1d[1]->set(k) = SIN_EVEN ;
629 index.
set(0) = j ; index.
set(1) = k ;
630 base.
bases_1d[0]->set(index) = LEG_ODD ;
643 base.
bases_1d[2]->set(0) = COS_EVEN ;
645 base.
bases_1d[1]->set(k) = SIN_ODD ;
647 index.
set(0) = j ; index.
set(1) = k ;
648 base.
bases_1d[0]->set(index) = LEG_ODD ;
661 base.
bases_1d[2]->set(0) = SIN_ODD ;
663 base.
bases_1d[1]->set(k) = SIN_ODD ;
665 index.
set(0) = j ; index.
set(1) = k ;
666 base.
bases_1d[0]->set(index) = LEG_ODD ;
679 base.
bases_1d[2]->set(0) = COS_ODD ;
681 base.
bases_1d[1]->set(k) = SIN_ODD ;
683 index.
set(0) = j ; index.
set(1) = k ;
684 base.
bases_1d[0]->set(index) = LEG_ODD ;
698 base.
bases_1d[2]->set(0) = SIN_EVEN ;
700 base.
bases_1d[1]->set(k) = COS_ODD ;
702 index.
set(0) = j ; index.
set(1) = k ;
703 base.
bases_1d[0]->set(index) = LEG_ODD ;
717 base.
bases_1d[2]->set(0) = COS_EVEN ;
719 base.
bases_1d[1]->set(k) = COS_EVEN ;
721 index.
set(0) = j ; index.
set(1) = k ;
722 base.
bases_1d[0]->set(index) = LEG_ODD ;
735 base.
bases_1d[2]->set(0) = COS_ODD ;
737 base.
bases_1d[1]->set(k) = SIN_ODD ;
739 index.
set(0) = j ; index.
set(1) = k ;
740 base.
bases_1d[0]->set(index) = LEG_ODD ;
753 base.
bases_1d[2]->set(0) = SIN_ODD ;
755 base.
bases_1d[1]->set(k) = SIN_ODD ;
757 index.
set(0) = j ; index.
set(1) = k ;
758 base.
bases_1d[0]->set(index) = LEG_ODD ;
772 base.
bases_1d[2]->set(0) = COS_EVEN ;
774 base.
bases_1d[1]->set(k) = COS_ODD ;
776 index.
set(0) = j ; index.
set(1) = k ;
777 base.
bases_1d[0]->set(index) = LEG_ODD ;
804 if (strcmp(p,
"R ")==0)
806 if (strcmp(p,
"T ")==0)
808 if (strcmp(p,
"P ")==0)
815 if (bound != OUTER_BC)
817 cerr <<
"Domain_nucleus_symphi::integ only defined for r=rmax yet..." << endl ;
821 int baset((*so.base.bases_1d[1])(0));
822 if (baset != COS_EVEN)
828 Index pos(get_nbr_coefs());
830 for (
int j(0) ; j < nbr_coefs(1) ; ++j)
833 double fact_tet(2.0/(1.0 - 4.0*j*j));
835 for (
int i(0) ; i < nbr_coefs(0) ; ++i)
838 res += fact_tet*(*so.cf)(pos);
851 if (!a.
def) res_def =
false;
852 if (!b.
def) res_def =
false;
861 int basea((*a.
bases_1d[2])(index_2));
862 int baseb((*b.
bases_1d[2])(index_2));
864 }
while(index_2.
inc());
872 int basea((*a.
bases_1d[1])(index_1));
873 int baseb((*b.
bases_1d[1])(index_1));
875 }
while(index_1.
inc());
882 int basea((*a.
bases_1d[0])(index_0));
883 int baseb((*b.
bases_1d[0])(index_0));
885 }
while(index_0.
inc());
888 for (
int dim(0) ; dim < a.
ndim ; ++dim)
901 if ( basea == COS_EVEN and baseb == COS_EVEN) res = COS_EVEN;
902 if ((basea == COS_EVEN and baseb == COS_ODD) or (basea == COS_ODD and baseb == COS_EVEN)) res = COS_ODD;
903 if ((basea == COS_EVEN and baseb == SIN_EVEN) or (basea == SIN_EVEN and baseb == COS_EVEN)) res = SIN_EVEN;
904 if ((basea == COS_EVEN and baseb == SIN_ODD) or (basea == SIN_ODD and baseb == COS_EVEN)) res = SIN_ODD;
905 if ( basea == COS_ODD and baseb == COS_ODD) res = COS_EVEN;
906 if ((basea == COS_ODD and baseb == SIN_EVEN) or (basea == SIN_EVEN and baseb == COS_ODD)) res = SIN_ODD;
907 if ((basea == COS_ODD and baseb == SIN_ODD) or (basea == SIN_ODD and baseb == COS_ODD)) res = SIN_EVEN;
908 if ( basea == SIN_EVEN and baseb == SIN_EVEN) res = COS_EVEN;
909 if ((basea == SIN_EVEN and baseb == SIN_ODD) or (basea == SIN_ODD and baseb == SIN_EVEN)) res = COS_ODD;
910 if ( basea == SIN_ODD and baseb == SIN_ODD) res = COS_EVEN;
917 if ( basea == CHEB_EVEN and baseb == CHEB_EVEN) res = CHEB_EVEN;
918 if ((basea == CHEB_EVEN and baseb == CHEB_ODD) or (basea == CHEB_ODD and baseb == CHEB_EVEN)) res = CHEB_ODD;
919 if ( basea == CHEB_ODD and baseb == CHEB_ODD) res = CHEB_EVEN;
Class for storing the basis of decompositions of a field.
Bases_container bases_1d
Arrays containing the various basis of decomposition.
void allocate(const Dim_array &nbr_coefs)
Allocates the various arrays, for a given number of coefficients.
bool def
true if the Base_spectral is defined and false otherwise.
int ndim
Number of dimensions.
Class for storing the dimensions of an array.
int get_ndim() const
Returns the number of dimensions.
void save(FILE *) const
Save function.
Class for a spherical domain containing the origin a symmetry with respect to the plane and an quadr...
void set_cheb_base_fory_cart(Base_spectral &so) const
Sets the base to the standard one for Chebyshev polynomials for a field like the component of a vect...
virtual void set_cheb_base_x_cart(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
virtual void do_coloc()
Computes the colocation points.
virtual void do_cart() const
Computes the Cartesian coordinates.
virtual void set_cheb_base_yz_cart(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
virtual void do_cart_surr() const
Computes the Cartesian coordinates over the radius.
virtual void set_legendre_base_y_cart(Base_spectral &) const
Gives the base using Legendre polynomials, for the component of a vector.
virtual Val_domain mult_sin_theta(const Val_domain &) const
Multiplication by .
virtual void do_der_abs_from_der_var(const Val_domain *const *const der_var, Val_domain **const der_abs) const
Computes the derivative with respect to the absolute Cartesian coordinates from the derivative with r...
virtual double integ(const Val_domain &so, int bound) const
Surface integral on a given boundary.
virtual int give_place_var(char *) const
Translates a name of a coordinate into its corresponding numerical name.
virtual Val_domain div_x(const Val_domain &) const
Division by .
virtual void set_cheb_base_xz_cart(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
int mult_base_r_int(int basea, int baseb) const
Multiply two radial basis.
virtual Val_domain mult_cos_theta(const Val_domain &) const
Multiplication by .
virtual Val_domain mult_sin_phi(const Val_domain &) const
Multiplication by .
virtual ostream & print(ostream &o) const
Delegate function to virtualize the << operator.
virtual void save(FILE *) const
Saving function.
virtual const Point absol_to_num_bound(const Point &, int) const
Computes the numerical coordinates from the physical ones for a point lying on a boundary.
virtual void set_legendre_base_p_spher(Base_spectral &) const
Gives the base using Legendre polynomials, for the component of a vector.
virtual const Point absol_to_num(const Point &xxx) const
Computes the numerical coordinates from the physical ones.
void set_cheb_base_forr(Base_spectral &so) const
Sets the base to the standard one for Chebyshev polynomials for a field like the radius .
void set_cheb_base_forx_cart(Base_spectral &so) const
Sets the base to the standard one for Chebyshev polynomials for a field like the component of a vect...
int mult_base_angle_int(int basea, int baseb) const
Multiply two angular basis.
virtual void set_cheb_base_y_cart(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
virtual Val_domain mult_cos_phi(const Val_domain &) const
Multiplication by .
double alpha
Relates the numerical to the physical radii.
virtual void set_legendre_base_x_cart(Base_spectral &) const
Gives the base using Legendre polynomials, for the component of a vector.
virtual void do_absol() const
Computes the absolute coordinates.
virtual bool is_in(const Point &xx, double prec=1e-13) const
Check whether a point lies inside Domain.
virtual void set_legendre_base_t_spher(Base_spectral &) const
Gives the base using Legendre polynomials, for the component of a vector.
virtual void set_cheb_base_p_spher(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
virtual void set_legendre_base(Base_spectral &) const
Gives the standard base for Legendre polynomials.
virtual void set_cheb_base_r_spher(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the radial component of a vector.
void set_legendre_base_forr(Base_spectral &so) const
Sets the base to the standard one for Legendre polynomials for a field like the radius .
Domain_nucleus_symphi(int num, int ttype, double radius, const Point &cr, const Dim_array &nbr)
Standard constructor :
void set_legendre_base_forx_cart(Base_spectral &so) const
Sets the base to the standard one for Legendre polynomials for a field like the component of a vecto...
virtual void set_cheb_base_rp_spher(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a 2-tensor.
void set_legendre_base_fory_cart(Base_spectral &so) const
Sets the base to the standard one for Legendre polynomials for a field like the component of a vecto...
virtual void set_cheb_base_xy_cart(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
virtual void set_cheb_base_z_cart(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
virtual void set_legendre_base_z_cart(Base_spectral &) const
Gives the base using Legendre polynomials, for the component of a vector.
virtual void set_cheb_base(Base_spectral &) const
Gives the standard base for Chebyshev polynomials.
virtual void set_cheb_base_t_spher(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a vector.
Point center
Absolute coordinates of the center.
virtual Val_domain der_normal(const Val_domain &, int) const
Normal derivative with respect to a given surface.
virtual void set_cheb_base_tp_spher(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a 2-tensor.
virtual Base_spectral mult(const Base_spectral &, const Base_spectral &) const
Method for the multiplication of two Base_spectral.
virtual void set_legendre_base_r_spher(Base_spectral &) const
Gives the base using Legendre polynomials, for the radial component of a vector.
virtual void set_cheb_base_rt_spher(Base_spectral &) const
Gives the base using Chebyshev polynomials, for the component of a 2-tensor.
virtual void do_radius() const
Computes the generalized radius.
void set_cheb_base_forz_cart(Base_spectral &so) const
Sets the base to the standard one for Chebyshev polynomials for a field like the component of a vect...
void set_legendre_base_forz_cart(Base_spectral &so) const
Sets the base to the standard one for Legendre polynomials for a field like the component of a vecto...
Abstract class that implements the fonctionnalities common to all the type of domains.
virtual void del_deriv()
Destroys the derivated members (like coloc, cart and radius), when changing the type of colocation po...
Val_domain * radius
The generalized radius.
Memory_mapped_array< Val_domain * > cart
Cartesian coordinates.
Memory_mapped_array< Val_domain * > absol
Asbolute coordinates (if defined ; usually Cartesian-like)
int ndim
Number of dimensions.
Memory_mapped_array< Val_domain * > cart_surr
Cartesian coordinates divided by the radius.
Dim_array nbr_coefs
Number of coefficients.
Dim_array nbr_points
Number of colocation points.
int type_base
Type of colocation point :
Memory_mapped_array< Array< double > * > coloc
Colocation points in each dimension (stored in ndim 1d- arrays)
Class that gives the position inside a multi-dimensional Array.
int & set(int i)
Read/write of the position in a given dimension.
bool inc(int increm, int var=0)
Increments the position of the Index.
The class Point is used to store the coordinates of a point.
void save(FILE *) const
Saving function.
const int & get_ndim() const
Returns the number of dimensions.
double & set(int i)
Read/write of a coordinate.
Class for storing the basis of decompositions of a field and its values on both the configuration and...
Val_domain mult_sin_phi() const
Multiplication by .
Val_domain mult_sin_theta() const
Multiplication by .
Val_domain mult_cos_phi() const
Multiplication by .
double & set(const Index &pos)
Read/write the value of the field in the configuration space.
Val_domain div_sin_theta() const
Division by .
Val_domain mult_cos_theta() const
Multiplication by .
Val_domain der_var(int i) const
Computes the derivative with respect to a numerical coordinate.
Base_spectral & set_base()
Sets the basis of decomposition.
void allocate_conf()
Allocates the values in the configuration space and destroys the values in the coefficients space.