28char map_af_integ_C[] =
"$Header: /cvsroot/Lorene/C++/Source/Map/map_af_integ.C,v 1.9 2014/10/13 08:53:02 j_novak Exp $" ;
83 static double* cx_tcp = 0 ;
84 static double* cx_rrp = 0 ;
85 static double* cx_rf_x2 = 0 ;
86 static double* cx_rf_x = 0 ;
87 static double* cx_rf = 0 ;
89 static int nt_cp_pre = 0 ;
90 static int nr_p_pre = 0 ;
91 static int nr_f_pre = 0 ;
93 assert(ci.
get_etat() != ETATNONDEF) ;
104 assert( ci.
get_etat() == ETATQCQ ) ;
112 assert(p_mti->
get_etat() == ETATQCQ) ;
118 for (
int l=0 ; l<nz ; l++) {
120 const Tbl* p_tbi = p_mti->
t[l] ;
122 if ( p_tbi->
get_etat() == ETATZERO ) {
127 assert( p_tbi->
get_etat() == ETATQCQ ) ;
132 int base = (p_mti->
base).get_b(l) ;
133 int base_r = base &
MSQ_R ;
134 int base_t = base &
MSQ_T ;
135 int base_p = base &
MSQ_P ;
141 double* s_tr =
new double[nr] ;
142 double* x_spec = p_tbi->
t ;
147 if (nt > nt_cp_pre) {
150 =
static_cast<double*
>(realloc(cx_tcp, nt*
sizeof(
double))) ;
151 for (
int j=0 ; j<nt ; j++) {
152 cx_tcp[j] = 2./(1. - 4.*j*j) ;
157 for (
int i=0 ; i<nr ; i++) s_tr[i] = 0 ;
158 for (
int j=0 ; j<nt ; j++) {
159 for (
int i=0 ; i<nr ; i++) {
160 s_tr[i] += cx_tcp[j] * x_spec[i] ;
168 for (
int i=0 ; i<nr ; i++) s_tr[i] = 0 ;
169 for (
int j=0 ; j<nt ; j++) {
171 for (
int i=0 ; i<nr ; i++) {
172 s_tr[i] += (2. / (1.-j*j)) * x_spec[i] ;
179 cout <<
"Map_af::integrale: unknown theta basis ! " << endl ;
197 assert(
beta[l] == 0) ;
200 cx_rrp =
static_cast<double*
>(realloc(cx_rrp, nr*
sizeof(
double))) ;
201 for (
int i=0 ; i<nr ; i++) {
202 cx_rrp[i] = (3. - 4.*i*i) /
203 (9. - 40. * i*i + 16. * i*i*i*i) ;
207 for (
int i=0 ; i<nr ; i++) {
208 som += cx_rrp[i] * s_tr[i] ;
210 double rmax =
alpha[l] ;
211 som *= rmax*rmax*rmax ;
218 cx_rf_x2 =
static_cast<double*
>(realloc(cx_rf_x2, nr*
sizeof(
double))) ;
219 cx_rf_x =
static_cast<double*
>(realloc(cx_rf_x, nr*
sizeof(
double))) ;
220 cx_rf =
static_cast<double*
>(realloc(cx_rf, nr*
sizeof(
double))) ;
221 for (
int i=0 ; i<nr ; i +=2 ) {
222 cx_rf_x2[i] = 2.*(3. - i*i)/(9. - 10. * i*i + i*i*i*i) ;
224 cx_rf[i] = 2./(1. - i*i) ;
226 for (
int i=1 ; i<nr ; i +=2 ) {
228 cx_rf_x[i] = 2./(4. - i*i) ;
233 for (
int i=0 ; i<nr ; i +=2 ) {
234 som_x2 += cx_rf_x2[i] * s_tr[i] ;
236 for (
int i=1 ; i<nr ; i +=2 ) {
237 som_x += cx_rf_x[i] * s_tr[i] ;
239 for (
int i=0 ; i<nr ; i +=2 ) {
240 som_c += cx_rf[i] * s_tr[i] ;
242 double a =
alpha[l] ;
244 som = a*a*a * som_x2 + 2.*a*a*b * som_x + a*b*b * som_c ;
251 cx_rf_x2 =
static_cast<double*
>(realloc(cx_rf_x2, nr*
sizeof(
double))) ;
252 cx_rf_x =
static_cast<double*
>(realloc(cx_rf_x, nr*
sizeof(
double))) ;
253 cx_rf =
static_cast<double*
>(realloc(cx_rf, nr*
sizeof(
double))) ;
255 for (
int i=0 ; i<nr ; i +=1 ) {
257 cx_rf_x[i] = 2*signe/double(i+1)/double(i+2);
261 cx_rf_x2[0] = double(8)/double(3) ;
264 for (
int i=0 ; i<nr ; i +=1 ) {
265 som_x2 += cx_rf_x2[i] * s_tr[i] ;
267 for (
int i=1 ; i<nr ; i +=1 ) {
268 som_x += cx_rf_x[i] * s_tr[i] ;
270 for (
int i=0 ; i<nr ; i +=1 ) {
271 som_c += cx_rf[i] * s_tr[i] ;
273 double a =
alpha[l] ;
276 som = a*a*a * som_x2 + 2.*a*a*(b-a) * som_x + a*(b-a)*(b-a) * som_c ;
284 cx_rf =
static_cast<double*
>(realloc(cx_rf, nr*
sizeof(
double))) ;
285 for (
int i=0 ; i<nr ; i +=2 ) {
286 cx_rf[i] = 2./(1. - i*i) ;
288 for (
int i=1 ; i<nr ; i +=2 ) {
293 for (
int i=0 ; i<nr ; i +=2 ) {
294 som_c += cx_rf[i] * s_tr[i] ;
296 som = -
alpha[l] * som_c ;
302 cout <<
"Map_af::integrale: unknown r basis ! " << endl ;
326 cout <<
"Map_af::integrale: unknown phi basis ! " << endl ;
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
int get_etat() const
Returns the logical state.
Valeur va
The numerical value of the Cmp
bool check_dzpuis(int dzi) const
Returns false if the last domain is compactified and *this is not zero in this domain and dzpuis is n...
double * beta
Array (size: mg->nzone ) of the values of in each domain.
double * alpha
Array (size: mg->nzone ) of the values of in each domain.
virtual Tbl * integrale(const Cmp &) const
Computes the integral over all space of a Cmp.
const Mg3d * mg
Pointer on the multi-grid Mgd3 on which this is defined
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.
Coefficients storage for the multi-domain spectral method.
Base_val base
Bases of the spectral expansions.
int get_etat() const
Returns the logical state.
Tbl ** t
Array (size nzone ) of pointers on the Tbl 's which contain the spectral coefficients in each domain.
int get_etat() const
Gives the logical state.
void annule_hard()
Sets the Tbl to zero in a hard way.
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
double * t
The array of double.
#define P_COSSIN_P
dev. sur Phi = 2*phi, freq. paires
#define R_CHEBU
base de Chebychev ordinaire (fin), dev. en 1/r
#define R_JACO02
base de Jacobi(0,2) ordinaire (finjac)
#define MSQ_R
Extraction de l'info sur R.
#define T_COS_P
dev. cos seulement, harmoniques paires
#define P_COSSIN
dev. standart
#define R_CHEBPIM_P
Cheb. pair-impair suivant m, pair pour m=0.
#define MSQ_T
Extraction de l'info sur Theta.
#define T_COSSIN_CP
cos pair-sin impair alternes, cos pour m=0
#define R_CHEB
base de Chebychev ordinaire (fin)
#define T_COS
dev. cos seulement
#define R_CHEBP
base de Cheb. paire (rare) seulement
#define MSQ_P
Extraction de l'info sur Phi.
#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.