21char param_elliptic_C[] =
"$Header: /cvsroot/Lorene/C++/Source/Param_elliptic/param_elliptic.C,v 1.21 2014/10/13 08:53:37 j_novak Exp $" ;
110#include "ope_elementary.h"
111#include "param_elliptic.h"
128 cout <<
"Unknown mapping in Param_elliptic" <<
endl ;
134double Param_elliptic::get_alpha(
int l)
const {
138 return mp_af->get_alpha()[
l] ;
144 cout <<
"Unknown mapping in Param_elliptic" <<
endl ;
150double Param_elliptic::get_beta(
int l)
const {
154 return mp_af->get_beta()[l] ;
157 return mp_log->get_beta(l) ;
160 cout <<
"Unknown mapping in Param_elliptic" << endl ;
166int Param_elliptic::get_type(
int l)
const {
173 return mp_log->get_type(l) ;
176 cout <<
"Unknown mapping in Param_elliptic" << endl ;
184 done_F (
so.get_mp().get_mg()->get_nzone(),
185 so.get_mp().get_mg()->get_np(0) + 1,
186 so.get_mp().get_mg()->get_nt(0)),
187 done_G (
so.get_mp().get_mg()->get_nzone()),
188 val_F_plus (
so.get_mp().get_mg()->get_nzone(),
189 so.get_mp().get_mg()->get_np(0) + 1,
190 so.get_mp().get_mg()->get_nt(0)),
191 val_F_minus (
so.get_mp().get_mg()->get_nzone(),
192 so.get_mp().get_mg()->get_np(0) + 1,
193 so.get_mp().get_mg()->get_nt(0)),
194 val_dF_plus (
so.get_mp().get_mg()->get_nzone(),
195 so.get_mp().get_mg()->get_np(0) + 1,
196 so.get_mp().get_mg()->get_nt(0)),
197 val_dF_minus (
so.get_mp().get_mg()->get_nzone(),
198 so.get_mp().get_mg()->get_np(0) + 1,
199 so.get_mp().get_mg()->get_nt(0)),
200 val_G_plus (
so.get_mp().get_mg()->get_nzone()),
201 val_G_minus (
so.get_mp().get_mg()->get_nzone()),
202 val_dG_plus (
so.get_mp().get_mg()->get_nzone()),
203 val_dG_minus (
so.get_mp().get_mg()->get_nzone())
210 auxi.set_spectral_va().coef() ;
211 auxi.set_spectral_va().ylm() ;
214 int dzpuis = (
auxi.dz_nonzero() ?
auxi.get_dzpuis() : 4) ;
222 cout <<
"Param_elliptic not yet defined on this type of mapping" <<
endl ;
237 int nz =
get_mp().get_mg()->get_nzone() ;
239 for (
int l=0 ;
l<nz ;
l++)
245 for (
int l=0 ;
l<
nbr ;
l++)
252 for (
int l=0 ;
l<nz ;
l++) {
254 nr =
get_mp().get_mg()->get_nr(
l) ;
256 for (
int k=0 ;
k<
get_mp().get_mg()->get_np(
l)+1 ;
k++)
257 for (
int j=0 ;
j<
get_mp().get_mg()->get_nt(
l) ;
j++) {
259 so.get_spectral_va().base.give_quant_numbers
266 get_beta(
l), l_quant, dzpuis) ;
269 if (
mp_log->get_type(
l) == AFFINE)
272 get_beta(
l), l_quant, dzpuis) ;
279 cout <<
"Unknown mapping in Param_elliptic::Param_elliptic"
313 for (
int l=0 ;
l<
get_mp().get_mg()->get_nzone() ;
l++)
316 for (
int i=0 ;
i<
nbr ;
i++)
328 for (
int l=0 ;
l<
get_mp().get_mg()->get_nzone() ;
l++) {
330 np =
get_mp().get_mg()->get_np(
l) ;
331 nt =
get_mp().get_mg()->get_nt(
l) ;
333 for (
int k=0 ;
k<np+1 ;
k++)
334 for (
int j=0 ;
j<nt ;
j++) {
345 source.set_spectral_va().coef() ;
346 source.set_spectral_va().ylm() ;
348 int nz =
get_mp().get_mg()->get_nzone() ;
357 for (
int l=0 ;
l<nz ;
l++) {
359 nr =
get_mp().get_mg()->get_nr(
l) ;
361 for (
int k=0 ;
k<
get_mp().get_mg()->get_np(
l)+1 ;
k++)
362 for (
int j=0 ;
j<
get_mp().get_mg()->get_nt(
l) ;
j++) {
366 source.get_spectral_va().base.give_quant_numbers
369 get_beta(
l) , masse) ;
379 cout <<
"Operator not define with LOG mapping..." <<
endl ;
383 for (
int l=0 ;
l<nz ;
l++) {
385 nr =
get_mp().get_mg()->get_nr(
l) ;
387 for (
int k=0 ;
k<
get_mp().get_mg()->get_np(
l)+1 ;
k++)
388 for (
int j=0 ;
j<
get_mp().get_mg()->get_nt(
l) ;
j++) {
392 source.get_spectral_va().base.give_quant_numbers
395 get_alpha(
l), get_beta(
l), masse) ;
405 cout <<
"Unkown mapping in set_helmhotz_minus" <<
endl ;
413 source.set_spectral_va().coef() ;
414 source.set_spectral_va().ylm() ;
416 int nz =
get_mp().get_mg()->get_nzone() ;
425 for (
int l=0 ;
l<nz ;
l++) {
427 nr =
get_mp().get_mg()->get_nr(
l) ;
429 for (
int k=0 ;
k<
get_mp().get_mg()->get_np(
l)+1 ;
k++)
430 for (
int j=0 ;
j<
get_mp().get_mg()->get_nt(
l) ;
j++) {
437 source.get_spectral_va().base.give_quant_numbers
440 get_beta(
l) , masse) ;
451 cout <<
"Operator not define with LOG mapping..." <<
endl ;
455 for (
int l=0 ;
l<nz ;
l++) {
457 nr =
get_mp().get_mg()->get_nr(
l) ;
459 for (
int k=0 ;
k<
get_mp().get_mg()->get_np(
l)+1 ;
k++)
460 for (
int j=0 ;
j<
get_mp().get_mg()->get_nt(
l) ;
j++) {
467 source.get_spectral_va().base.give_quant_numbers
470 get_alpha(
l), get_beta(
l), masse) ;
481 cout <<
"Unkown mapping in set_helmhotz_minus" <<
endl ;
490 cout <<
"set_poisson_vect_r only defined for an affine mapping..." <<
endl ;
495 int nz =
get_mp().get_mg()->get_nzone() ;
500 for (
int l=0 ;
l<nz ;
l++) {
502 nr =
get_mp().get_mg()->get_nr(
l) ;
504 for (
int k=0 ;
k<
get_mp().get_mg()->get_np(
l)+1 ;
k++)
505 for (
int j=0 ;
j<
get_mp().get_mg()->get_nt(
l) ;
j++) {
535 int nz =
get_mp().get_mg()->get_nzone() ;
538 for (
int l=0 ;
l<nz ;
l++) {
540 bool ced = (
get_mp().get_mg()->get_type_r(
l) == UNSURR ) ;
541 for (
int k=0 ;
k<
get_mp().get_mg()->get_np(
l)+1 ;
k++)
542 for (
int j=0 ;
j<
get_mp().get_mg()->get_nt(
l) ;
j++) {
563 cout <<
"set_poisson_tens_rr only defined for an affine mapping..."
569 int nz =
get_mp().get_mg()->get_nzone() ;
574 for (
int l=0 ;
l<nz ;
l++) {
576 nr =
get_mp().get_mg()->get_nr(
l) ;
578 for (
int k=0 ;
k<
get_mp().get_mg()->get_np(
l)+1 ;
k++)
579 for (
int j=0 ;
j<
get_mp().get_mg()->get_nt(
l) ;
j++) {
605 int nz =
get_mp().get_mg()->get_nzone() ;
613 for (
int l=0 ;
l<nz ;
l++) {
615 nr =
get_mp().get_mg()->get_nr(
l) ;
617 for (
int k=0 ;
k<
get_mp().get_mg()->get_np(
l)+1 ;
k++)
618 for (
int j=0 ;
j<
get_mp().get_mg()->get_nt(
l) ;
j++) {
626 get_beta(
l),
a, b, c) ;
636 cout <<
"Operator not define with LOG mapping..." <<
endl ;
640 for (
int l=0 ;
l<nz ;
l++) {
642 nr =
get_mp().get_mg()->get_nr(
l) ;
644 for (
int k=0 ;
k<
get_mp().get_mg()->get_np(
l)+1 ;
k++)
645 for (
int j=0 ;
j<
get_mp().get_mg()->get_nt(
l) ;
j++) {
653 get_beta(
l),
a, b, c) ;
663 cout <<
"Unkown mapping in set_sec_order_r2" <<
endl ;
672 cout <<
"set_sec_order only defined for a log mapping" <<
endl ;
677 int nz =
get_mp().get_mg()->get_nzone() ;
682 for (
int l=0 ;
l<nz ;
l++) {
684 nr =
get_mp().get_mg()->get_nr(
l) ;
686 for (
int k=0 ;
k<
get_mp().get_mg()->get_np(
l)+1 ;
k++)
687 for (
int j=0 ;
j<
get_mp().get_mg()->get_nt(
l) ;
j++) {
696 get_beta(
l),
a, b, c) ;
707 source.set_spectral_va().coef() ;
708 source.set_spectral_va().ylm() ;
710 int nz =
get_mp().get_mg()->get_nzone() ;
714 int dzpuis =
source.get_dzpuis() ;
720 for (
int l=0 ;
l<nz ;
l++) {
721 int dz = (
l==nz-1) ? dzpuis : 0 ;
722 nr =
get_mp().get_mg()->get_nr(
l) ;
724 for (
int k=0 ;
k<
get_mp().get_mg()->get_np(
l)+1 ;
k++)
725 for (
int j=0 ;
j<
get_mp().get_mg()->get_nt(
l) ;
j++) {
729 source.get_spectral_va().base.give_quant_numbers
732 get_beta(
l), l_quant,
dz) ;
742 cout <<
"Operator not define with LOG mapping..." <<
endl ;
746 for (
int l=0 ;
l<nz ;
l++) {
747 int dz = (
l==nz-1) ? dzpuis : 0 ;
748 nr =
get_mp().get_mg()->get_nr(
l) ;
750 for (
int k=0 ;
k<
get_mp().get_mg()->get_np(
l)+1 ;
k++)
751 for (
int j=0 ;
j<
get_mp().get_mg()->get_nt(
l) ;
j++) {
755 source.get_spectral_va().base.give_quant_numbers
758 get_alpha(
l), get_beta(
l), l_quant,
dz) ;
768 cout <<
"Unkown mapping in set_ope_vorton" <<
endl ;
776 assert (
so.get_etat() != ETATNONDEF) ;
785 assert (
so.get_etat() != ETATNONDEF) ;
Bases of the spectral expansions.
Time evolution with partial storage (*** under development ***).
void annule_hard()
Sets the Itbl to zero in a hard way.
Logarithmic radial mapping.
Base class for pure radial mappings.
Basic class for elementary elliptic operators.
int get_base_r() const
Returns base_r}.
virtual void inc_l_quant()=0
Increases the quatum number l by one unit.
Class for the Helmholtz operator ( ).
Class for the Helmholtz operator (m > 0).
Class for the operator of the rr component of the divergence-free tensor Poisson equation.
Class for the operator of the r component of the vector Poisson equation.
Class for the operator of the Poisson equation (i.e.
Class for operator of the type .
Class for operator of the type .
Class for the operator appearing for the vortons.
Scalar var_G
Multiplicative variable change that must be sphericaly symetric !
Itbl done_F
Stores what has been computed for F.
Tbl val_dG_minus
Values of the derivative of G at the inner boundaries of the various domains.
void set_sec_order_r2(int zone, double a, double b, double c)
Set the operator to in one domain (only in the shells).
void inc_l_quant(int zone)
Increases the quantum number l in the domain zone .
void set_helmholtz_minus(int zone, double mas, Scalar &so)
Set the operator to in one domain (not in the nucleus).
Tbl val_dF_minus
Values of the derivative of F at the inner boundaries of the various domains.
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...
const Map_log * mp_log
The mapping if log type.
Param_elliptic(const Scalar &)
Standard constructor from a Scalar.
const Map_radial & get_mp() const
Returns the mapping.
Tbl val_F_minus
Values of F at the inner boundaries of the various domains.
Tbl val_F_plus
Values of F at the outer boundaries of the various domains.
~Param_elliptic()
Destructor.
void set_sec_order(int zone, double a, double b, double c)
Set the operator to in one domain (only in the shells).
Tbl val_G_minus
Values of G at the inner boundaries of the various domains.
Tbl val_dG_plus
Values of the derivative of G at the outer boundaries of the various domains.
Ope_elementary ** operateurs
Array on the elementary operators.
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.
const Map_af * mp_af
The mapping, if affine.
void set_poisson_tens_rr(int zone)
Sets the operator to in all domains, for only.
void set_variable_F(const Scalar &)
Changes the variable function F.
void set_helmholtz_plus(int zone, double mas, Scalar &so)
Set the operator to in one domain (only in the shells).
void set_variable_G(const Scalar &)
Changes the variable function G.
int type_map
Type of mapping either MAP_AFF or MAP_LOG.
Tbl val_dF_plus
Values of the derivative of F at the outer boundaries of the various domains.
Tbl val_G_plus
Values of G at the outer boundaries of the various domains.
void set_ope_vorton(int zone, Scalar &so)
Set the operator to in one domain (not implemented in the nucleus).
Itbl done_G
Stores what has been computed for G.
Scalar var_F
Additive variable change function.
Tensor field of valence 0 (or component of a tensorial field).
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.
Valeur & set_spectral_va()
Returns va (read/write version)
void annule_hard()
Sets the Scalar to zero in a hard way.
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
void set_base(const Base_val &)
Sets the bases for spectral expansions (member base )
#define R_CHEBI
base de Cheb. impaire (rare) seulement