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
27char map_log_deriv_C[] = "$Header: /cvsroot/Lorene/C++/Source/Map/map_log_deriv.C,v 1.3 2014/10/13 08:53:05 j_novak Exp $" ;
28
29/*
30 * $Id: map_log_deriv.C,v 1.3 2014/10/13 08:53:05 j_novak Exp $
31 * $Log: map_log_deriv.C,v $
32 * Revision 1.3 2014/10/13 08:53:05 j_novak
33 * Lorene classes and functions now belong to the namespace Lorene.
34 *
35 * Revision 1.2 2012/01/17 10:33:43 j_penner
36 * added a derivative with respect to the computational coordinate xi
37 *
38 * Revision 1.1 2004/06/22 08:49:58 p_grandclement
39 * Addition of everything needed for using the logarithmic mapping
40 *
41 *
42 * $Header: /cvsroot/Lorene/C++/Source/Map/map_log_deriv.C,v 1.3 2014/10/13 08:53:05 j_novak Exp $
43 *
44 */
45
46// Header Lorene
47#include "map.h"
48#include "tensor.h"
49
50 //---------------------------------------------------
51 // d/d\xi
52 //---------------------------------------------------
53namespace Lorene {
54void Map_log::dsdxi(const Scalar& uu, Scalar& resu) const {
55
56 assert (uu.get_etat() != ETATNONDEF) ;
57 assert (uu.get_mp().get_mg() == mg) ;
58
59
60 if (uu.get_etat() == ETATZERO) {
61 resu.set_etat_zero() ;
62 }
63 else {
64 assert( uu.get_etat() == ETATQCQ ) ;
65
66 const Valeur& uuva = uu.get_spectral_va() ;
67
68 uuva.coef() ; // (uu.va).c_cf is up to date
69
70 int nz = mg->get_nzone() ;
71 int nzm1 = nz - 1 ;
72
73 if ( uu.get_dzpuis() == 0 ) {
74 resu = uuva.dsdx() ; // dsdx == d/d\xi
75
76 if (mg->get_type_r(nzm1) == UNSURR) {
77 resu.set_dzpuis(2) ; // r^2 d/dr has been computed in the
78 // external domain
79 }
80 }
81 else {
82 assert(mg->get_type_r(nzm1) == UNSURR) ;
83
84 int dzp = uu.get_dzpuis() ;
85
86 resu = uuva.dsdx() ;
87 resu.annule_domain(nzm1) ; // zero in the CED
88
89 // Special treatment in the CED
90 Valeur tmp_ced = - uuva.dsdx() ;
91 tmp_ced.annule(0, nz-2) ; // only non zero in the CED
92 tmp_ced.mult_xm1_zec() ;
93 tmp_ced.set(nzm1) -= dzp * uuva(nzm1) ;
94
95 // Recombination shells + CED :
96 resu.set_spectral_va() += tmp_ced ;
97
98 resu.set_dzpuis(dzp+1) ;
99
100 }
101 resu.set_spectral_base( uuva.dsdx().get_base() ) ; // same basis as d/dxi
102 }
103}
104
105 //---------------------------------------------------
106 // Derivee standard par rapport au "vrai" rayon .....
107 //---------------------------------------------------
108void Map_log::dsdr(const Scalar& uu, Scalar& resu) const {
109
110 assert (uu.get_etat() != ETATNONDEF) ;
111 assert (uu.get_mp().get_mg() == mg) ;
112
113
114 if (uu.get_etat() == ETATZERO) {
115 resu.set_etat_zero() ;
116 }
117 else {
118 assert( uu.get_etat() == ETATQCQ ) ;
119
120 const Valeur& uuva = uu.get_spectral_va() ;
121
122 uuva.coef() ; // (uu.va).c_cf is up to date
123
124 int nz = mg->get_nzone() ;
125 int nzm1 = nz - 1 ;
126
127 if ( uu.get_dzpuis() == 0 ) {
128 resu = uuva.dsdx() * dxdr ; // dxdr = dxi/dR, - dxi/dU (ZEC)
129
130 if (mg->get_type_r(nzm1) == UNSURR) {
131 resu.set_dzpuis(2) ; // r^2 d/dr has been computed in the
132 // external domain
133 }
134 }
135 else {
136 assert(mg->get_type_r(nzm1) == UNSURR) ;
137
138 int dzp = uu.get_dzpuis() ;
139
140 resu = uuva.dsdx() * dxdr ;
141 resu.annule_domain(nzm1) ; // zero in the CED
142
143 // Special treatment in the CED
144 Valeur tmp_ced = - uuva.dsdx() ;
145 tmp_ced.annule(0, nz-2) ; // only non zero in the CED
146 tmp_ced.mult_xm1_zec() ;
147 tmp_ced.set(nzm1) -= dzp * uuva(nzm1) ;
148
149 // Recombination shells + CED :
150 resu.set_spectral_va() += tmp_ced ;
151
152 resu.set_dzpuis(dzp+1) ;
153
154 }
155 resu.set_spectral_base( uuva.dsdx().get_base() ) ; // same basis as d/dxi
156 }
157}
158
159 //---------------------------------------------------
160 // Derivee par rapport au rayon numerique (r ou lnr)
161 //---------------------------------------------------
162void Map_log::dsdradial (const Scalar& uu, Scalar& resu) const {
163
164 assert (uu.get_etat() != ETATNONDEF) ;
165 assert (uu.get_mp().get_mg() == mg) ;
166
167
168 if (uu.get_etat() == ETATZERO) {
169 resu.set_etat_zero() ;
170 }
171 else {
172 assert( uu.get_etat() == ETATQCQ ) ;
173
174 const Valeur& uuva = uu.get_spectral_va() ;
175
176 uuva.coef() ; // (uu.va).c_cf is up to date
177
178 int nz = mg->get_nzone() ;
179 int nzm1 = nz - 1 ;
180
181 if ( uu.get_dzpuis() == 0 ) {
182 resu = uuva.dsdx() * dxdlnr ; // dxdr = dxi/dR, - dxi/dU (ZEC)
183
184 if (mg->get_type_r(nzm1) == UNSURR) {
185 resu.set_dzpuis(2) ; // r^2 d/dr has been computed in the
186 // external domain
187 }
188 }
189 else {
190 assert(mg->get_type_r(nzm1) == UNSURR) ;
191
192 int dzp = uu.get_dzpuis() ;
193
194 resu = uuva.dsdx() * dxdlnr ;
195 resu.annule_domain(nzm1) ; // zero in the CED
196
197 // Special treatment in the CED
198 Valeur tmp_ced = - uuva.dsdx() ;
199 tmp_ced.annule(0, nz-2) ; // only non zero in the CED
200 tmp_ced.mult_xm1_zec() ;
201 tmp_ced.set(nzm1) -= dzp * uuva(nzm1) ;
202
203 // Recombination shells + CED :
204 resu.set_spectral_va() += tmp_ced ;
205
206 resu.set_dzpuis(dzp+1) ;
207
208 }
209 resu.set_spectral_base( uuva.dsdx().get_base() ) ; // same basis as d/dxi
210 }
211}
212}
virtual void dsdxi(const Scalar &ci, Scalar &resu) const
Computes of a Scalar.
virtual void dsdradial(const Scalar &uu, Scalar &resu) const
Computes of a Scalar if the description is affine and if it is logarithmic.
Coord dxdlnr
Same as dxdr if the domains where the description is affine and where it is logarithmic.
Definition map.h:3603
virtual void dsdr(const Scalar &ci, Scalar &resu) const
Computes of a Scalar.
Coord dxdr
in the nucleus and in the non-compactified shells; \ in the compactified outer domain.
Definition map.h:1560
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition map.h:765
const Mg3d * mg
Pointer on the multi-grid Mgd3 on which this is defined
Definition map.h:676
int get_nzone() const
Returns the number of domains.
Definition grilles.h:448
int get_type_r(int l) const
Returns the type of sampling in the radial direction in domain no.
Definition grilles.h:474
Tensor field of valence 0 (or component of a tensorial field).
Definition scalar.h:387
int get_dzpuis() const
Returns dzpuis.
Definition scalar.h:557
virtual void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition scalar.C:324
Valeur & set_spectral_va()
Returns va (read/write version)
Definition scalar.h:604
const Valeur & get_spectral_va() const
Returns va (read only version)
Definition scalar.h:601
int get_etat() const
Returns the logical state ETATNONDEF (undefined), ETATZERO (null) or ETATQCQ (ordinary).
Definition scalar.h:554
void set_dzpuis(int)
Modifies the dzpuis flag.
Definition scalar.C:808
void set_spectral_base(const Base_val &)
Sets the spectral bases of the Valeur va
Definition scalar.C:797
Values and coefficients of a (real-value) function.
Definition valeur.h:287
void mult_xm1_zec()
Applies the following operator to *this : \ Id (r sampling = RARE, FIN ) \ (r -sampling = UNSURR )
const Base_val & get_base() const
Return the bases for spectral expansions (member base )
Definition valeur.h:480
void annule(int l)
Sets the Valeur to zero in a given domain.
Definition valeur.C:744
Tbl & set(int l)
Read/write of the value in a given domain (configuration space).
Definition valeur.h:363
const Valeur & dsdx() const
Returns of *this.
void coef() const
Computes the coeffcients of *this.
const Map & get_mp() const
Returns the mapping.
Definition tensor.h:861
void annule_domain(int l)
Sets the Tensor to zero in a given domain.
Definition tensor.C:666
Lorene prototypes.
Definition app_hor.h:64