23 char bhole_kij_C[] =
"$Header: /cvsroot/Lorene/C++/Source/Bhole_binaire/bhole_kij.C,v 1.6 2014/10/13 08:52:40 j_novak Exp $" ;
83 #include "utilitaires.h"
84 #include "graphique.h"
110 assert ((orientation_un==0) || (orientation_un==M_PI)) ;
112 assert ((orientation_deux==0) || (orientation_deux==M_PI)) ;
113 int same_orient = (orientation_un == orientation_deux) ? 1 : -1 ;
160 ntot_un = division_xpun (ntot_un, 0) ;
164 ntot_deux = division_xpun (ntot_deux, 0) ;
168 for (
int lig = 0 ; lig<3 ; lig++)
169 for (
int col = lig ; col<3 ; col++) {
177 if (same_orient == 1)
183 auxi_un = division_xpun (auxi_un, 0) ;
184 auxi_un = auxi_un / ntot_un ;
190 auxi_deux = division_xpun (auxi_deux, 0) ;
191 auxi_deux = auxi_deux / ntot_deux ;
213 double xabs, yabs, zabs, air, theta, phi ;
216 for (
int l=2 ; l<nz_un ; l++) {
226 for (
int k=0 ; k<np ; k++)
227 for (
int j=0 ; j<nt ; j++)
228 for (
int i=0 ; i<nr ; i++) {
230 xabs = xabs_un (l, k, j, i) ;
231 yabs = yabs_un (l, k, j, i) ;
232 zabs = zabs_un (l, k, j, i) ;
236 (xabs, yabs, zabs, air, theta, phi) ;
240 auxi_un.
set(l, k, j, i) =
241 copie_un(l, k, j, i) / ntot_un (l, k, j, i)/2. ;
244 auxi_un.
set(l, k, j, i) =
245 ind * auxi_deux.
val_point (air, theta, phi) ;
250 for (
int k=0 ; k<np ; k++)
251 for (
int j=0 ; j<nt ; j++)
252 auxi_un.
set(nz_un-1, k, j, nr-1) = 0 ;
256 for (
int l=2 ; l<nz_deux ; l++) {
266 for (
int k=0 ; k<np ; k++)
267 for (
int j=0 ; j<nt ; j++)
268 for (
int i=0 ; i<nr ; i++) {
270 xabs = xabs_deux (l, k, j, i) ;
271 yabs = yabs_deux (l, k, j, i) ;
272 zabs = zabs_deux (l, k, j, i) ;
276 (xabs, yabs, zabs, air, theta, phi) ;
280 auxi_deux.
set(l, k, j, i) =
281 copie_deux(l, k, j, i) / ntot_deux (l, k, j, i) /2.;
284 auxi_deux.
set(l, k, j, i) =
285 ind * auxi_un.
val_point (air, theta, phi) ;
289 for (
int k=0 ; k<np ; k++)
290 for (
int j=0 ; j<nt ; j++)
291 auxi_deux.
set(nz_deux-1, k, j, nr-1) = 0 ;
304 for (
int lig=0 ; lig<3 ; lig++)
305 for (
int col=lig ; col<3 ; col++) {
321 double lim_un = distance/2. ;
322 double lim_deux = distance/2. ;
323 double int_un = distance/6. ;
324 double int_deux = distance/6. ;
328 fonction_f_un = 0.5*
pow(
329 cos((
hole1.
mp.
r-int_un)*M_PI/2./(lim_un-int_un)), 2.)+0.5 ;
333 fonction_g_un = 0.5*
pow
334 (
sin((
hole1.
mp.
r-int_un)*M_PI/2./(lim_un-int_un)), 2.) ;
338 fonction_f_deux = 0.5*
pow(
339 cos((
hole2.
mp.
r-int_deux)*M_PI/2./(lim_deux-int_deux)), 2.)+0.5 ;
343 fonction_g_deux = 0.5*
pow(
344 sin((
hole2.
mp.
r-int_deux)*M_PI/2./(lim_un-int_deux)), 2.) ;
361 double xabs, yabs, zabs, air_un, air_deux, theta, phi ;
364 for (
int l=0 ; l<nz_un ; l++) {
373 for (
int k=0 ; k<np ; k++)
374 for (
int j=0 ; j<nt ; j++)
375 for (
int i=0 ; i<nr ; i++) {
377 xabs = xabs_un (l, k, j, i) ;
378 yabs = yabs_un (l, k, j, i) ;
379 zabs = zabs_un (l, k, j, i) ;
383 (xabs, yabs, zabs, air_un, theta, phi) ;
385 (xabs, yabs, zabs, air_deux, theta, phi) ;
387 if (air_un <= lim_un)
389 decouple_un.
set(l, k, j, i) = 1 ;
392 decouple_un.
set(l, k, j, i) =
393 fonction_f_un (l, k, j, i) ;
395 if (air_deux <= lim_deux)
396 if (air_deux < int_deux)
397 decouple_un.
set(l, k, j, i) = 0 ;
400 decouple_un.
set(l, k, j, i) =
401 fonction_g_deux.
val_point (air_deux, theta, phi) ;
405 decouple_un.
set(l, k, j, i) = 0.5 ;
410 for (
int k=0 ; k<np ; k++)
411 for (
int j=0 ; j<nt ; j++)
412 decouple_un.
set(nz_un-1, k, j, nr) = 0.5 ;
415 for (
int l=0 ; l<nz_deux ; l++) {
425 for (
int k=0 ; k<np ; k++)
426 for (
int j=0 ; j<nt ; j++)
427 for (
int i=0 ; i<nr ; i++) {
429 xabs = xabs_deux (l, k, j, i) ;
430 yabs = yabs_deux (l, k, j, i) ;
431 zabs = zabs_deux (l, k, j, i) ;
435 (xabs, yabs, zabs, air_un, theta, phi) ;
437 (xabs, yabs, zabs, air_deux, theta, phi) ;
439 if (air_deux <= lim_deux)
440 if (air_deux < int_deux)
441 decouple_deux.
set(l, k, j, i) = 1 ;
444 decouple_deux.
set(l, k, j, i) =
445 fonction_f_deux (l, k, j, i) ;
447 if (air_un <= lim_un)
449 decouple_deux.
set(l, k, j, i) = 0 ;
452 decouple_deux.
set(l, k, j, i) =
453 fonction_g_un.
val_point (air_un, theta, phi) ;
457 decouple_deux.
set(l, k, j, i) = 0.5 ;
462 for (
int k=0 ; k<np ; k++)
463 for (
int j=0 ; j<nt ; j++)
464 decouple_deux.
set(nz_deux-1, k, j, nr) = 0.5 ;