28 char vector_poisson_C[] =
"$Header: /cvsroot/Lorene/C++/Source/Tensor/vector_poisson.C,v 1.24 2014/10/13 08:53:45 j_novak Exp $" ;
126 #include "param_elliptic.h" 133 bool nullite = true ;
134 for (
int i=0; i<3; i++) {
135 assert(
cmp[i]->check_dzpuis(4)) ;
136 if (
cmp[i]->get_etat() != ETATZERO) nullite = false ;
138 assert ((method>=0) && (method<7)) ;
150 if (fabs(lambda+1) < 1.e-8)
166 if (fabs(lambda+1) < 1.e-8)
169 divf = (
potential(met_f) / (lambda + 1)) ;
181 if (div_lnot0.
get_etat() != ETATZERO) {
184 for (
int lz=0; lz<nz; lz++) {
188 if (va_div.
c_cf->operator()(lz).get_etat() != ETATZERO)
189 for (
int k=0; k<np+1; k++)
190 for (
int j=0; j<nt; j++) {
191 int l_q, m_q, base_r ;
192 if (nullite_plm(j, nt, k, np, va_div.
base) == 1) {
193 donne_lm(nz, lz, j, k, va_div.
base, m_q, l_q, base_r) ;
195 for (
int i=0; i<nr; i++)
196 va_div.
c_cf->
set(lz, k, j, i) = 0. ;
211 source_r += *(
cmp[0]) - lambda*divf.
dsdr() ;
213 if (source_r.
get_etat() != ETATZERO) {
219 for (
int lz=0; lz<nz; lz++)
230 source_eta -= 2*f_r ;
241 for (
int i=0; i<3; i++) {
242 source_p.set(i) =
Cmp(*
cmp[i]) ;
250 Tenseur resu_p(source_p.poisson_vect(lambda, vect_auxi, scal_auxi)) ;
253 for (
int i=1; i<=3; i++)
254 resu.
set(i) = resu_p(i-1) ;
263 for (
int i=0; i<3; i++) {
264 source_p.set(i) =
Cmp(*
cmp[i]) ;
270 Tenseur resu_p(source_p.poisson_vect_oohara(lambda, scal_auxi)) ;
273 for (
int i=1; i<=3; i++)
274 resu.
set(i) = resu_p(i-1) ;
281 if (fabs(lambda+1) < 1.e-8)
284 divf = (
potential(met_f) / (lambda + 1)) ;
296 if (div_tmp.
get_etat() != ETATZERO) {
298 for (
int lz=0; lz<nz; lz++) {
302 if (va_div.
c_cf->operator()(lz).get_etat() != ETATZERO)
303 for (
int k=0; k<np+1; k++)
304 for (
int j=0; j<nt; j++) {
305 int l_q, m_q, base_r ;
306 if (nullite_plm(j, nt, k, np, va_div.
base) == 1) {
307 donne_lm(nz, lz, j, k, va_div.
base, m_q, l_q, base_r) ;
309 for (
int i=0; i<nr; i++)
310 va_div.
c_cf->
set(lz, k, j, i) = 0. ;
325 source_r += *(
cmp[0]) - lambda*divf.
dsdr() ;
327 if (source_r.
get_etat() != ETATZERO) {
333 for (
int lz=0; lz<nz; lz++)
346 source_eta = source_eta.
dsdt() ;
348 source_eta = (source_eta +
cmp[2]->
stdsdp()).poisson_angu() ;
366 for (
int lz=0; lz<nz; lz++) {
371 if (va_eta.
c_cf->operator()(lz).get_etat() != ETATZERO)
372 for (
int k=0; k<np+1; k++)
373 for (
int j=0; j<nt; j++) {
374 int l_q, m_q, base_r ;
375 if (nullite_plm(j, nt, k, np, va_eta.
base) == 1) {
376 donne_lm(nz, lz, j, k, va_eta.
base, m_q, l_q, base_r) ;
378 for (
int i=0; i<nr; i++) {
379 if (va_div.
c_cf->operator()(lz).get_etat() != ETATZERO)
381 -= (lambda + 2. / double(ced ? -l_q : (l_q+1) ))
382 * va_div.
c_cf->operator()(lz, k, j, i) ;
383 if (va_fsr.
c_cf->operator()(lz).get_etat() != ETATZERO) {
385 += 2. / double(ced ? -l_q : (l_q+1) )
386 * va_dfr.
c_cf->operator()(lz, k, j, i) ;
388 -= 2.*( ced ? double(l_q+2)/double(l_q)
389 : double(l_q-1)/double(l_q+1) )
390 * va_fsr.
c_cf->
set(lz, k, j, i) ;
396 if (va_eta.
c != 0x0) {
406 for (
int lz=0; lz<nz; lz++)
417 if (fabs(lambda+1) < 1.e-8)
420 divf = (
potential(met_f) / (lambda + 1)) ;
432 if (div_lnot0.
get_etat() != ETATZERO) {
435 for (
int lz=0; lz<nz; lz++) {
439 if (va_div.
c_cf->operator()(lz).get_etat() != ETATZERO)
440 for (
int k=0; k<np+1; k++)
441 for (
int j=0; j<nt; j++) {
442 int l_q, m_q, base_r ;
443 if (nullite_plm(j, nt, k, np, va_div.
base) == 1) {
444 donne_lm(nz, lz, j, k, va_div.
base, m_q, l_q, base_r) ;
446 for (
int i=0; i<nr; i++)
447 va_div.
c_cf->
set(lz, k, j, i) = 0. ;
462 source_r += *(
cmp[0]) - lambda*divf.
dsdr() ;
464 if (source_r.
get_etat() != ETATZERO) {
471 for (
int lz=0; lz<nz; lz++)
481 source_eta -= (lambda+2.)*divf ;
487 tmp = (1.+lambda)*divf ;
490 source_eta = source_eta.
primr() ;
499 poisson_block(lambda, resu) ;
504 cout <<
"Vector::poisson : unexpected type of method !" << endl
505 <<
" method = " << method << endl ;
523 assert ((tspher != 0x0) || (tcart != 0x0)) ;
527 assert (tcart == 0x0) ;
532 assert (tspher == 0x0) ;
536 return (
poisson(lambda, *met_f, method) );
546 for (
int i=0; i<3; i++)
547 assert(
cmp[i]->check_dzpuis(4)) ;
549 assert ((method==0) || (method==2)) ;
558 int nitermax = par.
get_int(0) ;
574 if (fabs(lambda+1) < 1.e-8)
586 par_free.
add_int(nitermax, 0) ;
603 for (
int i=0; i<3; i++) {
604 source_p.set(i) =
Cmp(*
cmp[i]) ;
610 for (
int i=0; i<3 ;i++){
611 vect_auxi.set(i) = 0. ;
620 source_p.poisson_vect(lambda, par, resu_p, vect_auxi, scal_auxi) ;
623 for (
int i=1; i<=3; i++)
624 resu.
set(i) = resu_p(i-1) ;
630 cout <<
"Vector::poisson : unexpected type of method !" << endl
631 <<
" method = " << method << endl ;
void set_poisson_vect_eta(int zone)
Sets the operator to be a regular elliptic operator to solve for the component of the vector Poisson...
Mtbl_cf * c_cf
Coefficients of the spectral expansion of the function.
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
void add_int(const int &n, int position=0)
Adds the address of a new int to the list.
void ylm_i()
Inverse of ylm()
const Scalar & lapang() const
Returns the angular Laplacian of *this , where .
int get_np(int l) const
Returns the number of points in the azimuthal direction ( ) in domain no. l.
void coef() const
Computes the coeffcients of *this.
virtual void set_etat_zero()
Sets the logical state to ETATZERO (zero).
const Vector_divfree & div_free(const Metric &) const
Returns the div-free vector in the Helmholtz decomposition.
void set_std_base()
Set the standard spectal basis of decomposition for each component.
const Base_vect_spher & get_bvect_spher() const
Returns the orthonormal vectorial basis associated with the coordinates of the mapping.
Tbl & set(int l)
Read/write of the Tbl containing the coefficients in a given domain.
void ylm()
Computes the coefficients of *this.
const Scalar & dsdt() const
Returns of *this .
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Flat metric for tensor calculation.
Tensor field of valence 0 (or component of a tensorial field).
Scalar poisson() const
Solves the scalar Poisson equation with *this as a source.
Scalar primr(bool null_infty=true) const
Computes the radial primitive which vanishes for .
void mult_sint()
Multiplication by .
void annule_hard()
Sets the Cmp to zero in a hard way.
const Base_vect * triad
Vectorial basis (triad) with respect to which the tensor components are defined.
Values and coefficients of a (real-value) function.
int get_etat() const
Returns the logical state ETATNONDEF (undefined), ETATZERO (null) or ETATQCQ (ordinary).
virtual void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
const Vector & derive_con(const Metric &gam) const
Returns the "contravariant" derivative of *this with respect to some metric , by raising the index of...
Tensor field of valence 1.
Scalar poisson_angu(double lambda=0) const
Solves the (generalized) angular Poisson equation with *this as source.
const Metric_flat & flat_met_cart() const
Returns the flat metric associated with the Cartesian coordinates and with components expressed in th...
void set_dzpuis(int)
Modifies the dzpuis flag.
virtual void dec_dzpuis(int dec=1)
Decreases by dec units the value of dzpuis and changes accordingly the values in the compactified ext...
void annule_hard()
Sets the Scalar to zero in a hard way.
Cmp & get_cmp_mod(int position=0) const
Returns the reference of a modifiable Cmp stored in the list.
Cmp & set()
Read/write for a scalar (see also operator=(const Cmp&) ).
const int & get_int(int position=0) const
Returns the reference of a int stored in the list.
void change_triad(const Base_vect &new_triad)
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.
virtual void inc_dzpuis(int inc=1)
Increases by inc units the value of dzpuis and changes accordingly the values of the Scalar in the co...
Vector poisson(double lambda, int method=6) const
Solves the vector Poisson equation with *this as a source.
Base_val base
Bases on which the spectral expansion is performed.
void set_vr_eta_mu(const Scalar &vr_i, const Scalar &eta_i, const Scalar &mu_i)
Defines the components through potentials and (see members p_eta and p_mu ), as well as the compon...
Scalar ** cmp
Array of size n_comp of pointers onto the components.
Scalar sol_elliptic(Param_elliptic ¶ms) const
Resolution of a general elliptic equation, putting zero at infinity.
Mtbl * c
Values of the function at the points of the multi-grid.
int get_nzone() const
Returns the number of domains.
int & get_int_mod(int position=0) const
Returns the reference of a modifiable int stored in the list.
const Scalar & stdsdp() const
Returns of *this .
Vector_divfree poisson() const
Computes the solution of a vectorial Poisson equation with *this as a source: .
This class contains the parameters needed to call the general elliptic solver.
Cartesian vectorial bases (triads).
Spherical orthonormal vectorial bases (triads).
int get_nr(int l) const
Returns the number of points in the radial direction ( ) in domain no. l.
const Base_vect_cart & get_bvect_cart() const
Returns the Cartesian basis associated with the coordinates (x,y,z) of the mapping, i.e.
const Scalar & potential(const Metric &) const
Returns the potential in the Helmholtz decomposition.
void add_double(const double &x, int position=0)
Adds the the address of a new double to the list.
const Scalar & dsdr() const
Returns of *this .
virtual const Scalar & mu() const
Gives the field such that the angular components of the vector are written: .
void div_sint()
Division by .
virtual void set_etat_zero()
Sets the logical state of all components to ETATZERO (zero state).
const double & get_double(int position=0) const
Returns the reference of a double stored in the list.
int get_nt(int l) const
Returns the number of points in the co-latitude direction ( ) in domain no. l.
Valeur & set_spectral_va()
Returns va (read/write version)
Scalar & set(int)
Read/write access to a component.
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
void add_cmp_mod(Cmp &ti, int position=0)
Adds the address of a new modifiable Cmp to the list.
int get_type_r(int l) const
Returns the type of sampling in the radial direction in domain no.
const Map *const mp
Mapping on which the numerical values at the grid points are defined.
void div_r_dzpuis(int ced_mult_r)
Division by r everywhere but with the output flag dzpuis set to ced_mult_r .
virtual void dec_dzpuis(int dec=1)
Decreases by dec units the value of dzpuis and changes accordingly the values of the Scalar in the co...
const Metric_flat & flat_met_spher() const
Returns the flat metric associated with the spherical coordinates and with components expressed in th...
void mult_r_dzpuis(int ced_mult_r)
Multiplication by r everywhere but with the output flag dzpuis set to ced_mult_r .
Tensor handling *** DEPRECATED : use class Tensor instead ***.
void add_int_mod(int &n, int position=0)
Adds the address of a new modifiable int to the list.
void set_poisson_vect_r(int zone, bool only_l_zero=false)
Sets the operator to in all domains, for ; and to in all domains otherwise.