23char map_af_dalembert_C[] =
"$Header: /cvsroot/Lorene/C++/Source/Map/map_af_dalembert.C,v 1.18 2014/10/13 08:53:02 j_novak Exp $" ;
123 const Scalar& source)
const {
126 assert(source.
get_etat() != ETATNONDEF) ;
128 assert(fj.
get_etat() != ETATNONDEF) ;
130 assert(fjm1.
get_etat() != ETATNONDEF) ;
142 int nz0 = (ced ? nz - 1 : nz) ;
146 Scalar fjm1_local = fjm1 ;
151 Scalar sigma = 2*fj_local - fjm1_local ;
158 coeff =
new Tbl(12,nz);
164 Tbl a1(nz) ; a1 = 1 ;
165 Tbl a2(nz) ; a2 = 0 ;
166 Tbl a3(nz) ; a3 = 0 ;
171 assert(metri != 0x0) ;
172 assert (metri->
get_etat() == ETATQCQ) ;
176 double* bornes =
new double[nz+1] ;
177 bornes[0] =
beta[0] ;
178 for (
int i=0; i<nz; i++) bornes[i+1] =
alpha[i] +
beta[i] ;
185 assert (tmap != 0x0) ;
196 for (
int lz=0; lz<nz0; lz++)
197 for (
int ir=0; ir<
mg->
get_nr(lz); ir++)
200 if (
mg->
get_nt(0) != 1) xmetr = xmetr /
sqrt(
double(2)) ;
203 const Mtbl& erre = this->
r ;
210 for (
int lz=0; lz<nz0; lz++) {
212 double r1 = erre(lz, 0, 0, nr-1) ;
213 double rm1 = erre(lz, 0, 0, 0) ;
220 a3.
set(lz) = (x1 - a1(lz)) / (r1 * r1);
224 double r0 = erre(lz, 0, 0, i0) ;
226 double p1 = (r1 - rm1)*(r1 - r0) ;
227 double pm1 = (r0 - rm1)*(r1 - rm1) ;
228 double p0 = (r0 - rm1)*(r1 - r0) ;
229 a1.
set(lz) = xm1*r1*r0/pm1 + x1*rm1*r0/p1 - x0*rm1*r1/p0 ;
230 a2.
set(lz) = x0*(rm1 + r1)/p0 - xm1*(r1 + r0)/pm1
232 a3.
set(lz) = xm1/pm1+x1/p1-x0/p0 ;
235 for (
int k=0; k<
mg->
get_np(lz)+2; k++)
236 for (
int j=0; j<
mg->
get_nt(lz); j++)
237 for (
int i=0; i<nr; i++)
239 (a2(lz) + erre(lz, 0, 0, i)*a3(lz)) ;
242 double offset =
max(diff) ;
243 a1.
set(lz) += offset ;
249 sigma += (dt*dt)*(source + reste) ;
251 sigma += (0.5*dt*dt)*mime*fjm1_local.
laplacian() ;
254 sigma += (dt*dt) * source ;
256 sigma += (0.5*dt*dt)*fjm1_local.
laplacian() ;
265 int l_q, m_q, baser ;
267 for (
int lz=0; lz<nz-1; lz++) {
270 for (
int k=0; k<np+2; k++)
271 for (
int j=0; j<nt; j++) {
273 if ((nullite_plm(j, nt, k, np, base) == 1) && (l_q+dl >= l_min) ) {
274 for (
int i=0; i<
mg->
get_nr(lz); i++) {
276 0.5*dt*dt*dl*(2*l_q + dl +1)
294 for (
int i=0; i<nz; i++) {
295 coeff->
set(1,i) = a1(i) ;
296 coeff->
set(2,i) = a2(i) ;
297 coeff->
set(3,i) = a3(i) ;
298 coeff->
set(4,i) = 0. ;
299 coeff->
set(5,i) = 0. ;
300 coeff->
set(6,i) = 0. ;
301 coeff->
set(7,i) = 0. ;
302 coeff->
set(8,i) = 0. ;
303 coeff->
set(9,i) = 0. ;
310 double R = this->
val_r(nz0-1, 1., 0., 0.) ;
326 tbc3 =
new Tbl(np2,nt) ;
334 phijm1 =
new Tbl(np2,nt) ;
335 phij =
new Tbl(np2,nt) ;
363 if (bound3.
get_etat() == ETATZERO) {
369 if (nz0>1) bound3.
annule(0,nz0-2) ;
379 for (
int k=0; k<np2; k++)
380 for (
int j=0; j<nt; j++) {
382 for (
int i=0; i<nr; i++)
383 val += (*bound3.
c_cf)(nz0-1,k,j,i) ;
384 tbc3->
set(k,j) = val ;
397 if (nz0>1) souphi.
annule(0,nz0-2) ;
401 bool zero = (souphi.
get_etat() == ETATZERO) ;
409 int l_s, m_s, base_r ;
411 int dl = (par.
get_n_int() > 1) ? -1 : 0 ;
412 for (
int k=0; k<np2; k++) {
413 for (
int j=0; j<nt; j++) {
414 donne_lm(nz, nz0-1, j, k, souphi.
base, m_s, l_s, base_r) ;
419 double multi = 8*R*R + dt*dt*(6+3*l_s*(l_s+1)) + 12*R*dt ;
420 val = ( 16*R*R*(*phij)(k,j) -
421 (multi-24*R*dt)*(*phijm1)(k,j)
423 phijm1->
set(k,j) = (*phij)(k,j) ;
424 phij->
set(k,j) = val ;
431 if (bound3.
get_etat() == ETATZERO) *tbc3 = 0 ;
433 if (nz0 > 1) bound3.
annule(0,nz0-2) ;
437 for (
int k=0; k<np2; k++)
438 for (
int j=0; j<nt; j++) {
440 for (
int i=0; i<nr; i++)
441 val += (*bound3.
c_cf)(nz0-1,k,j,i) ;
442 tbc3->
set(k,j) = val + 2*R*dt*(*phij)(k,j);
448 cout <<
"ERROR: Map_af::dalembert" << endl ;
449 cout <<
"The boundary condition par.get_int(0) = "<< par.
get_int(0)
450 <<
" is unknown!" << endl ;
464 assert(sourva.
get_etat() == ETATQCQ) ;
Bases of the spectral expansions.
void give_quant_numbers(int, int, int, int &, int &, int &) const
Computes the various quantum numbers and 1d radial base.
void ylm()
The basis is transformed as with a transformation to basis.
void dsdx()
The basis is transformed as with a operation.
double * beta
Array (size: mg->nzone ) of the values of in each domain.
virtual void dalembert(Param &par, Scalar &fJp1, const Scalar &fJ, const Scalar &fJm1, const Scalar &source) const
Performs one time-step integration of the d'Alembert scalar equation.
double * alpha
Array (size: mg->nzone ) of the values of in each domain.
virtual double val_r(int l, double xi, double theta, double pphi) const
Returns the value of the radial coordinate r for a given in a given domain.
Coord r
r coordinate centered on the grid
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
const Mg3d * mg
Pointer on the multi-grid Mgd3 on which this 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.
const Mg3d * get_radial() const
Returns the pointer on the associated radial grid.
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.
Base_val std_base_scal() const
Returns the standard spectral bases for a scalar.
Coefficients storage for the multi-domain spectral method.
Tbl & set(int l)
Read/write of the Tbl containing the coefficients in a given domain.
void annule_hard()
Sets the Mtbl_cf to zero in a hard way.
double val_out_bound_jk(int l, int j, int k) const
Computes the angular coefficient of index j,k of the field represented by *this at by means of the s...
Tensor & get_tensor_mod(int position=0) const
Returns the reference of a modifiable Tensor stored in the list.
void add_map(const Map &mi, int position=0)
Adds the address of a new Map to the list.
const int & get_int(int position=0) const
Returns the reference of a int stored in the list.
const double & get_double(int position=0) const
Returns the reference of a double stored in the list.
int get_n_tbl_mod() const
Returns the number of modifiable Tbl 's addresses in the list.
int get_n_int_mod() const
Returns the number of modifiable int 's addresses in the list.
Tbl & get_tbl_mod(int position=0) const
Returns the reference of a modifiable Tbl stored in the list.
const Map & get_map(int position=0) const
Returns the reference of a Map stored in the list.
int get_n_tensor_mod() const
Returns the number of modifiable Tensor 's addresses in the list.
void add_double_mod(double &x, int position=0)
Adds the address of a new modifiable double to the list.
void add_tbl_mod(Tbl &ti, int position=0)
Adds the address of a new modifiable Tbl to the list.
int get_n_double() const
Returns the number of stored double 's addresses.
int get_n_int() const
Returns the number of stored int 's addresses.
int & get_int_mod(int position=0) const
Returns the reference of a modifiable int stored in the list.
double & get_double_mod(int position=0) const
Returns the reference of a stored modifiable double .
Tensor field of valence 0 (or component of a tensorial field).
int get_dzpuis() const
Returns dzpuis.
const Scalar & laplacian(int ced_mult_r=4) const
Returns the Laplacian of *this.
virtual void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
virtual void std_spectral_base()
Sets the spectral bases of the Valeur va to the standard ones for a scalar field.
double val_grid_point(int l, int k, int j, int i) const
Returns the value of the field at a specified grid point.
void div_r()
Division by r everywhere; dzpuis is not changed.
virtual void set_etat_zero()
Sets the logical state to ETATZERO (zero).
const Scalar & dsdr() const
Returns of *this .
Tbl & set_domain(int l)
Read/write of the value in a given domain.
Valeur & set_spectral_va()
Returns va (read/write version)
const Valeur & get_spectral_va() const
Returns va (read only version)
void annule_hard()
Sets the Scalar to zero in a hard way.
int get_etat() const
Returns the logical state ETATNONDEF (undefined), ETATZERO (null) or ETATQCQ (ordinary).
const Base_val & get_spectral_base() const
Returns the spectral bases of the Valeur va
void set_dzpuis(int)
Modifies the dzpuis flag.
double & set_grid_point(int l, int k, int j, int i)
Setting the value of the field at a given grid point.
const Tbl & domain(int l) const
Read-only of the value in a given domain.
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 & set(int i)
Read/write of a particular element (index i) (1D case)
Values and coefficients of a (real-value) function.
void set_etat_cf_qcq()
Sets the logical state to ETATQCQ (ordinary state) for values in the configuration space (Mtbl_cf c_c...
int get_etat() const
Returns the logical state.
void set_base(const Base_val &)
Sets the bases for spectral expansions (member base )
void ylm()
Computes the coefficients of *this.
Mtbl * c
Values of the function at the points of the multi-grid
void annule(int l)
Sets the Valeur to zero in a given domain.
void coef_i() const
Computes the physical value of *this.
Mtbl_cf * c_cf
Coefficients of the spectral expansion of the function.
void coef() const
Computes the coeffcients of *this.
void ylm_i()
Inverse of ylm()
Base_val base
Bases on which the spectral expansion is performed.
void set_base_t(int base_t)
Sets the expansion basis for functions in all domains.
Cmp sqrt(const Cmp &)
Square root.
Tbl max(const Cmp &)
Maximum values of a Cmp in each domain.
#define T_LEG_PP
fct. de Legendre associees paires avec m pair
const Map & get_mp() const
Returns the mapping.
void annule_domain(int l)
Sets the Tensor to zero in a given domain.