23char raccord_c1_C[] =
"$Header: /cvsroot/Lorene/C++/Source/Non_class_members/Utilities/raccord_c1.C,v 1.4 2014/10/13 08:53:32 j_novak Exp $" ;
53Cmp raccord_c1(
const Cmp& uu,
int l1) {
55 const Map_radial* mpi =
dynamic_cast<const Map_radial*
>( uu.get_mp() ) ;
59 "raccord_c1 : The mapping does not belong to the class Map_radial !"
64 const Map_radial& mp = *mpi ;
66 const Mg3d& mg = *(mp.get_mg()) ;
69 int nz = mg.get_nzone() ;
79 Mtbl r2 = mp.r * mp.r ;
83 Valeur& va = resu.va ;
88 va.c->set_etat_qcq() ;
89 va.c->t[l1]->set_etat_qcq() ;
91 int np = mg.get_np(l1) ;
92 int nt = mg.get_nt(l1) ;
93 assert( mg.get_np(l0) == np ) ;
94 assert( mg.get_nt(l0) == nt ) ;
95 assert( mg.get_np(l2) == np ) ;
96 assert( mg.get_nt(l2) == nt ) ;
98 int nr0 = mg.get_nr(l0) ;
99 int nr1 = mg.get_nr(l1) ;
101 for (
int k=0; k<np; k++) {
102 for (
int j=0; j<nt; j++) {
103 double lam0 = uu(l0, k, j, nr0-1) ;
104 double lam1 = uu.dsdr()(l0, k, j, nr0-1) / dxdr(l1, k, j, 0) ;
105 double mu0 = uu(l2, k, j, 0) ;
106 double mu1 = uu.dsdr()(l2, k, j, 0) / dxdr(l1, k, j, nr1-1) ;
108 if ( mg.get_type_r(l2) == UNSURR ) {
109 mu1 /= r2(l2, k, j, 0) ;
112 double s0 = 0.25 * (mu0 + lam0) ;
113 double s1 = 0.25 * (mu1 + lam1) ;
114 double d0 = 0.25 * (mu0 - lam0) ;
115 double d1 = 0.25 * (mu1 - lam1) ;
117 double a0 = 2. * s0 - d1 ;
118 double a1 = 3. * d0 - s1 ;
120 double a3 = s1 - d0 ;
122 for (
int i=0; i<nr1; i++) {
123 double x = mg.get_grille3d(l1)->x[i] ;
124 va.set(l1, k, j, i) = a0 + x * ( a1 + x * ( a2 + x * a3 ) ) ;
137Scalar raccord_c1(
const Scalar& uu,
int l1) {
140 return Scalar( raccord_c1(cuu, l1) ) ;