5 #include "utilitaires.h"
24 const Mg3d* mgrid = (*map).get_mg();
30 const Coord& rr = (*map).r ;
35 Scalar source_coq = source ;
52 int l_q, m_q, base_r ;
55 for (
int k=0; k < np; k++)
56 for (
int j=0; j<nt; j++) {
58 if (nullite_plm(j, nt, k, np, base) == 1) {
59 for (
int ii=0 ; ii<nr ; ii++){
60 sol_hom1.
set(lz, k, j, ii) = 0 ;
61 sol_part.
set(lz, k, j, ii) = 0 ;
70 double alpha = (*map).get_alpha()[lz] ;
71 double beta = (*map).get_beta()[lz] ;
72 double ech = beta / alpha ;
73 for (
int k=0; k < np; k++)
74 for (
int j=0; j<nt; j++) {
76 if (nullite_plm(j, nt, k, np, base) == 1) {
105 ope = mx2dx2 + 2*ech*mxdx2 + ech*ech*mdx2 + 2*(mxdx + ech*mdx)
106 - l_q*(l_q+1)*mid -2*(l_q+1)*mid - ((mx2dx2 + 2*ech*mxdx2 + ech*ech*mdx2) -(fit)*(mx2dx2 + 2*ech*mxdx2 + ech*ech*mdx2) + (fit)*beta*(mx2dx2 + 2*ech*mxdx2 + ech*ech*mdx2)
107 + (fit)*alpha*(mx3dx2 + 2*ech*mx2dx2 + ech*ech*mxdx2) + fit2*alpha*alpha*(mx4dx2 + 2*ech*mx3dx2 + ech*ech*mx2dx2)
108 + fit2*alpha*(beta-1.)*(mx3dx2 + 2*ech*mx2dx2 + ech*ech*mxdx2) + fit2*alpha*(beta- rrr.
val_grid_point(1,0,0, nr-1))*(mx3dx2 + 2*ech*mx2dx2 + ech*ech*mxdx2)+ fit2*(beta-1.)*(beta- rrr.
val_grid_point(1,0,0,nr -1))*(mx2dx2 + 2*ech*mxdx2 + ech*ech*mdx2));
113 for (
int col=0; col<nr; col++)
114 ope.
set(nr-1, col) = 0 ;
115 ope.
set(nr-1, 0) = 1 ;
120 for (
int i=0; i<nr; i++)
127 for (
int i=0; i<nr; i++)
128 sol_part.
set(1, k, j, i) = sol(i) ;
135 for (
int i=0; i<nr; i++)
136 sol_hom1.
set(1, k, j, i) = sol(i) ;
150 double alpha = (*map).get_alpha()[lz] ;
151 double beta = (*map).get_beta()[lz] ;
152 double ech = beta / alpha ;
154 for (
int k=0; k < np; k++)
155 for (
int j=0; j<nt; j++) {
157 if (nullite_plm(j, nt, k, np, base) == 1) {
170 ope = mx2dx2 + 2*ech*mxdx2 + ech*ech*mdx2 + 2*(mxdx + ech*mdx)
171 - l_q*(l_q+1)*mid -2*(l_q+1)*mid - fit0d2*(mx2dx2 + 2*ech*mxdx2 + ech*ech*mdx2) + (fit1d2)*(rrr.
val_grid_point(lz, 0, 0, 0))*(mx2dx2 + 2*ech*mxdx2 + ech*ech*mdx2) - (fit1d2)*beta*(mx2dx2 + 2*ech*mxdx2 + ech*ech*mdx2) - (fit1d2)*alpha*(mx3dx2 + 2*ech*mx2dx2 + ech*ech*mxdx2) ;
178 for (
int col=0; col<nr; col++)
179 ope.
set(nr-1, col) = 0 ;
180 ope.
set(nr-1, 0) = 1 ;
181 for (
int col=0; col<nr; col++) {
182 ope.
set(nr-2, col) = 0 ;
184 ope.
set(nr-2, 1) = 1 ;
188 for (
int i=0; i<nr; i++)
195 for (
int i=0; i<nr; i++)
196 sol_part.
set(lz, k, j, i) = sol(i) ;
201 for (
int i=0; i<nr; i++)
202 sol_hom1.
set(lz, k, j, i) = sol(i) ;
207 for (
int i=0; i<nr; i++)
208 sol_hom2.
set(lz, k, j, i) = sol(i) ;
219 double alpha = (*map).get_alpha()[lz] ;
220 double beta = (*map).get_beta()[lz] ;
221 double ech = beta / alpha ;
223 for (
int k=0; k < np; k++)
224 for (
int j=0; j<nt; j++) {
226 if (nullite_plm(j, nt, k, np, base) == 1) {
239 ope = mx2dx2 + 2*ech*mxdx2 + ech*ech*mdx2 + 2*(mxdx + ech*mdx)
240 - l_q*(l_q+1)*mid -2*(l_q+1)*mid - fit0d3*(mx2dx2 + 2*ech*mxdx2 + ech*ech*mdx2) + (fit1d3)*(rrr.
val_grid_point(lz, 0, 0, 0))*(mx2dx2 + 2*ech*mxdx2 + ech*ech*mdx2) - (fit1d3)*beta*(mx2dx2 + 2*ech*mxdx2 + ech*ech*mdx2) - (fit1d3)*alpha*(mx3dx2 + 2*ech*mx2dx2 + ech*ech*mxdx2) ;
246 for (
int col=0; col<nr; col++)
247 ope.
set(nr-1, col) = 0 ;
248 ope.
set(nr-1, 0) = 1 ;
249 for (
int col=0; col<nr; col++) {
250 ope.
set(nr-2, col) = 0 ;
252 ope.
set(nr-2, 1) = 1 ;
256 for (
int i=0; i<nr; i++)
263 for (
int i=0; i<nr; i++)
264 sol_part.
set(lz, k, j, i) = sol(i) ;
269 for (
int i=0; i<nr; i++)
270 sol_hom1.
set(lz, k, j, i) = sol(i) ;
275 for (
int i=0; i<nr; i++)
276 sol_hom2.
set(lz, k, j, i) = sol(i) ;
287 for (
int lz=4; lz<nz-1; lz++) {
288 double ech = (*map).get_beta()[lz] / (*map).get_alpha()[lz] ;
289 for (
int k=0; k < np; k++)
290 for (
int j=0; j<nt; j++) {
292 if (nullite_plm(j, nt, k, np, base) == 1) {
303 ope = mx2dx2 + 2*ech*mxdx2 + ech*ech*mdx2 + 2*(mxdx + ech*mdx)
304 - l_q*(l_q+1)*mid -2*(l_q+1)*mid;
306 for (
int col=0; col<nr; col++)
307 ope.
set(nr-1, col) = 0 ;
308 ope.
set(nr-1, 0) = 1 ;
309 for (
int col=0; col<nr; col++) {
310 ope.
set(nr-2, col) = 0 ;
312 ope.
set(nr-2, 1) = 1 ;
316 for (
int i=0; i<nr; i++)
323 for (
int i=0; i<nr; i++)
324 sol_part.
set(lz, k, j, i) = sol(i) ;
329 for (
int i=0; i<nr; i++)
330 sol_hom1.
set(lz, k, j, i) = sol(i) ;
335 for (
int i=0; i<nr; i++)
336 sol_hom2.
set(lz, k, j, i) = sol(i) ;
343 double alpha = (*map).get_alpha()[lz] ;
344 for (
int k=0; k < np; k++)
345 for (
int j=0; j<nt; j++) {
347 if (nullite_plm(j, nt, k, np, base) == 1) {
354 ope = (mdx2 - l_q*(l_q+1)*ms2 -2*(l_q+1)*ms2)/(alpha*alpha) ;
356 for (
int i=0; i<nr; i++)
357 ope.
set(nr-1, i) = 0 ;
358 ope.
set(nr-1, 0) = 1 ;
360 for (
int i=0; i<nr; i++) {
361 ope.
set(nr-2, i) = 1 ;
365 for (
int i=0; i<nr; i++) {
366 ope.
set(nr-3, i) = i*i ;
373 for (
int i=0; i<nr; i++)
375 if (l_q>0) rhs.
set(nr-3) = 0 ;
381 for (
int i=0; i<nr; i++)
382 sol_part.
set(lz, k, j, i) = sol(i) ;
387 for (
int i=0; i<nr; i++)
388 sol_hom2.
set(lz, k, j, i) = sol(i) ;
411 for (
int k=0; k < np; k++)
412 for (
int j=0; j<nt; j++) {
414 if (nullite_plm(j, nt, k, np, base) == 1) {
415 Matrice systeme(2*nz-4, 2*nz-4) ;
424 double alpha = (*map).get_alpha()[1] ;
435 for (
int lz=2; lz<nz-1; lz++) {
436 alpha = (*map).get_alpha()[lz] ;
460 alpha = (*map).get_alpha()[nz-1] ;
481 for (
int i=0; i<(*mgrid).get_nr(1); i++)
482 sol_coef.
set(1, k, j, i) = sol_part(1, k, j, i)
483 +coef(indice)*sol_hom1(1, k, j, i) ;
487 for (
int lz=2; lz<nz-1; lz++) {
488 for (
int i=0; i<(*mgrid).get_nr(lz); i++)
489 sol_coef.
set(lz, k, j, i) = sol_part(lz, k, j, i)
490 +coef(indice)*sol_hom1(lz, k, j, i)
491 +coef(indice+1)*sol_hom2(lz, k, j, i) ;
494 for (
int i=0; i<(*mgrid).get_nr(nz-1); i++)
495 sol_coef.
set(nz-1, k, j, i) = sol_part(nz-1, k, j, i)
496 +coef(indice)*sol_hom2(nz-1, k, j, i) ;
510 for (
int lz=0; lz<nz; lz++) {
511 for (
int i=0; i<(*mgrid).get_nr(lz); i++)
513 sol_coef.
set(lz,0,0,i) = 0.;
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.
Active physical coordinates and mapping derivatives.
Class for the elementary differential operator (see the base class Diff ).
virtual const Matrice & get_matrice() const
Returns the matrix associated with the operator.
Class for the elementary differential operator (see the base class Diff ).
virtual const Matrice & get_matrice() const
Returns the matrix associated with the operator.
Class for the elementary differential operator Identity (see the base class Diff ).
Class for the elementary differential operator division by (see the base class Diff ).
virtual const Matrice & get_matrice() const
Returns the matrix associated with the operator.
Class for the elementary differential operator (see the base class Diff ).
virtual const Matrice & get_matrice() const
Returns the matrix associated with the operator.
Class for the elementary differential operator (see the base class Diff ).
virtual const Matrice & get_matrice() const
Returns the matrix associated with the operator.
Class for the elementary differential operator (see the base class Diff ).
virtual const Matrice & get_matrice() const
Returns the matrix associated with the operator.
Class for the elementary differential operator (see the base class Diff ).
virtual const Matrice & get_matrice() const
Returns the matrix associated with the operator.
Class for the elementary differential operator (see the base class Diff ).
virtual const Matrice & get_matrice() const
Returns the matrix associated with the operator.
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
double & set(int j, int i)
Read/write of a particuliar element.
Tbl inverse(const Tbl &sec_membre) const
Solves the linear system represented by the matrix.
void annule_hard()
Sets the logical state to ETATQCQ (undefined state).
void set_lu() const
Calculate the LU-representation, assuming the band-storage has been done.
Coefficients storage for the multi-domain spectral method.
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...
double val_in_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...
Tbl & set(int l)
Read/write of the Tbl containing the coefficients in a given domain.
Tensor field of valence 0 (or component of a tensorial field).
const Base_val & get_spectral_base() const
Returns the spectral bases of the Valeur va
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 annule_hard()
Sets the Scalar to zero in a hard way.
const Valeur & get_spectral_va() const
Returns va (read only version)
Valeur & set_spectral_va()
Returns va (read/write version)
void mult_r()
Multiplication by r everywhere; dzpuis is not changed.
void annule_hard()
Sets the Tbl to zero in a hard way.
double & set(int i)
Read/write of a particular element (index i) (1D case)
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
Mtbl_cf * c_cf
Coefficients of the spectral expansion of the function.
Base_val base
Bases on which the spectral expansion is performed.
const Map & get_mp() const
Returns the mapping.
void annule_domain(int l)
Sets the Tensor to zero in a given domain.
void tensorellipticCt(Scalar source, Scalar &resu, double fitd1, double fit2d1, double fit0d2, double fit1d2, double fit0d3, double fit1d3)