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