LORENE
param_elliptic_val_lim.C
1 /*
2  * Copyright (c) 2004 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 version 2
8  * as published by the Free Software Foundation.
9  *
10  * LORENE is distributed in the hope that it will be useful,
11  * but WITHOUT ANY WARRANTY; without even the implied warranty of
12  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13  * GNU General Public License for more details.
14  *
15  * You should have received a copy of the GNU General Public License
16  * along with LORENE; if not, write to the Free Software
17  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
18  *
19  */
20 
21 
22 
23 /*
24  * $Id: param_elliptic_val_lim.C,v 1.4 2016/12/05 16:18:14 j_novak Exp $
25  * $Log: param_elliptic_val_lim.C,v $
26  * Revision 1.4 2016/12/05 16:18:14 j_novak
27  * Suppression of some global variables (file names, loch, ...) to prevent redefinitions
28  *
29  * Revision 1.3 2014/10/13 08:53:37 j_novak
30  * Lorene classes and functions now belong to the namespace Lorene.
31  *
32  * Revision 1.2 2014/10/06 15:13:16 j_novak
33  * Modified #include directives to use c++ syntax.
34  *
35  * Revision 1.1 2004/08/24 09:14:49 p_grandclement
36  * Addition of some new operators, like Poisson in 2d... It now requieres the
37  * GSL library to work.
38  *
39  * Also, the way a variable change is stored by a Param_elliptic is changed and
40  * no longer uses Change_var but rather 2 Scalars. The codes using that feature
41  * will requiere some modification. (It should concern only the ones about monopoles)
42  *
43 
44  * $Header: /cvsroot/Lorene/C++/Source/Param_elliptic/param_elliptic_val_lim.C,v 1.4 2016/12/05 16:18:14 j_novak Exp $
45  *
46  */
47 
48 #include "headcpp.h"
49 
50 #include <cmath>
51 #include <cstdlib>
52 
53 #include "base_val.h"
54 #include "param_elliptic.h"
55 #include "proto.h"
56 
57 namespace Lorene {
58 double Param_elliptic::F_plus (int zone, int k, int j) const {
59 
60  if (done_F(zone, k, j) == 0)
61  compute_val_F(zone, k, j) ;
62 
63  return val_F_plus (zone, k, j) ;
64 }
65 
66 double Param_elliptic::F_minus (int zone, int k, int j) const {
67 
68  if (done_F(zone, k, j) == 0)
69  compute_val_F(zone, k, j) ;
70 
71  return val_F_minus (zone, k, j) ;
72 }
73 
74 double Param_elliptic::dF_plus (int zone, int k, int j) const {
75 
76  if (done_F(zone, k, j) == 0)
77  compute_val_F(zone, k, j) ;
78 
79  return val_dF_plus (zone, k, j) ;
80 }
81 
82 double Param_elliptic::dF_minus (int zone, int k, int j) const {
83 
84  if (done_F(zone, k, j) == 0)
85  compute_val_F(zone, k, j) ;
86 
87  return val_dF_minus (zone, k, j) ;
88 }
89 
90 
91 double Param_elliptic::G_plus (int zone) const {
92 
93  if (done_G(zone) == 0)
94  compute_val_G(zone) ;
95 
96  return val_G_plus (zone) ;
97 }
98 
99 double Param_elliptic::G_minus (int zone) const {
100 
101  if (done_G(zone) == 0)
102  compute_val_G(zone) ;
103 
104  return val_G_minus (zone) ;
105 }
106 
107 double Param_elliptic::dG_plus (int zone) const {
108 
109  if (done_G(zone) == 0)
110  compute_val_G(zone) ;
111 
112  return val_dG_plus (zone) ;
113 }
114 
115 double Param_elliptic::dG_minus (int zone) const {
116 
117  if (done_G(zone) == 0)
118  compute_val_G(zone) ;
119 
120  return val_dG_minus (zone) ;
121 }
122 
123 
124 void Param_elliptic::compute_val_F (int zone, int k, int j) const {
125 
126  int nr = get_mp().get_mg()->get_nr(zone) ;
127  Tbl coefs (nr) ;
128  coefs.set_etat_qcq() ;
129 
130  bool zero ;
131 
132  if (var_F.get_spectral_va().c_cf->get_etat() == ETATZERO)
133  zero = true ;
134  else
135  if ((*var_F.get_spectral_va().c_cf)(zone).get_etat() == ETATZERO)
136  zero = true ;
137  else
138  zero = false ;
139 
140  if (zero)
141  coefs.annule_hard() ;
142  else
143  for (int i=0 ; i<nr ; i++)
144  coefs.set(i) = (*var_F.get_spectral_va().c_cf)(zone, k, j, i) ;
145 
146  int lq, mq ;
147  int base_r ;
148  var_F.get_spectral_va().base.give_quant_numbers (zone, k,j, lq, mq, base_r) ;
149 
150  double alpha = get_alpha(zone) ;
151 
152  Tbl output (val_solp(coefs, alpha, base_r)) ;
153 
154  // On range :
155  val_F_plus.set(zone, k, j) = output(0) ;
156  val_F_minus.set(zone, k, j) = output(1) ;
157  val_dF_plus.set(zone, k, j) = output(2) ;
158  val_dF_minus.set(zone, k, j) = output(3) ;
159  done_F.set(zone, k, j) = 1 ;
160 
161 }
162 
163 void Param_elliptic::compute_val_G (int zone) const {
164 
165  int nr = get_mp().get_mg()->get_nr(zone) ;
166  Tbl coefs (nr) ;
167  coefs.set_etat_qcq() ;
168 
169 
170  bool zero ;
171  if (var_G.get_spectral_va().c_cf->get_etat() == ETATZERO)
172  zero = true ;
173  else
174  if ((*var_G.get_spectral_va().c_cf)(zone).get_etat() == ETATZERO)
175  zero = true ;
176  else
177  zero = false ;
178  if (zero)
179  coefs.annule_hard() ;
180  else
181  for (int i=0 ; i<nr ; i++)
182  coefs.set(i) = (*var_G.get_spectral_va().c_cf)(zone, 0, 0, i) ;
183 
184  int lq, mq ;
185  int base_r ;
186  var_G.get_spectral_va().base.give_quant_numbers (zone, 0, 0, lq, mq, base_r) ;
187 
188 
189  double alpha = get_alpha(zone) ;
190 
191  Tbl output (val_solp(coefs, alpha, base_r)) ;
192 
193  // On range :
194  val_G_plus.set(zone) = output(0) ;
195  val_G_minus.set(zone) = output(1) ;
196  val_dG_plus.set(zone) = output(2) ;
197  val_dG_minus.set(zone) = output(3) ;
198  done_G.set(zone) = 1 ;
199 }
200 }
Mtbl_cf * c_cf
Coefficients of the spectral expansion of the function.
Definition: valeur.h:312
int & set(int i)
Read/write of a particular element (index i ) (1D case)
Definition: itbl.h:247
void give_quant_numbers(int, int, int, int &, int &, int &) const
Computes the various quantum numbers and 1d radial base.
double dF_plus(int zone, int k, int j) const
Returns the value of the radial derivative of F, for a given angular point, at the outer boundary of ...
double dG_plus(int zone) const
Returns the value of the radial derivative of G, for a given angular point, at the outer boundary of ...
Lorene prototypes.
Definition: app_hor.h:67
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:777
double & set(int i)
Read/write of a particular element (index i) (1D case)
Definition: tbl.h:301
Itbl done_G
Stores what has been computed for G.
void compute_val_F(int, int, int) const
Computes the various values of F.
int get_etat() const
Returns the logical state.
Definition: mtbl_cf.h:466
Tbl val_F_minus
Values of F at the inner boundaries of the various domains.
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: tbl.C:364
double G_plus(int zone) const
Returns the value of G, for a given angular point, at the outer boundary of the domain zone ;...
Tbl val_G_plus
Values of G at the outer boundaries of the various domains.
Scalar var_F
Additive variable change function.
double dF_minus(int zone, int k, int j) const
Returns the value of the radial derivative of F, for a given angular point, at the inner boundary of ...
Base_val base
Bases on which the spectral expansion is performed.
Definition: valeur.h:315
Tbl val_G_minus
Values of G at the inner boundaries of the various domains.
Tbl val_F_plus
Values of F at the outer boundaries of the various domains.
double F_minus(int zone, int k, int j) const
Returns the value of F, for a given angular point, at the inner boundary of the domain zone ;...
const Map_radial & get_mp() const
Returns the mapping.
Scalar var_G
Multiplicative variable change that must be sphericaly symetric !
double G_minus(int zone) const
Returns the value of G, for a given angular point, at the inner boundary of the domain zone ;...
Tbl val_dG_minus
Values of the derivative of G at the inner boundaries of the various domains.
Tbl val_dG_plus
Values of the derivative of G at the outer boundaries of the various domains.
int get_nr(int l) const
Returns the number of points in the radial direction ( ) in domain no. l.
Definition: grilles.h:469
double dG_minus(int zone) const
Returns the value of the radial derivative of G, for a given angular point, at the inner boundary of ...
Tbl val_dF_minus
Values of the derivative of F at the inner boundaries of the various domains.
Itbl done_F
Stores what has been computed for F.
Basic array class.
Definition: tbl.h:164
void compute_val_G(int) const
Computes the various values of G.
void annule_hard()
Sets the Tbl to zero in a hard way.
Definition: tbl.C:375
Tbl val_dF_plus
Values of the derivative of F at the outer boundaries of the various domains.
double F_plus(int zone, int k, int j) const
Returns the value of F, for a given angular point, at the outer boundary of the domain zone ;...
const Valeur & get_spectral_va() const
Returns va (read only version)
Definition: scalar.h:607