26char gval_from_spectral_C[] =
"$Header: /cvsroot/Lorene/C++/Source/Valencia/gval_from_spectral.C,v 1.14 2014/10/13 08:53:48 j_novak Exp $" ;
81#include "grille_val.h"
94 cout <<
"Gval_spher::somme_spectral2():\n" ;
95 cout <<
"grid size incompatible with array size... exiting!" <<
endl ;
105 resu[
i] =
meudon.get_spectral_va().val_point_jk(
l, xi, 0, 0) ;
117 cout <<
"Gval_spher::somme_spectral2():\n" ;
118 cout <<
"grid size incompatible with array size... exiting!" <<
endl ;
167 cout <<
"Gval_spher::somme_spectral2():\n" ;
168 cout <<
"grid size incompatible with array size... exiting!" <<
endl ;
238 cout <<
"Gval_spher::somme_spectral2():\n" ;
239 cout <<
"grid size incompatible with array size... exiting!" <<
endl ;
285 double*
resu =
new double[taille] ;
329 double*
resu =
new double[taille] ;
369 meudon.get_spectral_va().coef() ;
383 cout <<
"Gval_spher::somme_spectral3():\n" ;
384 cout <<
"grid size incompatible with array size... exiting!" <<
endl ;
397 for (
int lz=1;
lz<nz;
lz++) {
405 double* alpha =
new double[
nrv0*(
npm+2)*
ntm] ;
412 double**
coefm =
new double*[nz] ;
413 for (
int lz=0;
lz<nz;
lz++) {
435 *
p_coef += (*tbcf)*(*p_func) ;
445 for (
int lz=0;
lz<nz;
lz++) {
454 double*
tetlj = 0x0 ;
455 initialize_spectral_theta(mp,
meudon.get_spectral_va().get_base(),
tetlj) ;
466 *
p_coef += (*p_interm) * (*p_func) ;
486 double*
expmk = 0x0 ;
487 initialize_spectral_phi(mp,
meudon.get_spectral_va().get_base(),
expmk) ;
513 *
p_coef += (*p_interm) * (*p_func) ;
547void Gval_spher::initialize_spectral_r(
const Map& mp,
const Base_val& base,
557 double* xi =
new double[
nrv0] ;
565 assert (chebnri == 0x0) ;
566 chebnri =
new double[(npm+2)*ntm*nrmax] ;
567 double* p_out = chebnri ;
568 for (
int irv=0; irv<nrv0; irv++) {
569 bool nucleus = (mg->
get_type_r(idom[irv]) == RARE) ;
570 int nmax = (nucleus ? 2*mg->
get_nr(idom[irv]) + 1
571 : mg->
get_nr(idom[irv])) ;
572 double* cheb =
new double[nmax] ;
575 for (
int ir=2; ir<nmax; ir++) {
576 cheb[ir] = 2*xi[irv]*cheb[ir-1] - cheb[ir-2] ;
581 for (
int ip=0; ip<npm+2; ip++) {
582 for (
int it=0; it<ntm; it++) {
613 par = 1 - ((ip/2) % 2) ;
618 cout <<
"Gval_spher::initialize_spectral_r : " <<
'\n'
619 <<
"Unexpected radial base !" <<
'\n'
620 <<
"Base : " << base_r << endl ;
626 for (
int ir=0; ir<mg->
get_nr(idom[irv]); ir++) {
627 *p_out = cheb[fact*ir+par] ;
641void Gval_spher::initialize_spectral_theta(
const Map& mp,
const Base_val& base,
642 double*& tetlj)
const {
645 const Mg3d* mg = mp.get_mg() ;
647 int ntm = mg->get_nt(0) ;
648 int base_t = base.get_base_t(0) ;
650 assert (tetlj == 0x0) ;
651 tetlj =
new double[(npm+2)*ntv0*ntm] ;
652 double* p_out = tetlj ;
654 for (
int jtv=0; jtv<ntv0; jtv++) {
656 for (
int mpm=0; mpm < npm+2; mpm++) {
657 for (
int ltm=0; ltm<ntm; ltm++) {
660 *p_out =
cos(ltm*teta) ;
664 *p_out =
sin(ltm*teta) ;
668 *p_out =
cos(2*ltm*teta) ;
672 *p_out =
cos((2*ltm+1)*teta) ;
676 *p_out =
sin(2*ltm*teta) ;
680 *p_out =
sin((2*ltm+1)*teta) ;
684 *p_out = ( ((mpm/2) % 2 == 0) ?
cos(2*ltm*teta)
685 :
sin((2*ltm+1)*teta)) ;
689 *p_out = ( ((mpm/2) % 2 == 0) ?
cos((2*ltm+1)*teta)
694 *p_out = ( ((mpm/2) % 2 == 0) ?
sin(2*ltm*teta)
695 :
cos((2*ltm+1)*teta)) ;
699 *p_out = ( ((mpm/2) % 2 == 0) ?
sin((2*ltm+1)*teta)
704 *p_out = ( ((mpm/2) % 2 == 0) ?
cos(ltm*teta)
709 *p_out = ( ((mpm/2) % 2 == 0) ?
sin(ltm*teta)
714 cout <<
"Gval_spher::initialize_spectral_theta : " <<
'\n'
715 <<
"Unexpected theta base !" <<
'\n'
716 <<
"Base : " << base_t << endl ;
735void Gval_spher::initialize_spectral_phi(
const Map& mp,
const Base_val& base,
736 double*& expmk)
const {
739 const Mg3d* mg = mp.get_mg() ;
740 int npm = mg->get_np(0) ;
741 int base_p = base.get_base_p(0) ;
743 assert (expmk == 0x0) ;
744 expmk =
new double[(npm+2)*npv0] ;
745 double* p_out = expmk ;
747 for (
int kpv=0; kpv<npv0; kpv++) {
749 for (
int mpm=0; mpm < npm+2; mpm++) {
753 *p_out = ( (mpm%2 == 0) ?
cos(m*fi) :
sin(m*fi) ) ;
758 *p_out = ( (mpm%2 == 0) ?
cos(2*m*fi) :
sin(2*m*fi) ) ;
763 *p_out = ( (mpm%2 == 0) ?
cos((2*m+1)*fi) :
sin((2*m+1)*fi) ) ;
767 cout <<
"Gval_spher::initialize_spectral_phi : " <<
'\n'
768 <<
"Unexpected phi base !" <<
'\n'
769 <<
"Base : " << base_p << endl ;
Bases of the spectral expansions.
int get_base_r(int l) const
Returns the expansion basis for r ( ) functions in the domain of index l (e.g.
int * dim
Array of dimensions (size: ndim).
int ndim
Number of dimensions of the Tbl: can be 1, 2 or 3.
Time evolution with partial storage (*** under development ***).
int nfantome
The number of hidden cells (same on each side)
void somme_spectrale1(const Scalar &meudon, double *t, int taille) const
Makes the sommation of the spectral basis functions to know the values of the function described by t...
Dim_tbl dim
The dimensions of the grid.
Tbl * zr
Arrays containing the values of coordinate z (or r) on the nodes
Tbl * zri
Arrays containing the values of coordinate z (or r) on the interfaces.
virtual void somme_spectrale2(const Scalar &meudon, double *t, int taille) const
Makes the sommation of the spectral basis functions to know the values of the function described by t...
virtual void somme_spectrale3(const Scalar &meudon, double *t, int taille) const
Same as before but for the 3D case.
Tbl * x
Arrays containing the values of coordinate x on the nodes.
virtual void somme_spectrale3(const Scalar &meudon, double *t, int taille) const
Same as before but for the 3D case.
virtual void somme_spectrale2(const Scalar &meudon, double *t, int taille) const
Makes the sommation of the spectral basis functions to know the values of the function described by t...
Tbl * phi
Arrays containing the values of coordinate on the nodes.
Tbl * teti
Arrays containing the values of coordinate on the interfaces.
double * somme_spectrale2ri(const Scalar &meudon) const
Same as before but at radial grid interfaces.
double * somme_spectrale2ti(const Scalar &meudon) const
Same as before but at angular grid interfaces.
Tbl * tet
Arrays containing the values of coordinate on the nodes.
Base class for coordinate mappings.
virtual void val_lx(double rr, double theta, double pphi, int &l, double &xi) const =0
Computes the domain index l and the value of corresponding to a point given by its physical coordina...
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
int get_np(int l) const
Returns the number of points in the azimuthal direction ( ) in domain no. l.
int get_nt(int l) const
Returns the number of points in the co-latitude direction ( ) in domain no. l.
int get_nzone() const
Returns the number of domains.
int get_nr(int l) const
Returns the number of points in the radial direction ( ) in domain no. l.
int get_type_r(int l) const
Returns the type of sampling in the radial direction in domain no.
Coefficients storage for the multi-domain spectral method.
Tensor field of valence 0 (or component of a tensorial field).
double * t
The array of double.
Cmp sqrt(const Cmp &)
Square root.
Cmp sin(const Cmp &)
Sine.
Cmp acos(const Cmp &)
Arccosine.
Cmp cos(const Cmp &)
Cosine.
#define P_COSSIN_P
dev. sur Phi = 2*phi, freq. paires
#define T_COSSIN_SP
sin pair-cos impair alternes, sin pour m=0
#define R_CHEBI
base de Cheb. impaire (rare) seulement
#define T_SIN_P
dev. sin seulement, harmoniques paires
#define T_COSSIN_S
dev. cos-sin alternes, sin pour m=0
#define R_CHEBPIM_I
Cheb. pair-impair suivant m, impair pour m=0.
#define T_COSSIN_SI
sin impair-cos pair alternes, sin pour m=0
#define R_CHEBPI_I
Cheb. pair-impair suivant l impair pour l=0.
#define T_COS_P
dev. cos seulement, harmoniques paires
#define T_COSSIN_CI
cos impair-sin pair alternes, cos pour m=0
#define P_COSSIN
dev. standart
#define P_COSSIN_I
dev. sur Phi = 2*phi, freq. impaires
#define R_CHEBPIM_P
Cheb. pair-impair suivant m, pair pour m=0.
#define T_COSSIN_CP
cos pair-sin impair alternes, cos pour m=0
#define T_SIN_I
dev. sin seulement, harmoniques impaires
#define T_COS
dev. cos seulement
#define R_CHEBP
base de Cheb. paire (rare) seulement
#define T_SIN
dev. sin seulement
#define T_COS_I
dev. cos seulement, harmoniques impaires
#define T_COSSIN_C
dev. cos-sin alternes, cos pour m=0
#define R_CHEBPI_P
Cheb. pair-impair suivant l pair pour l=0.