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