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