LORENE
map_log_deriv.C
1 /*
2  * Computations of Scalar partial derivatives for a Map_log mapping
3  */
4 
5 /*
6  * Copyright (c) 2004 Philippe Grandclement
7  *
8  * This file is part of LORENE.
9  *
10  * LORENE is free software; you can redistribute it and/or modify
11  * it under the terms of the GNU General Public License as published by
12  * the Free Software Foundation; either version 2 of the License, or
13  * (at your option) any later version.
14  *
15  * LORENE is distributed in the hope that it will be useful,
16  * but WITHOUT ANY WARRANTY; without even the implied warranty of
17  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18  * GNU General Public License for more details.
19  *
20  * You should have received a copy of the GNU General Public License
21  * along with LORENE; if not, write to the Free Software
22  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
23  *
24  */
25 
26 
27 
28 
29 /*
30  * $Id: map_log_deriv.C,v 1.4 2016/12/05 16:17:58 j_novak Exp $
31  * $Log: map_log_deriv.C,v $
32  * Revision 1.4 2016/12/05 16:17:58 j_novak
33  * Suppression of some global variables (file names, loch, ...) to prevent redefinitions
34  *
35  * Revision 1.3 2014/10/13 08:53:05 j_novak
36  * Lorene classes and functions now belong to the namespace Lorene.
37  *
38  * Revision 1.2 2012/01/17 10:33:43 j_penner
39  * added a derivative with respect to the computational coordinate xi
40  *
41  * Revision 1.1 2004/06/22 08:49:58 p_grandclement
42  * Addition of everything needed for using the logarithmic mapping
43  *
44  *
45  * $Header: /cvsroot/Lorene/C++/Source/Map/map_log_deriv.C,v 1.4 2016/12/05 16:17:58 j_novak Exp $
46  *
47  */
48 
49 // Header Lorene
50 #include "map.h"
51 #include "tensor.h"
52 
53  //---------------------------------------------------
54  // d/d\xi
55  //---------------------------------------------------
56 namespace Lorene {
57 void Map_log::dsdxi(const Scalar& uu, Scalar& resu) const {
58 
59  assert (uu.get_etat() != ETATNONDEF) ;
60  assert (uu.get_mp().get_mg() == mg) ;
61 
62 
63  if (uu.get_etat() == ETATZERO) {
64  resu.set_etat_zero() ;
65  }
66  else {
67  assert( uu.get_etat() == ETATQCQ ) ;
68 
69  const Valeur& uuva = uu.get_spectral_va() ;
70 
71  uuva.coef() ; // (uu.va).c_cf is up to date
72 
73  int nz = mg->get_nzone() ;
74  int nzm1 = nz - 1 ;
75 
76  if ( uu.get_dzpuis() == 0 ) {
77  resu = uuva.dsdx() ; // dsdx == d/d\xi
78 
79  if (mg->get_type_r(nzm1) == UNSURR) {
80  resu.set_dzpuis(2) ; // r^2 d/dr has been computed in the
81  // external domain
82  }
83  }
84  else {
85  assert(mg->get_type_r(nzm1) == UNSURR) ;
86 
87  int dzp = uu.get_dzpuis() ;
88 
89  resu = uuva.dsdx() ;
90  resu.annule_domain(nzm1) ; // zero in the CED
91 
92  // Special treatment in the CED
93  Valeur tmp_ced = - uuva.dsdx() ;
94  tmp_ced.annule(0, nz-2) ; // only non zero in the CED
95  tmp_ced.mult_xm1_zec() ;
96  tmp_ced.set(nzm1) -= dzp * uuva(nzm1) ;
97 
98  // Recombination shells + CED :
99  resu.set_spectral_va() += tmp_ced ;
100 
101  resu.set_dzpuis(dzp+1) ;
102 
103  }
104  resu.set_spectral_base( uuva.dsdx().get_base() ) ; // same basis as d/dxi
105  }
106 }
107 
108  //---------------------------------------------------
109  // Derivee standard par rapport au "vrai" rayon .....
110  //---------------------------------------------------
111 void Map_log::dsdr(const Scalar& uu, Scalar& resu) const {
112 
113  assert (uu.get_etat() != ETATNONDEF) ;
114  assert (uu.get_mp().get_mg() == mg) ;
115 
116 
117  if (uu.get_etat() == ETATZERO) {
118  resu.set_etat_zero() ;
119  }
120  else {
121  assert( uu.get_etat() == ETATQCQ ) ;
122 
123  const Valeur& uuva = uu.get_spectral_va() ;
124 
125  uuva.coef() ; // (uu.va).c_cf is up to date
126 
127  int nz = mg->get_nzone() ;
128  int nzm1 = nz - 1 ;
129 
130  if ( uu.get_dzpuis() == 0 ) {
131  resu = uuva.dsdx() * dxdr ; // dxdr = dxi/dR, - dxi/dU (ZEC)
132 
133  if (mg->get_type_r(nzm1) == UNSURR) {
134  resu.set_dzpuis(2) ; // r^2 d/dr has been computed in the
135  // external domain
136  }
137  }
138  else {
139  assert(mg->get_type_r(nzm1) == UNSURR) ;
140 
141  int dzp = uu.get_dzpuis() ;
142 
143  resu = uuva.dsdx() * dxdr ;
144  resu.annule_domain(nzm1) ; // zero in the CED
145 
146  // Special treatment in the CED
147  Valeur tmp_ced = - uuva.dsdx() ;
148  tmp_ced.annule(0, nz-2) ; // only non zero in the CED
149  tmp_ced.mult_xm1_zec() ;
150  tmp_ced.set(nzm1) -= dzp * uuva(nzm1) ;
151 
152  // Recombination shells + CED :
153  resu.set_spectral_va() += tmp_ced ;
154 
155  resu.set_dzpuis(dzp+1) ;
156 
157  }
158  resu.set_spectral_base( uuva.dsdx().get_base() ) ; // same basis as d/dxi
159  }
160 }
161 
162  //---------------------------------------------------
163  // Derivee par rapport au rayon numerique (r ou lnr)
164  //---------------------------------------------------
165 void Map_log::dsdradial (const Scalar& uu, Scalar& resu) const {
166 
167  assert (uu.get_etat() != ETATNONDEF) ;
168  assert (uu.get_mp().get_mg() == mg) ;
169 
170 
171  if (uu.get_etat() == ETATZERO) {
172  resu.set_etat_zero() ;
173  }
174  else {
175  assert( uu.get_etat() == ETATQCQ ) ;
176 
177  const Valeur& uuva = uu.get_spectral_va() ;
178 
179  uuva.coef() ; // (uu.va).c_cf is up to date
180 
181  int nz = mg->get_nzone() ;
182  int nzm1 = nz - 1 ;
183 
184  if ( uu.get_dzpuis() == 0 ) {
185  resu = uuva.dsdx() * dxdlnr ; // dxdr = dxi/dR, - dxi/dU (ZEC)
186 
187  if (mg->get_type_r(nzm1) == UNSURR) {
188  resu.set_dzpuis(2) ; // r^2 d/dr has been computed in the
189  // external domain
190  }
191  }
192  else {
193  assert(mg->get_type_r(nzm1) == UNSURR) ;
194 
195  int dzp = uu.get_dzpuis() ;
196 
197  resu = uuva.dsdx() * dxdlnr ;
198  resu.annule_domain(nzm1) ; // zero in the CED
199 
200  // Special treatment in the CED
201  Valeur tmp_ced = - uuva.dsdx() ;
202  tmp_ced.annule(0, nz-2) ; // only non zero in the CED
203  tmp_ced.mult_xm1_zec() ;
204  tmp_ced.set(nzm1) -= dzp * uuva(nzm1) ;
205 
206  // Recombination shells + CED :
207  resu.set_spectral_va() += tmp_ced ;
208 
209  resu.set_dzpuis(dzp+1) ;
210 
211  }
212  resu.set_spectral_base( uuva.dsdx().get_base() ) ; // same basis as d/dxi
213  }
214 }
215 }
void annule_domain(int l)
Sets the Tensor to zero in a given domain.
Definition: tensor.C:675
const Valeur & dsdx() const
Returns of *this.
Definition: valeur_dsdx.C:114
void coef() const
Computes the coeffcients of *this.
Definition: valeur_coef.C:151
virtual void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition: scalar.C:330
Lorene prototypes.
Definition: app_hor.h:67
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:783
Tensor field of valence 0 (or component of a tensorial field).
Definition: scalar.h:393
void mult_xm1_zec()
Applies the following operator to *this : \ Id (r sampling = RARE, FIN ) \ (r -sampling = UNSURR ) ...
const Scalar & dsdx() const
Returns of *this , where .
Definition: scalar_deriv.C:266
Values and coefficients of a (real-value) function.
Definition: valeur.h:297
int get_etat() const
Returns the logical state ETATNONDEF (undefined), ETATZERO (null) or ETATQCQ (ordinary).
Definition: scalar.h:560
void annule(int l)
Sets the Valeur to zero in a given domain.
Definition: valeur.C:747
void set_dzpuis(int)
Modifies the dzpuis flag.
Definition: scalar.C:814
Coord dxdlnr
Same as dxdr if the domains where the description is affine and where it is logarithmic.
Definition: map.h:3636
Coord dxdr
in the nucleus and in the non-compactified shells; \ in the compactified outer domain.
Definition: map.h:1581
int get_dzpuis() const
Returns dzpuis.
Definition: scalar.h:563
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:465
virtual void dsdr(const Scalar &ci, Scalar &resu) const
Computes of a Scalar.
void set_spectral_base(const Base_val &)
Sets the spectral bases of the Valeur va
Definition: scalar.C:803
const Mg3d * mg
Pointer on the multi-grid Mgd3 on which this is defined.
Definition: map.h:694
virtual void dsdradial(const Scalar &uu, Scalar &resu) const
Computes of a Scalar if the description is affine and if it is logarithmic.
virtual void dsdxi(const Scalar &ci, Scalar &resu) const
Computes of a Scalar.
Definition: map_log_deriv.C:57
Valeur & set_spectral_va()
Returns va (read/write version)
Definition: scalar.h:610
int get_type_r(int l) const
Returns the type of sampling in the radial direction in domain no.
Definition: grilles.h:491
const Map & get_mp() const
Returns the mapping.
Definition: tensor.h:874
const Valeur & get_spectral_va() const
Returns va (read only version)
Definition: scalar.h:607
Tbl & set(int l)
Read/write of the value in a given domain (configuration space).
Definition: valeur.h:373