LORENE
scalar_pde.C
1 /*
2  * Methods of the class Scalar for various partial differential equations
3  *
4  * See file scalar.h for documentation.
5  */
6 
7 /*
8  * Copyright (c) 2003-2005 Eric Gourgoulhon & Jerome Novak
9  * Copyright (c) 2004 Philippe Grandclement
10  *
11  * Copyright (c) 1999-2001 Eric Gourgoulhon (for preceding class Cmp)
12  * Copyright (c) 1999-2001 Philippe Grandclement (for preceding class Cmp)
13  * Copyright (c) 2000-2001 Jerome Novak (for preceding class Cmp)
14  *
15  * This file is part of LORENE.
16  *
17  * LORENE is free software; you can redistribute it and/or modify
18  * it under the terms of the GNU General Public License as published by
19  * the Free Software Foundation; either version 2 of the License, or
20  * (at your option) any later version.
21  *
22  * LORENE is distributed in the hope that it will be useful,
23  * but WITHOUT ANY WARRANTY; without even the implied warranty of
24  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
25  * GNU General Public License for more details.
26  *
27  * You should have received a copy of the GNU General Public License
28  * along with LORENE; if not, write to the Free Software
29  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
30  *
31  */
32 
33 
34 char scalar_pde_C[] = "$Header: /cvsroot/Lorene/C++/Source/Tensor/Scalar/scalar_pde.C,v 1.20 2014/10/13 08:53:46 j_novak Exp $" ;
35 
36 /*
37  * $Id: scalar_pde.C,v 1.20 2014/10/13 08:53:46 j_novak Exp $
38  * $Log: scalar_pde.C,v $
39  * Revision 1.20 2014/10/13 08:53:46 j_novak
40  * Lorene classes and functions now belong to the namespace Lorene.
41  *
42  * Revision 1.19 2007/05/06 10:48:15 p_grandclement
43  * Modification of a few operators for the vorton project
44  *
45  * Revision 1.18 2007/01/16 15:10:00 n_vasset
46  * New function sol_elliptic_boundary, with Scalar on mono domain
47  * angular grid as boundary
48  *
49  * Revision 1.17 2005/11/30 11:09:09 p_grandclement
50  * Changes for the Bin_ns_bh project
51  *
52  * Revision 1.16 2005/08/26 14:02:41 p_grandclement
53  * Modification of the elliptic solver that matches with an oscillatory exterior solution
54  * small correction in Poisson tau also...
55  *
56  * Revision 1.15 2005/08/25 12:14:10 p_grandclement
57  * Addition of a new method to solve the scalar Poisson equation, based on a multi-domain Tau-method
58  *
59  * Revision 1.14 2005/06/09 08:00:10 f_limousin
60  * Implement a new function sol_elliptic_boundary() and
61  * Vector::poisson_boundary(...) which solve the vectorial poisson
62  * equation (method 6) with an inner boundary condition.
63  *
64  * Revision 1.13 2005/04/04 21:34:44 e_gourgoulhon
65  * Added argument lambda to method poisson_angu
66  * to deal with the generalized angular Poisson equation:
67  * Lap_ang u + lambda u = source.
68  *
69  * Revision 1.12 2004/08/24 09:14:52 p_grandclement
70  * Addition of some new operators, like Poisson in 2d... It now requieres the
71  * GSL library to work.
72  *
73  * Also, the way a variable change is stored by a Param_elliptic is changed and
74  * no longer uses Change_var but rather 2 Scalars. The codes using that feature
75  * will requiere some modification. (It should concern only the ones about monopoles)
76  *
77  * Revision 1.11 2004/06/22 08:50:00 p_grandclement
78  * Addition of everything needed for using the logarithmic mapping
79  *
80  * Revision 1.10 2004/05/25 14:30:48 f_limousin
81  * Minor modif.
82  *
83  * Revision 1.9 2004/03/17 15:58:50 p_grandclement
84  * Slight modification of sol_elliptic_no_zec
85  *
86  * Revision 1.8 2004/03/01 09:57:04 j_novak
87  * the wave equation is solved with Scalars. It now accepts a grid with a
88  * compactified external domain, which the solver ignores and where it copies
89  * the values of the field from one time-step to the next.
90  *
91  * Revision 1.7 2004/02/11 09:47:46 p_grandclement
92  * Addition of a new elliptic solver, matching with the homogeneous solution
93  * at the outer shell and not solving in the external domain (more details
94  * coming soon ; check your local Lorene dealer...)
95  *
96  * Revision 1.6 2004/01/28 16:59:14 p_grandclement
97  * *** empty log message ***
98  *
99  * Revision 1.5 2004/01/28 16:46:24 p_grandclement
100  * Addition of the sol_elliptic_fixe_der_zero stuff
101  *
102  * Revision 1.4 2004/01/14 10:11:51 f_limousin
103  * Corrected bug in poisson with parameters.
104  *
105  * Revision 1.3 2003/12/11 14:48:51 p_grandclement
106  * Addition of ALL (and that is a lot !) the files needed for the general elliptic solver ... UNDER DEVELOPEMENT...
107  *
108  * Revision 1.2 2003/10/15 21:14:23 e_gourgoulhon
109  * Added method poisson_angu().
110  *
111  * Revision 1.1 2003/09/25 08:06:56 e_gourgoulhon
112  * First versions (use Cmp as intermediate quantities).
113  *
114  *
115  * $Header: /cvsroot/Lorene/C++/Source/Tensor/Scalar/scalar_pde.C,v 1.20 2014/10/13 08:53:46 j_novak Exp $
116  *
117  */
118 
119 // Header Lorene:
120 #include "map.h"
121 #include "scalar.h"
122 #include "tensor.h"
123 #include "param.h"
124 #include "cmp.h"
125 #include "param_elliptic.h"
126 
127 
128  //-----------------------------------//
129  // Scalar Poisson equation //
130  //-----------------------------------//
131 
132 // Version without parameters
133 // --------------------------
134 
135 namespace Lorene {
137 
138  Param bidon ;
139  Cmp csource(*this) ;
140  Cmp cresu(mp) ;
141  cresu = 0. ;
142 
143  mp->poisson(csource, bidon, cresu) ;
144 
145  Scalar resu(cresu) ;
146  return resu ;
147 }
148 
149 // Version with parameters
150 // -----------------------
151 
152 void Scalar::poisson(Param& par, Scalar& uu) const {
153 
154  Cmp csource(*this) ;
155  Cmp cuu(uu) ;
156 
157  mp->poisson(csource, par, cuu) ;
158 
159  uu = cuu ;
160 }
161 
162  //-----------------------------------------------//
163  // Scalar Poisson equation (TAU method) //
164  //----------------------------------------------//
165 
166  // without parameters
167  // --------------------------
168 
170 
171  Param bidon ;
172  Cmp csource(*this) ;
173  Cmp cresu(mp) ;
174  cresu = 0. ;
175 
176  mp->poisson_tau(csource, bidon, cresu) ;
177 
178  Scalar resu(cresu) ;
179  return resu ;
180 }
181 
182  // Version with parameters
183  // -----------------------
184 void Scalar::poisson_tau (Param& par, Scalar& uu) const {
185 
186  Cmp csource(*this) ;
187  Cmp cuu(uu) ;
188 
189  mp->poisson_tau(csource, par, cuu) ;
190 
191  uu = cuu ;
192 }
193 
194 
195  //-----------------------------------//
196  // Angular Poisson equation //
197  //-----------------------------------//
198 
199 
200 Scalar Scalar::poisson_angu(double lambda) const {
201 
202  Param bidon ;
203 
204  Scalar resu(*mp) ;
205  resu = 0. ;
206 
207  mp->poisson_angu(*this, bidon, resu, lambda) ;
208 
209  return resu ;
210 }
211 
212 
213  //-----------------------------------//
214  // Scalar d'Alembert equation //
215  //-----------------------------------//
216 
218  const Scalar& source) const {
219 
220  Scalar fjp1(*mp) ;
221 
222  mp->dalembert(par, fjp1, *this, fjm1, source) ;
223 
224  return fjp1 ;
225 
226 }
227 
228 
229  //-----------------------------------//
230  // General elliptic equation //
231  //-----------------------------------//
232 
233 
235 
236  // Right now, only applicable with affine mapping or log one
237  const Map_af* map_affine = dynamic_cast <const Map_af*> (mp) ;
238  const Map_log* map_log = dynamic_cast <const Map_log*> (mp) ;
239 
240  if ((map_affine == 0x0) && (map_log == 0x0)) {
241  cout << "sol_elliptic only defined for affine or log mapping" << endl ;
242  abort() ;
243  }
244 
245  Scalar res (*mp) ;
246  res.set_etat_qcq() ;
247 
248  if (map_affine != 0x0)
249  map_affine->sol_elliptic (ope_var, *this, res) ;
250  else
251  map_log->sol_elliptic (ope_var, *this, res) ;
252 
253  return (res) ;
254 }
255 
257 double fact_dir, double fact_neu) const {
258 
259  // Right now, only applicable with affine mapping or log one
260  const Map_af* map_affine = dynamic_cast <const Map_af*> (mp) ;
261  const Map_log* map_log = dynamic_cast <const Map_log*> (mp) ;
262 
263  if ((map_affine == 0x0) && (map_log == 0x0)) {
264  cout << "sol_elliptic only defined for affine or log mapping" << endl ;
265  abort() ;
266  }
267 
268  Scalar res (*mp) ;
269  res.set_etat_qcq() ;
270 
271  if (map_affine != 0x0)
272  map_affine->sol_elliptic_boundary (ope_var, *this, res, bound,
273 fact_dir, fact_neu ) ;
274  else
275  map_log->sol_elliptic_boundary (ope_var, *this, res, bound,
276 fact_dir, fact_neu ) ;
277 
278  return (res) ;
279 }
280 
281 
283 double fact_dir, double fact_neu) const {
284 
285  // Right now, only applicable with affine mapping or log one
286  const Map_af* map_affine = dynamic_cast <const Map_af*> (mp) ;
287  const Map_log* map_log = dynamic_cast <const Map_log*> (mp) ;
288 
289  if ((map_affine == 0x0) && (map_log == 0x0)) {
290  cout << "sol_elliptic only defined for affine or log mapping" << endl ;
291  abort() ;
292  }
293 
294  Scalar res (*mp) ;
295  res.set_etat_qcq() ;
296 
297  if (map_affine != 0x0)
298  map_affine->sol_elliptic_boundary (ope_var, *this, res, bound,
299 fact_dir, fact_neu ) ;
300  else
301  map_log->sol_elliptic_boundary (ope_var, *this, res, bound,
302 fact_dir, fact_neu ) ;
303 
304  return (res) ;
305 }
306 
307 
308 
309  //-----------------------------------//
310  // General elliptic equation //
311  // with no ZEC //
312  //-----------------------------------//
313 
314 Scalar Scalar::sol_elliptic_no_zec(Param_elliptic& ope_var, double val) const {
315 
316  // Right now, only applicable with affine mapping
317  const Map_af* map_affine = dynamic_cast <const Map_af*> (mp) ;
318  const Map_log* map_log = dynamic_cast <const Map_log*> (mp) ;
319 
320  if ((map_affine == 0x0) && (map_log == 0x0)) {
321  cout << "sol_elliptic_no_zec only defined for affine or log mapping" << endl ;
322  abort() ;
323  }
324 
325  Scalar res (*mp) ;
326  res.set_etat_qcq() ;
327 
328  if (map_affine != 0x0)
329  map_affine->sol_elliptic_no_zec (ope_var, *this, res, val) ;
330  else
331  map_log->sol_elliptic_no_zec (ope_var, *this, res, val) ;
332 
333  return (res) ;
334 }
335 
336  //-----------------------------------//
337  // General elliptic equation //
338  // with no ZEC //
339  //-----------------------------------//
340 
342 
343  // Right now, only applicable with affine mapping
344  const Map_af* map_affine = dynamic_cast <const Map_af*> (mp) ;
345 
346  if (map_affine == 0x0) {
347  cout << "sol_elliptic_no_zec only defined for affine or log mapping" << endl ;
348  abort() ;
349  }
350 
351  Scalar res (*mp) ;
352  res.set_etat_qcq() ;
353 
354  map_affine->sol_elliptic_only_zec (ope_var, *this, res, val) ;
355  return (res) ;
356 }
357  //-----------------------------------//
358  // General elliptic equation //
359  // with no ZEC and a //
360  // matching with sin(r)/r //
361  //-----------------------------------//
362 
363 Scalar Scalar::sol_elliptic_sin_zec(Param_elliptic& ope_var, double* amplis, double* phases)
364  const {
365 
366  // Right now, only applicable with affine mapping
367  const Map_af* map_affine = dynamic_cast <const Map_af*> (mp) ;
368 
369  if (map_affine == 0x0) {
370  cout << "sol_elliptic_sin_zec only defined for affine mapping" << endl ;
371  abort() ;
372  }
373 
374  Scalar res (*mp) ;
375  res.set_etat_qcq() ;
376 
377  map_affine->sol_elliptic_sin_zec (ope_var, *this, res, amplis, phases) ;
378 
379  return (res) ;
380 }
381  //-----------------------------------//
382  // General elliptic equation //
383  // fixing the radial derivative //
384  //-----------------------------------//
385 
387  Param_elliptic& ope_var) const {
388 
389  // Right now, only applicable with affine mapping
390  const Map_af* map_affine = dynamic_cast <const Map_af*> (mp) ;
391 
392  if (map_affine == 0x0) {
393  cout << "sol_elliptic_no_zec only defined for affine mapping" << endl ;
394  abort() ;
395  }
396 
397  Scalar res (*mp) ;
398  res.set_etat_qcq() ;
399 
400  map_affine->sol_elliptic_fixe_der_zero (valeur, ope_var, *this, res) ;
401 
402  return (res) ;
403 }
404 
405  //-----------------------------------//
406  // Two-dimensional Poisson eq. //
407  //-----------------------------------//
408 
410 
411  // Right now, only applicable with affine mapping
412  const Map_af* map_affine = dynamic_cast <const Map_af*> (mp) ;
413 
414  if (map_affine == 0x0) {
415  cout << "Poisson 2D only defined for affine mapping" << endl ;
416  abort() ;
417  }
418 
419  Scalar res (*mp) ;
420  res.set_etat_qcq() ;
421 
422  map_affine->sol_elliptic_2d(ope_var, *this, res) ;
423 
424  return (res) ;
425 }
426  //-----------------------------------//
427  // Pseudo-1dimensional eq. //
428  //-----------------------------------//
429 
431 
432  // Right now, only applicable with affine mapping
433  const Map_af* map_affine = dynamic_cast <const Map_af*> (mp) ;
434 
435  if (map_affine == 0x0) {
436  cout << "Pseudo_1d only defined for affine mapping" << endl ;
437  abort() ;
438  }
439 
440  Scalar res (*mp) ;
441  res.set_etat_qcq() ;
442 
443  map_affine->sol_elliptic_pseudo_1d(ope_var, *this, res) ;
444 
445  return (res) ;
446 }
447 }
Lorene::Map::poisson_angu
virtual void poisson_angu(const Scalar &source, Param &par, Scalar &uu, double lambda=0) const =0
Computes the solution of the generalized angular Poisson equation.
Lorene::Scalar::sol_elliptic
Scalar sol_elliptic(Param_elliptic &params) const
Resolution of a general elliptic equation, putting zero at infinity.
Definition: scalar_pde.C:234
Lorene::Map_af::sol_elliptic_sin_zec
void sol_elliptic_sin_zec(Param_elliptic &params, const Scalar &so, Scalar &uu, double *coefs, double *) const
General elliptic solver.
Definition: map_af_elliptic.C:426
Lorene
Lorene prototypes.
Definition: app_hor.h:64
Lorene::Scalar::sol_elliptic_sin_zec
Scalar sol_elliptic_sin_zec(Param_elliptic &params, double *coefs, double *phases) const
General elliptic solver.
Definition: scalar_pde.C:363
Lorene::Map_log::sol_elliptic_boundary
void sol_elliptic_boundary(Param_elliptic &params, const Scalar &so, Scalar &uu, const Mtbl_cf &bound, double fact_dir, double fact_neu) const
General elliptic solver including inner boundary conditions.
Definition: map_log_elliptic.C:132
Lorene::Map_log::sol_elliptic_no_zec
void sol_elliptic_no_zec(Param_elliptic &params, const Scalar &so, Scalar &uu, double) const
General elliptic solver.
Definition: map_log_elliptic.C:275
Lorene::Scalar
Tensor field of valence 0 (or component of a tensorial field).
Definition: scalar.h:387
Lorene::Param_elliptic
This class contains the parameters needed to call the general elliptic solver.
Definition: param_elliptic.h:135
Lorene::Scalar::sol_elliptic_2d
Scalar sol_elliptic_2d(Param_elliptic &) const
Solves the scalar 2-dimensional elliptic equation with *this as a source.
Definition: scalar_pde.C:409
Lorene::Map_log::sol_elliptic
void sol_elliptic(Param_elliptic &params, const Scalar &so, Scalar &uu) const
General elliptic solver.
Definition: map_log_elliptic.C:75
Lorene::Scalar::sol_elliptic_pseudo_1d
Scalar sol_elliptic_pseudo_1d(Param_elliptic &) const
Solves a pseudo-1d elliptic equation with *this as a source.
Definition: scalar_pde.C:430
Lorene::Map_af::sol_elliptic
void sol_elliptic(Param_elliptic &params, const Scalar &so, Scalar &uu) const
General elliptic solver.
Definition: map_af_elliptic.C:100
Lorene::Map::poisson_tau
virtual void poisson_tau(const Cmp &source, Param &par, Cmp &uu) const =0
Computes the solution of a scalar Poisson equationwith a Tau method.
Lorene::Map_af::sol_elliptic_only_zec
void sol_elliptic_only_zec(Param_elliptic &params, const Scalar &so, Scalar &uu, double val) const
General elliptic solver.
Definition: map_af_elliptic.C:368
Lorene::Map_af::sol_elliptic_2d
void sol_elliptic_2d(Param_elliptic &, const Scalar &, Scalar &) const
General elliptic solver in a 2D case.
Definition: map_af_elliptic_2d.C:58
Lorene::Tensor::mp
const Map *const mp
Mapping on which the numerical values at the grid points are defined.
Definition: tensor.h:295
Lorene::Map::dalembert
virtual void dalembert(Param &par, Scalar &fJp1, const Scalar &fJ, const Scalar &fJm1, const Scalar &source) const =0
Performs one time-step integration of the d'Alembert scalar equation.
Lorene::Cmp
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition: cmp.h:446
Lorene::Scalar::sol_elliptic_fixe_der_zero
Scalar sol_elliptic_fixe_der_zero(double val, Param_elliptic &params) const
Resolution of a general elliptic equation fixing the dericative at the origin and relaxing one contin...
Definition: scalar_pde.C:386
Lorene::Map_af
Affine radial mapping.
Definition: map.h:2027
Lorene::Scalar::poisson_tau
Scalar poisson_tau() const
Solves the scalar Poisson equation with *this as a source using a real Tau method The source of the ...
Definition: scalar_pde.C:169
Lorene::Scalar::set_etat_qcq
virtual void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: scalar.C:353
Lorene::Map_af::sol_elliptic_no_zec
void sol_elliptic_no_zec(Param_elliptic &params, const Scalar &so, Scalar &uu, double val) const
General elliptic solver.
Definition: map_af_elliptic.C:312
Lorene::Map::poisson
virtual void poisson(const Cmp &source, Param &par, Cmp &uu) const =0
Computes the solution of a scalar Poisson equation.
Lorene::Map_af::sol_elliptic_pseudo_1d
void sol_elliptic_pseudo_1d(Param_elliptic &, const Scalar &, Scalar &) const
General elliptic solver in a pseudo 1d case.
Definition: map_af_elliptic_pseudo_1d.C:59
Lorene::Scalar::avance_dalembert
Scalar avance_dalembert(Param &par, const Scalar &fJm1, const Scalar &source) const
Performs one time-step integration (from ) of the scalar d'Alembert equation with *this being the val...
Definition: scalar_pde.C:217
Lorene::Param
Parameter storage.
Definition: param.h:125
Lorene::Scalar::sol_elliptic_no_zec
Scalar sol_elliptic_no_zec(Param_elliptic &params, double val=0) const
Resolution of a general elliptic equation, putting a given value at the outermost shell and not solvi...
Definition: scalar_pde.C:314
Lorene::Scalar::poisson
Scalar poisson() const
Solves the scalar Poisson equation with *this as a source.
Definition: scalar_pde.C:136
Lorene::Scalar::poisson_angu
Scalar poisson_angu(double lambda=0) const
Solves the (generalized) angular Poisson equation with *this as source.
Definition: scalar_pde.C:200
Lorene::Scalar::sol_elliptic_only_zec
Scalar sol_elliptic_only_zec(Param_elliptic &params, double val) const
Resolution of a general elliptic equation solving in the compactified domain and putting a given valu...
Definition: scalar_pde.C:341
Lorene::Map_log
Logarithmic radial mapping.
Definition: map.h:3583
Lorene::Mtbl_cf
Coefficients storage for the multi-domain spectral method.
Definition: mtbl_cf.h:186
Lorene::Map_af::sol_elliptic_fixe_der_zero
void sol_elliptic_fixe_der_zero(double val, Param_elliptic &params, const Scalar &so, Scalar &uu) const
General elliptic solver fixing the derivative at the origin and relaxing the continuity of the first ...
Definition: map_af_elliptic.C:483
Lorene::Scalar::sol_elliptic_boundary
Scalar sol_elliptic_boundary(Param_elliptic &params, const Mtbl_cf &bound, double fact_dir, double fact_neu) const
Resolution of a general elliptic equation, putting zero at infinity and with inner boundary condition...
Definition: scalar_pde.C:256
Lorene::Map_af::sol_elliptic_boundary
void sol_elliptic_boundary(Param_elliptic &params, const Scalar &so, Scalar &uu, const Mtbl_cf &bound, double fact_dir, double fact_neu) const
General elliptic solver including inner boundary conditions.
Definition: map_af_elliptic.C:159