LORENE
cmp_pde_frontiere.C
1 /*
2  * Copyright (c) 2000-2001 Philippe Grandclement
3  *
4  * This file is part of LORENE.
5  *
6  * LORENE is free software; you can redistribute it and/or modify
7  * it under the terms of the GNU General Public License as published by
8  * the Free Software Foundation; either version 2 of the License, or
9  * (at your option) any later version.
10  *
11  * LORENE is distributed in the hope that it will be useful,
12  * but WITHOUT ANY WARRANTY; without even the implied warranty of
13  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  * GNU General Public License for more details.
15  *
16  * You should have received a copy of the GNU General Public License
17  * along with LORENE; if not, write to the Free Software
18  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
19  *
20  */
21 
22 
23 
24 
25 /*
26  * $Id: cmp_pde_frontiere.C,v 1.8 2016/12/05 16:17:49 j_novak Exp $
27  * $Log: cmp_pde_frontiere.C,v $
28  * Revision 1.8 2016/12/05 16:17:49 j_novak
29  * Suppression of some global variables (file names, loch, ...) to prevent redefinitions
30  *
31  * Revision 1.7 2014/10/13 08:52:48 j_novak
32  * Lorene classes and functions now belong to the namespace Lorene.
33  *
34  * Revision 1.6 2005/02/18 13:14:08 j_novak
35  * Changing of malloc/free to new/delete + suppression of some unused variables
36  * (trying to avoid compilation warnings).
37  *
38  * Revision 1.5 2004/11/23 12:49:58 f_limousin
39  * Intoduce function poisson_dir_neu(...) to solve a scalar poisson
40  * equation with a mixed boundary condition (Dirichlet + Neumann).
41  *
42  * Revision 1.4 2004/05/20 07:04:02 k_taniguchi
43  * Recovery of "->get_angu()" in the assertion of Map_af::poisson_frontiere
44  * because "limite" is the boundary value.
45  *
46  * Revision 1.3 2004/03/31 11:18:42 f_limousin
47  * Methods Map_et::poisson_interne, Map_af::poisson_interne and
48  * Cmp::poisson_neumann_interne have been implemented to solve the
49  * continuity equation for strange stars.
50  *
51  * Revision 1.2 2003/10/03 15:58:45 j_novak
52  * Cleaning of some headers
53  *
54  * Revision 1.1.1.1 2001/11/20 15:19:27 e_gourgoulhon
55  * LORENE
56  *
57  * Revision 2.6 2000/05/22 16:07:03 phil
58  * *** empty log message ***
59  *
60  * Revision 2.5 2000/05/22 16:03:48 phil
61  * ajout du cas dzpuis = 3
62  *
63  * Revision 2.4 2000/04/27 15:18:27 phil
64  * ajout des procedures relatives a la resolution dans une seule zone avec deux conditions limites.
65  *
66  * Revision 2.3 2000/03/31 15:59:54 phil
67  * gestion des cas ou la source est nulle.
68  *
69  * Revision 2.2 2000/03/20 13:08:53 phil
70  * *** empty log message ***
71  *
72  * Revision 2.1 2000/03/17 17:33:05 phil
73  * *** empty log message ***
74  *
75  * Revision 2.0 2000/03/17 17:25:08 phil
76  * *** empty log message ***
77  *
78  *
79  * $Header: /cvsroot/Lorene/C++/Source/Cmp/cmp_pde_frontiere.C,v 1.8 2016/12/05 16:17:49 j_novak Exp $
80  *
81  */
82 
83 // Header Lorene:
84 #include "scalar.h"
85 #include "param.h"
86 #include "cmp.h"
87 
88 namespace Lorene {
89 Mtbl_cf sol_poisson_frontiere(const Map_af&, const Mtbl_cf&, const Mtbl_cf&,
90  int, int, int, double = 0.,
91  double = 0.) ;
92 
93 Mtbl_cf sol_poisson_frontiere_double (const Map_af&, const Mtbl_cf&, const Mtbl_cf&,
94  const Mtbl_cf&, int) ;
95 
96 Mtbl_cf sol_poisson_interne(const Map_af&, const Mtbl_cf&, const Mtbl_cf&) ;
97 
98 Cmp Cmp::poisson_dirichlet (const Valeur& limite, int num_front) const {
99 
100  Cmp resu(*mp) ;
101  mp->poisson_frontiere (*this, limite, 1, num_front, resu) ;
102  return resu ;
103 }
104 
105 
106 Cmp Cmp::poisson_neumann (const Valeur& limite, int num_front) const {
107 
108  Cmp resu(*mp) ;
109  mp->poisson_frontiere (*this, limite, 2, num_front, resu) ;
110  return resu ;
111 }
112 
114  Param& par, Cmp& resu) const {
115 
116  mp->poisson_interne (*this, limite, par, resu) ;
117  return resu ;
118 }
119 
120 Cmp Cmp::poisson_frontiere_double (const Valeur& lim_func, const Valeur& lim_der,
121  int num_zone) const {
122  Cmp resu(*mp) ;
123  mp->poisson_frontiere_double (*this, lim_func, lim_der, num_zone, resu) ;
124  return resu ;
125 }
126 
127 void Map_et::poisson_frontiere(const Cmp&, const Valeur&, int, int, Cmp&, double, double) const {
128  cout << "Procedure non implantee ! " << endl ;
129  abort() ;
130 }
131 
132 void Map_et::poisson_frontiere_double (const Cmp&, const Valeur&, const Valeur&,
133  int, Cmp&) const {
134  cout << "Procedure non implantee ! " << endl ;
135  abort() ;
136 }
137 
138 void Map_af::poisson_frontiere(const Cmp& source, const Valeur& limite,
139  int type_raccord, int num_front,
140  Cmp& pot, double fact_dir, double fact_neu) const {
141 
142  assert(source.get_etat() != ETATNONDEF) ;
143  assert(source.get_mp()->get_mg() == mg) ;
144  assert(pot.get_mp()->get_mg() == mg) ;
145 
146  assert( source.check_dzpuis(2) || source.check_dzpuis(4)
147  || source.check_dzpuis(3)) ;
148  assert ((type_raccord == 1) || (type_raccord==2)|| (type_raccord==3)) ;
149  int dzpuis ;
150 
151  if (source.dz_nonzero()){
152  dzpuis = source.get_dzpuis() ;
153  }
154  else{
155  dzpuis = 4 ;
156  }
157 
158  // Spherical harmonic expansion of the source
159  // ------------------------------------------
160 
161  Valeur sourva = source.va ;
162  sourva.coef() ;
163  sourva.ylm() ;
164 
165  // Pour gerer le cas ou source est dans ETAT_ZERO...
166  Mtbl_cf so_cf (sourva.get_mg(), sourva.base) ;
167  if (sourva.get_etat() == ETATZERO) {
168  so_cf.set_etat_zero() ;
169  }
170  else
171  so_cf = *sourva.c_cf ;
172 
173 
174  Valeur conditions (limite) ;
175  conditions.coef() ;
176  conditions.ylm() ; // spherical harmonic transforms
177 
178  // Pour gerer le cas ou condition est dans ETAT_ZERO...
179  Mtbl_cf auxiliaire (conditions.get_mg(), conditions.base) ;
180  if (conditions.get_etat() == ETATZERO)
181  auxiliaire.set_etat_zero() ;
182  else
183  auxiliaire = *conditions.c_cf ;
184 
185  Mtbl_cf resu (sourva.get_mg(), sourva.base) ;
186 
187  if (type_raccord == 3){
188 
189  // Call to the Mtbl_cf version
190  // ---------------------------
191  resu = sol_poisson_frontiere(*this, so_cf, auxiliaire,
192  type_raccord, num_front, dzpuis,
193  fact_dir, fact_neu) ;
194  }
195  else{
196  resu = sol_poisson_frontiere(*this, so_cf, auxiliaire,
197  type_raccord, num_front, dzpuis) ;
198  }
199 
200  // Final result returned as a Cmp
201  // ------------------------------
202 
203  pot.set_etat_zero() ; // to call Cmp::del_t().
204  pot.set_etat_qcq() ;
205  pot.va = resu ;
206  (pot.va).ylm_i() ; // On repasse en base standard.
207  pot.set_dzpuis(0) ;
208 }
209 
210 
211 void Map_af::poisson_frontiere_double (const Cmp& source, const Valeur& lim_func,
212  const Valeur& lim_der, int num_zone, Cmp& pot) const {
213 
214  assert(source.get_etat() != ETATNONDEF) ;
215  assert(source.get_mp()->get_mg() == mg) ;
216  assert(pot.get_mp()->get_mg() == mg) ;
217  assert (source.get_mp()->get_mg()->get_angu() == lim_func.get_mg()) ;
218  assert (source.get_mp()->get_mg()->get_angu() == lim_der.get_mg()) ;
219 
220  // Spherical harmonic expansion of the source
221  // ------------------------------------------
222 
223  Valeur sourva = source.va ;
224  sourva.coef() ;
225  sourva.ylm() ;
226 
227  // Pour gerer le cas ou source est dans ETAT_ZERO...
228  Mtbl_cf so_cf (sourva.get_mg(), sourva.base) ;
229  if (sourva.get_etat() == ETATZERO) {
230  so_cf.set_etat_zero() ;
231  }
232  else
233  so_cf = *sourva.c_cf ;
234 
235 
236  Valeur cond_func (lim_func) ;
237  cond_func.coef() ;
238  cond_func.ylm() ; // spherical harmonic transforms
239 
240  // Pour gerer le cas ou condition est dans ETAT_ZERO...
241  Mtbl_cf auxi_func (cond_func.get_mg(), cond_func.base) ;
242  if (cond_func.get_etat() == ETATZERO)
243  auxi_func.set_etat_zero() ;
244  else
245  auxi_func = *cond_func.c_cf ;
246 
247  Valeur cond_der (lim_der) ;
248  cond_der.coef() ;
249  cond_der.ylm() ; // spherical harmonic transforms
250 
251  // Pour gerer le cas ou condition est dans ETAT_ZERO...
252  Mtbl_cf auxi_der (cond_der.get_mg(), cond_der.base) ;
253  if (cond_der.get_etat() == ETATZERO)
254  auxi_der.set_etat_zero() ;
255  else
256  auxi_der = *cond_der.c_cf ;
257 
258 
259 
260  // Call to the Mtbl_cf version
261  // ---------------------------
262 
263  Mtbl_cf resu = sol_poisson_frontiere_double (*this, so_cf, auxi_func,
264  auxi_der, num_zone) ;
265 
266  // Final result returned as a Cmp
267  // ------------------------------
268 
269  pot.set_etat_zero() ; // to call Cmp::del_t().
270  pot.set_etat_qcq() ;
271  pot.va = resu ;
272  (pot.va).ylm_i() ; // On repasse en base standard.
273 }
274 
275 
276 
277 void Map_et::poisson_interne(const Cmp& source, const Valeur& limite,
278  Param& par, Cmp& uu) const {
279 
280  assert(source.get_etat() != ETATNONDEF) ;
281  assert(source.get_mp() == this) ;
282 
283  assert(uu.get_mp() == this) ;
284 
285  int nz = mg->get_nzone() ;
286 
287  //-------------------------------
288  // Computation of the prefactor a ---> Cmp apre
289  //-------------------------------
290 
291  Mtbl unjj = 1 + srdrdt*srdrdt + srstdrdp*srstdrdp ;
292 
293  Mtbl apre1(*mg) ;
294  apre1.set_etat_qcq() ;
295  for (int l=0; l<nz; l++) {
296  *(apre1.t[l]) = alpha[l]*alpha[l] ;
297  }
298 
299  apre1 = apre1 * dxdr * dxdr * unjj ;
300 
301  Cmp apre(*this) ;
302  apre = apre1 ;
303 
304  Tbl amax0 = max(apre1) ; // maximum values in each domain
305 
306  // The maximum values of a in each domain are put in a Mtbl
307  Mtbl amax1(*mg) ;
308  amax1.set_etat_qcq() ;
309  for (int l=0; l<nz; l++) {
310  *(amax1.t[l]) = amax0(l) ;
311  }
312 
313  Cmp amax(*this) ;
314  amax = amax1 ;
315 
316 
317  //-------------------
318  // Initializations
319  //-------------------
320 
321  int nitermax = par.get_int() ;
322  int& niter = par.get_int_mod() ;
323  double lambda = par.get_double() ;
324  double unmlambda = 1. - lambda ;
325  double precis = par.get_double(1) ;
326 
327  Cmp& ssj = par.get_cmp_mod() ;
328 
329  Cmp ssjm1 = ssj ;
330  Cmp ssjm2 = ssjm1 ;
331 
332  Valeur& vuu = uu.va ;
333 
334  Valeur vuujm1(*mg) ;
335  if (uu.get_etat() == ETATZERO) {
336  vuujm1 = 1 ; // to take relative differences
337  vuujm1.set_base( vuu.base ) ;
338  }
339  else{
340  vuujm1 = vuu ;
341  }
342 
343  // Affine mapping for the Laplacian-tilde
344 
345  Map_af mpaff(*this) ;
346  Param par_nul ;
347 
348  cout << "Map_et::poisson : relat. diff. u^J <-> u^{J-1} : " << endl ;
349 
350 //==========================================================================
351 //==========================================================================
352 // Start of iteration
353 //==========================================================================
354 //==========================================================================
355 
356  Tbl tdiff(nz) ;
357  niter = 0 ;
358 
359  do {
360 
361  //====================================================================
362  // Computation of R(u) (the result is put in uu)
363  //====================================================================
364 
365 
366  //-----------------------
367  // First operations on uu
368  //-----------------------
369 
370  Valeur duudx = (uu.va).dsdx() ; // d/dx
371 
372  const Valeur& d2uudtdx = duudx.dsdt() ; // d^2/dxdtheta
373 
374  const Valeur& std2uudpdx = duudx.stdsdp() ; // 1/sin(theta) d^2/dxdphi
375 
376  //------------------
377  // Angular Laplacian
378  //------------------
379 
380  Valeur sxlapang = uu.va ;
381 
382  sxlapang.ylm() ;
383 
384  sxlapang = sxlapang.lapang() ;
385 
386  sxlapang = sxlapang.sx() ; // Lap_ang(uu) /x in the nucleus
387 
388  //---------------------------------------------------------------
389  // Computation of
390  // [ 2 /(dRdx) ( A - 1 ) duu/dx + 1/R (B - 1) Lap_ang(uu) ] / x
391  //
392  // with A = 1/(dRdx) R/(x+beta_l/alpha_l) unjj
393  // B = [1/(dRdx) R/(x+beta_l/alpha_l)]^2 unjj
394  //
395  // The result is put in uu (via vuu)
396  //---------------------------------------------------------------
397 
398  Valeur varduudx = duudx ;
399 
400  uu.set_etat_qcq() ;
401 
402  Base_val sauve_base = varduudx.base ;
403 
404  vuu = 2. * dxdr * ( rsxdxdr * unjj - 1.) * varduudx
405  + ( rsxdxdr*rsxdxdr * unjj - 1.) * xsr * sxlapang ;
406 
407  vuu.set_base(sauve_base) ;
408 
409  vuu = vuu.sx() ;
410 
411  //---------------------------------------
412  // Computation of R(u)
413  //
414  // The result is put in uu (via vuu)
415  //---------------------------------------
416 
417 
418  sauve_base = vuu.base ;
419 
420  vuu = xsr * vuu
421  + 2. * dxdr * ( sr2drdt * d2uudtdx
422  + sr2stdrdp * std2uudpdx ) ;
423 
424  vuu += dxdr * ( lapr_tp + dxdr * (
425  dxdr* unjj * d2rdx2
426  - 2. * ( sr2drdt * d2rdtdx + sr2stdrdp * sstd2rdpdx ) )
427  ) * duudx ;
428 
429  vuu.set_base(sauve_base) ;
430 
431  //====================================================================
432  // Computation of the effective source s^J of the ``affine''
433  // Poisson equation
434  //====================================================================
435 
436  ssj = lambda * ssjm1 + unmlambda * ssjm2 ;
437 
438  ssj = ( source + uu + (amax - apre) * ssj ) / amax ;
439 
440  (ssj.va).set_base((source.va).base) ;
441 
442  //====================================================================
443  // Resolution of the ``affine'' Poisson equation
444  //====================================================================
445 
446  mpaff.poisson_interne(ssj, limite, par_nul, uu) ;
447 
448  tdiff = diffrel(vuu, vuujm1) ;
449 
450  cout << " step " << niter << " : " ;
451  cout << tdiff(0) << " " ;
452 
453  cout << endl ;
454 
455  //=================================
456  // Updates for the next iteration
457  //=================================
458 
459  ssjm2 = ssjm1 ;
460  ssjm1 = ssj ;
461  vuujm1 = vuu ;
462 
463  niter++ ;
464 
465  } // End of iteration
466  while ( (tdiff(0) > precis) && (niter < nitermax) ) ;
467 
468 //==========================================================================
469 //==========================================================================
470 // End of iteration
471 //==========================================================================
472 //==========================================================================
473 
474 }
475 
476 
477 void Map_af::poisson_interne(const Cmp& source, const Valeur& limite,
478  Param& , Cmp& pot) const {
479 
480  assert(source.get_etat() != ETATNONDEF) ;
481  assert(source.get_mp()->get_mg() == mg) ;
482  assert(pot.get_mp()->get_mg() == mg) ;
483  assert (source.get_mp()->get_mg()->get_angu() == limite.get_mg()) ;
484 
485  // Spherical harmonic expansion of the source
486  // ------------------------------------------
487 
488  Valeur sourva = source.va ;
489  sourva.coef() ;
490  sourva.ylm() ;
491 
492  // Pour gerer le cas ou source est dans ETAT_ZERO...
493  Mtbl_cf so_cf (sourva.get_mg(), sourva.base) ;
494  if (sourva.get_etat() == ETATZERO) {
495  so_cf.set_etat_zero() ;
496  }
497  else
498  so_cf = *sourva.c_cf ;
499 
500  Valeur conditions (limite) ;
501  conditions.coef() ;
502  conditions.ylm() ; // spherical harmonic transforms
503 
504 
505  // Pour gerer le cas ou condition est dans ETAT_ZERO...
506  Mtbl_cf auxiliaire (conditions.get_mg(), conditions.base) ;
507  if (conditions.get_etat() == ETATZERO)
508  auxiliaire.set_etat_zero() ;
509  else
510  auxiliaire = *conditions.c_cf ;
511 
512  // Call to the Mtbl_cf version
513  // ---------------------------
514 
515  Mtbl_cf resu = sol_poisson_interne(*this, so_cf, auxiliaire) ;
516 
517  // Final result returned as a Cmp
518  // ------------------------------
519 
520  pot.set_etat_zero() ; // to call Cmp::del_t().
521  pot.set_etat_qcq() ;
522  pot.va = resu ;
523  (pot.va).ylm_i() ; // On repasse en base standard.
524  pot.set_dzpuis(0) ;
525 }
526 
527 }
const Map * get_mp() const
Returns the mapping.
Definition: cmp.h:901
Coord d2rdx2
in the nucleus and in the non-compactified shells; \ in the compactified outer domain.
Definition: map.h:1634
const Valeur & dsdt() const
Returns of *this.
Definition: valeur_dsdt.C:115
Coord sr2stdrdp
in the nucleus and in the non-compactified shells; \ in the compactified outer domain.
Definition: map.h:1623
Mtbl_cf * c_cf
Coefficients of the spectral expansion of the function.
Definition: valeur.h:312
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition: cmp.h:446
const Valeur & lapang() const
Returns the angular Laplacian of *this.
Definition: valeur_lapang.C:75
void coef() const
Computes the coeffcients of *this.
Definition: valeur_coef.C:151
Multi-domain array.
Definition: mtbl.h:118
double * alpha
Array (size: mg->nzone ) of the values of in each domain.
Definition: map.h:2776
Lorene prototypes.
Definition: app_hor.h:67
void ylm()
Computes the coefficients of *this.
Definition: valeur_ylm.C:141
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:777
int get_etat() const
Returns the logical state.
Definition: cmp.h:899
Coord sr2drdt
in the nucleus and in the non-compactified shells; \ in the compactified outer domain.
Definition: map.h:1615
virtual void poisson_frontiere_double(const Cmp &source, const Valeur &lim_func, const Valeur &lim_der, int num_zone, Cmp &pot) const
Solver of the Poisson equation with boundary condition for the affine mapping case, cases with boundary conditions of both Dirichlet and Neumann type (no condition imposed at infinity).
const Valeur & sx() const
Returns (r -sampling = RARE ) \ Id (r sampling = FIN ) \ (r -sampling = UNSURR ) ...
Definition: valeur_sx.C:113
Values and coefficients of a (real-value) function.
Definition: valeur.h:297
void set_base(const Base_val &)
Sets the bases for spectral expansions (member base )
Definition: valeur.C:813
Tbl diffrel(const Cmp &a, const Cmp &b)
Relative difference between two Cmp (norme version).
Definition: cmp_math.C:507
Coord srstdrdp
in the nucleus and in the non-compactified shells; \ in the compactified outer domain.
Definition: map.h:1607
const Mg3d * get_mg() const
Returns the Mg3d on which the this is defined.
Definition: valeur.h:763
Cmp & get_cmp_mod(int position=0) const
Returns the reference of a modifiable Cmp stored in the list.
Definition: param.C:1052
virtual void poisson_interne(const Cmp &source, const Valeur &limite, Param &par, Cmp &pot) const =0
Computes the solution of a Poisson equation in the shell, imposing a boundary condition at the surfac...
int get_etat() const
Returns the logical state.
Definition: valeur.h:760
virtual void poisson_frontiere(const Cmp &, const Valeur &, int, int, Cmp &, double=0., double=0.) const
Not yet implemented.
void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition: cmp.C:292
const int & get_int(int position=0) const
Returns the reference of a int stored in the list.
Definition: param.C:295
const Mg3d * get_angu() const
Returns the pointer on the associated angular grid.
Definition: mg3d.C:604
Coord dxdr
in the nucleus and in the non-compactified shells; \ in the compactified outer domain.
Definition: map.h:1575
Base_val base
Bases on which the spectral expansion is performed.
Definition: valeur.h:315
virtual void poisson_frontiere(const Cmp &, const Valeur &, int, int, Cmp &, double=0., double=0.) const
Solver of the Poisson equation with boundary condition for the affine mapping case.
bool dz_nonzero() const
Returns true if the last domain is compactified and *this is not zero in this domain.
Definition: cmp.C:663
Parameter storage.
Definition: param.h:125
const Valeur & stdsdp() const
Returns of *this.
Definition: valeur_stdsdp.C:63
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:465
int & get_int_mod(int position=0) const
Returns the reference of a modifiable int stored in the list.
Definition: param.C:433
Cmp poisson_dirichlet(const Valeur &limite, int num) const
Is identicall to Cmp::poisson() .
Tbl max(const Cmp &)
Maximum values of a Cmp in each domain.
Definition: cmp_math.C:438
Coord sstd2rdpdx
in the nucleus and in the non-compactified shells; \ in the compactified outer domain.
Definition: map.h:1663
Coord xsr
in the nucleus; \ 1/R in the non-compactified shells; \ in the compactified outer domain...
Definition: map.h:1564
virtual void poisson_interne(const Cmp &source, const Valeur &limite, Param &par, Cmp &pot) const
Computes the solution of a Poisson equation in the shell .
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: cmp.C:307
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: mtbl.C:302
Coord rsxdxdr
in the nucleus; \ in the shells; \ in the outermost compactified domain.
Definition: map.h:2852
void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition: mtbl_cf.C:291
const Mg3d * mg
Pointer on the multi-grid Mgd3 on which this is defined.
Definition: map.h:688
virtual void poisson_frontiere(const Cmp &source, const Valeur &limite, int raccord, int num_front, Cmp &pot, double=0., double=0.) const =0
Computes the solution of a Poisson equation from the domain num_front+1 .
Bases of the spectral expansions.
Definition: base_val.h:325
Affine radial mapping.
Definition: map.h:2042
const Map * mp
Reference mapping.
Definition: cmp.h:451
int get_dzpuis() const
Returns dzpuis.
Definition: cmp.h:903
Coord lapr_tp
in the nucleus and in the non-compactified shells; \ in the compactified outer domain.
Definition: map.h:1646
Coefficients storage for the multi-domain spectral method.
Definition: mtbl_cf.h:196
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...
Definition: cmp.C:718
virtual void poisson_interne(const Cmp &source, const Valeur &limite, Param &par, Cmp &pot) const
Computes the solution of a Poisson equation in the shell, imposing a boundary condition at the surfac...
Cmp poisson_neumann(const Valeur &, int) const
Idem as Cmp::poisson_dirichlet , the boundary condition being on the radial derivative of the solutio...
void set_dzpuis(int)
Set a value to dzpuis.
Definition: cmp.C:657
const double & get_double(int position=0) const
Returns the reference of a double stored in the list.
Definition: param.C:364
Basic array class.
Definition: tbl.h:164
Tbl ** t
Array (size nzone ) of pointers on the Tbl &#39;s.
Definition: mtbl.h:132
Coord srdrdt
in the nucleus and in the non-compactified shells; \ in the compactified outer domain.
Definition: map.h:1599
Valeur va
The numerical value of the Cmp.
Definition: cmp.h:464
Cmp poisson_neumann_interne(const Valeur &, Param &par, Cmp &resu) const
Idem as Cmp::poisson_neumann , the boundary condition is on the radial derivative of the solution...
Coord d2rdtdx
in the nucleus and in the non-compactified shells; \ in the compactified outer domain.
Definition: map.h:1655