23char separation_C[] =
"$Header: /cvsroot/Lorene/C++/Source/Bhole_binaire/separation.C,v 1.4 2014/10/13 08:52:41 j_novak Exp $" ;
74void separation (
const Cmp& c1,
const Cmp& c2, Cmp& res1, Cmp& res2,
int decrois,
75 int puiss,
int lmax,
double precision,
const double relax,
const int itemax,
const int flag) {
77 assert (c1.get_etat() != ETATNONDEF) ;
78 assert (c2.get_etat() != ETATNONDEF) ;
80 if ((c1.get_etat() == ETATZERO) && (c2.get_etat() == ETATZERO)) {
81 res1.set_etat_zero() ;
82 res2.set_etat_zero() ;
88 if (res1.get_etat() == ETATZERO) {
90 res1.std_base_scal() ;
93 res1.raccord_externe (decrois, puiss, lmax) ;
94 for (
int i=0 ; i<decrois ; i++)
98 if (res2.get_etat() == ETATZERO) {
100 res2.std_base_scal() ;
102 res2.raccord_externe (decrois, puiss, lmax) ;
103 for (
int i=0 ; i<decrois ; i++)
111 Cmp old_deux (res2) ;
114 Mtbl xa_mtbl_un (c1.get_mp()->xa) ;
115 Mtbl ya_mtbl_un (c1.get_mp()->ya) ;
116 Mtbl za_mtbl_un (c1.get_mp()->za) ;
118 Mtbl xa_mtbl_deux (c2.get_mp()->xa) ;
119 Mtbl ya_mtbl_deux (c2.get_mp()->ya) ;
120 Mtbl za_mtbl_deux (c2.get_mp()->za) ;
122 double xabs, yabs, zabs, air, theta, phi ;
126 int nz_un = c1.get_mp()->get_mg()->get_nzone() ;
127 for (
int l=1 ; l<nz_un-1 ; l++) {
129 np = c1.get_mp()->get_mg()->get_np(l) ;
130 nt = c1.get_mp()->get_mg()->get_nt(l) ;
131 nr = c1.get_mp()->get_mg()->get_nr(l) ;
133 for (
int k=0 ; k<np ; k++)
134 for (
int j=0 ; j<nt ; j++)
135 for (
int i=0 ; i<nr ; i++) {
136 xabs = xa_mtbl_un (l, k, j, i) ;
137 yabs = ya_mtbl_un (l, k, j, i) ;
138 zabs = za_mtbl_un (l, k, j, i) ;
140 c2.get_mp()->convert_absolute(xabs, yabs, zabs, air, theta, phi) ;
141 res1.set(l, k, j, i) =
142 (1-relax)*res1.set(l, k, j, i) +
143 relax*(c1(l, k, j, i) - old_deux.val_point(air, theta, phi)) ;
148 int nz_deux = c2.get_mp()->get_mg()->get_nzone() ;
149 for (
int l=1 ; l<nz_deux-1 ; l++) {
151 np = c2.get_mp()->get_mg()->get_np(l) ;
152 nt = c2.get_mp()->get_mg()->get_nt(l) ;
153 nr = c2.get_mp()->get_mg()->get_nr(l) ;
155 for (
int k=0 ; k<np ; k++)
156 for (
int j=0 ; j<nt ; j++)
157 for (
int i=0 ; i<nr ; i++) {
159 xabs = xa_mtbl_deux (l, k, j, i) ;
160 yabs = ya_mtbl_deux (l, k, j, i) ;
161 zabs = za_mtbl_deux (l, k, j, i) ;
163 c1.get_mp()->convert_absolute(xabs, yabs, zabs, air, theta, phi) ;
164 res2.set(l, k, j, i) =
165 (1-relax)*res2.set(l, k, j, i) +
166 relax*(c2(l, k, j, i) - old_un.val_point(air, theta, phi)) ;
171 res1.va.set_etat_c_qcq() ;
172 res2.va.set_etat_c_qcq() ;
174 res1.raccord_externe (decrois, puiss, lmax) ;
175 for (
int i=0 ; i<decrois ; i++)
179 res2.raccord_externe (decrois, puiss, lmax) ;
180 for (
int i=0 ; i<decrois ; i++)
189 for (
int i=1 ; i<nz_un-1 ; i++)
190 if (diff_un(i)>erreur)
191 erreur = diff_un(i) ;
194 for (
int i=1 ; i<nz_deux-1 ; i++)
195 if (diff_deux(i)>erreur)
196 erreur = diff_deux(i) ;
199 cout <<
"Pas " << conte <<
" : erreur = " << erreur << endl ;
200 if (erreur<=precision)
Tbl diffrelmax(const Cmp &a, const Cmp &b)
Relative difference between two Cmp (max version).