LORENE
scalar_raccord.C
1 /*
2  * Copyright (c) 2003 Eric Gourgoulhon & Jerome Novak
3  *
4  * Copyright (c) 2000-2001 Philippe Grandclement (for preceding Cmp version)
5  *
6  * This file is part of LORENE.
7  *
8  * LORENE is free software; you can redistribute it and/or modify
9  * it under the terms of the GNU General Public License as published by
10  * the Free Software Foundation; either version 2 of the License, or
11  * (at your option) any later version.
12  *
13  * LORENE is distributed in the hope that it will be useful,
14  * but WITHOUT ANY WARRANTY; without even the implied warranty of
15  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  * GNU General Public License for more details.
17  *
18  * You should have received a copy of the GNU General Public License
19  * along with LORENE; if not, write to the Free Software
20  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
21  *
22  */
23 
24 
25 char scalar_raccord_C[] = "$Header: /cvsroot/Lorene/C++/Source/Tensor/Scalar/scalar_raccord.C,v 1.4 2014/10/13 08:53:47 j_novak Exp $" ;
26 
27 /*
28  * $Id: scalar_raccord.C,v 1.4 2014/10/13 08:53:47 j_novak Exp $
29  * $Log: scalar_raccord.C,v $
30  * Revision 1.4 2014/10/13 08:53:47 j_novak
31  * Lorene classes and functions now belong to the namespace Lorene.
32  *
33  * Revision 1.3 2014/10/06 15:16:16 j_novak
34  * Modified #include directives to use c++ syntax.
35  *
36  * Revision 1.2 2003/10/01 13:04:44 e_gourgoulhon
37  * The method Tensor::get_mp() returns now a reference (and not
38  * a pointer) onto a mapping.
39  *
40  * Revision 1.1 2003/09/25 08:58:10 e_gourgoulhon
41  * First version.
42  *
43  *
44  * $Header: /cvsroot/Lorene/C++/Source/Tensor/Scalar/scalar_raccord.C,v 1.4 2014/10/13 08:53:47 j_novak Exp $
45  *
46  */
47 
48 //standard
49 #include <cstdlib>
50 
51 // LORENE
52 #include "tensor.h"
53 #include "proto.h"
54 #include "matrice.h"
55 
56 namespace Lorene {
57 Matrice matrice_raccord_pair (int cont, double alpha_kernel) ;
58 Matrice matrice_raccord_impair (int cont, double alpha_kernel) ;
59 Tbl sec_membre_raccord (Tbl coef, int cont, double alpha_shell) ;
60 Tbl regularise (Tbl coef, int nr, int base_r) ;
61 
62 void Scalar::raccord (int aux) {
63 
64  assert (etat != ETATNONDEF) ;
65 
66  assert (aux >=0) ;
67  int cont = aux+1 ;
68 
69  const Map_af* mapping = dynamic_cast<const Map_af*>( mp ) ;
70 
71  if (mapping == 0x0) {
72  cout <<
73  "Scalar::raccord : The mapping does not belong to the class Map_af !"
74  << endl ;
75  abort() ;
76  }
77 
78  assert (mapping->get_mg()->get_type_r(1) == FIN) ;
79  assert (mapping->get_mg()->get_type_r(0) == RARE) ;
80 
81  // On passe en Ylm et vire tout dans la zone interne...
82  va.coef() ;
83  va.ylm() ;
85  va.c_cf->t[0]->annule_hard() ;
86 
87  // Confort :
88  int nz = mapping->get_mg()->get_nzone() ;
89  int nbrer_kernel = mapping->get_mg()->get_nr(0) ;
90  int nbrer_shell = mapping->get_mg()->get_nr(1) ;
91 
92  int nbret_kernel = mapping->get_mg()->get_nt(0) ;
93  int nbret_shell = mapping->get_mg()->get_nt(1) ;
94 
95  int nbrep_kernel = mapping->get_mg()->get_np(0) ;
96  int nbrep_shell = mapping->get_mg()->get_np(1) ;
97 
98  double alpha_kernel = mapping->get_alpha()[0] ;
99  double alpha_shell = mapping->get_alpha()[1] ;
100 
101  int base_r, m_quant, l_quant ;
102 
103  for (int k=0 ; k<nbrep_kernel+1 ; k++)
104  for (int j=0 ; j<nbret_kernel ; j++)
105  if (nullite_plm(j, nbret_kernel, k,nbrep_kernel, va.base) == 1)
106  if (nullite_plm(j, nbret_shell, k, nbrep_shell, va.base) == 1)
107  {
108  // calcul des nombres quantiques :
109  donne_lm(nz, 0, j, k, va.base, m_quant, l_quant, base_r) ;
110  assert ((base_r == R_CHEBP) || (base_r == R_CHEBI)) ;
111 
112  Matrice systeme(cont, cont) ;
113 
114  Tbl facteur (nbrer_kernel) ;
115  facteur.annule_hard() ;
116  for (int i=0 ; i<nbrer_shell ; i++)
117  if (i<nbrer_kernel)
118  facteur.set(i) = (*va.c_cf)(1, k, j, i) ;
119 
120  Tbl sec_membre (sec_membre_raccord (facteur, cont, alpha_shell)) ;
121 
122  if (base_r == R_CHEBP)
123  systeme = matrice_raccord_pair (cont, alpha_kernel) ;
124  else
125  systeme = matrice_raccord_impair (cont, alpha_kernel) ;
126 
127  Tbl soluce (systeme.inverse(sec_membre)) ;
128  Tbl regulier (nbrer_kernel) ;
129 
130  if (l_quant == 0)
131  for (int i=0 ; i<cont ; i++)
132  va.c_cf->set(0, k, j, i) = soluce(i) ;
133  else {
134  if (l_quant %2 == 0)
135  regulier = regularise (soluce, nbrer_kernel, R_CHEBP) ;
136  else
137  regulier = regularise (soluce, nbrer_kernel, R_CHEBI) ;
138 
139  for (int i=0 ; i<nbrer_kernel ; i++)
140  va.c_cf->set(0, k, j, i) = regulier(i) ;
141  }
142  }
143  va.ylm_i() ;
144 }
145 }
Lorene::Scalar::raccord
void raccord(int n)
Performs the matching of the nucleus with respect to the first shell.
Definition: scalar_raccord.C:62
Lorene::Tbl::set
double & set(int i)
Read/write of a particular element (index i) (1D case)
Definition: tbl.h:281
Lorene::Mg3d::get_np
int get_np(int l) const
Returns the number of points in the azimuthal direction ( ) in domain no. l.
Definition: grilles.h:462
Lorene::Valeur::set_etat_cf_qcq
void set_etat_cf_qcq()
Sets the logical state to ETATQCQ (ordinary state) for values in the configuration space (Mtbl_cf c_c...
Definition: valeur.C:712
Lorene::Scalar::va
Valeur va
The numerical value of the Scalar
Definition: scalar.h:405
Lorene::Mtbl_cf::t
Tbl ** t
Array (size nzone ) of pointers on the Tbl 's which contain the spectral coefficients in each domain.
Definition: mtbl_cf.h:205
Lorene::Valeur::c_cf
Mtbl_cf * c_cf
Coefficients of the spectral expansion of the function.
Definition: valeur.h:302
Lorene::Mg3d::get_nt
int get_nt(int l) const
Returns the number of points in the co-latitude direction ( ) in domain no. l.
Definition: grilles.h:457
Lorene
Lorene prototypes.
Definition: app_hor.h:64
Lorene::Tbl::annule_hard
void annule_hard()
Sets the Tbl to zero in a hard way.
Definition: tbl.C:372
Lorene::Map::get_mg
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:765
Lorene::Mg3d::get_type_r
int get_type_r(int l) const
Returns the type of sampling in the radial direction in domain no.
Definition: grilles.h:474
Lorene::Mtbl_cf::set
Tbl & set(int l)
Read/write of the Tbl containing the coefficients in a given domain.
Definition: mtbl_cf.h:294
R_CHEBP
#define R_CHEBP
base de Cheb. paire (rare) seulement
Definition: type_parite.h:168
R_CHEBI
#define R_CHEBI
base de Cheb. impaire (rare) seulement
Definition: type_parite.h:170
Lorene::Valeur::coef
void coef() const
Computes the coeffcients of *this.
Definition: valeur_coef.C:148
Lorene::Scalar::etat
int etat
The logical state ETATNONDEF (undefined), ETATZERO (null), ETATUN (one), or ETATQCQ (ordinary).
Definition: scalar.h:396
Lorene::Matrice::inverse
Tbl inverse(const Tbl &sec_membre) const
Solves the linear system represented by the matrix.
Definition: matrice.C:424
Lorene::Tensor::mp
const Map *const mp
Mapping on which the numerical values at the grid points are defined.
Definition: tensor.h:295
Lorene::Tbl
Basic array class.
Definition: tbl.h:161
Lorene::Mg3d::get_nzone
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:448
Lorene::Mg3d::get_nr
int get_nr(int l) const
Returns the number of points in the radial direction ( ) in domain no. l.
Definition: grilles.h:452
Lorene::Map_af
Affine radial mapping.
Definition: map.h:2027
Lorene::Matrice
Matrix handling.
Definition: matrice.h:152
Lorene::Valeur::base
Base_val base
Bases on which the spectral expansion is performed.
Definition: valeur.h:305
Lorene::Valeur::ylm
void ylm()
Computes the coefficients of *this.
Definition: valeur_ylm.C:138
Lorene::Map_af::get_alpha
const double * get_alpha() const
Returns the pointer on the array alpha.
Definition: map_af.C:477
Lorene::Valeur::ylm_i
void ylm_i()
Inverse of ylm()
Definition: valeur_ylm_i.C:131