bh_boundary_n_Bfun.f90

Go to the documentation of this file.
00001 subroutine bh_boundary_n_Bfun(dsou_surf,potx,poty,potz)
00002   use phys_constant, only  : long
00003   use grid_parameter, only : nrg, ntg, npg, rgin
00004   use trigonometry_grav_theta, only : hsinthg, hcosthg
00005   use trigonometry_grav_phi,   only : hsinphig, hcosphig
00006   use def_binary_parameter,    only : sepa
00007   use def_metric
00008   implicit none
00009   real(long), pointer :: dsou_surf(:,:), potx(:,:,:), poty(:,:,:), potz(:,:,:)
00010   real(long) :: st, ct, sp, cp, xa,ya,za, rcm2, xycm2, rcm, tcm, pcm, xcm, ycm, zcm
00011   real(long) :: fx,fy,gx,gy
00012   integer    :: itg, ipg
00013 !
00014   do ipg = 1, npg
00015     do itg = 1, ntg
00016       st = hsinthg(itg)
00017       ct = hcosthg(itg)
00018       sp = hsinphig(ipg)
00019       cp = hcosphig(ipg)
00020 ! Coordinates wrt the center of black hole
00021       xa = rgin*st*cp
00022       ya = rgin*st*sp
00023       za = rgin*ct
00024 ! Coordinates wrt the CM
00025       xcm = xa - 0.5d0*sepa
00026       ycm = ya
00027       zcm = za
00028  
00029       gx = 0.25d0*(potx(0,itg,ipg)   + potx(0,itg-1,ipg)  &
00030                  + potx(0,itg,ipg-1) + potx(0,itg-1,ipg-1) )
00031       gy = 0.25d0*(poty(0,itg,ipg)   + poty(0,itg-1,ipg)  &
00032                  + poty(0,itg,ipg-1) + poty(0,itg-1,ipg-1) )
00033 
00034       fx = 0.25d0*(bvxd(0,itg,ipg)   + bvxd(0,itg-1,ipg)  &
00035                  + bvxd(0,itg,ipg-1) + bvxd(0,itg-1,ipg-1) )
00036       fy = 0.25d0*(bvyd(0,itg,ipg)   + bvyd(0,itg-1,ipg)  &
00037                  + bvyd(0,itg,ipg-1) + bvyd(0,itg-1,ipg-1) )
00038 
00039       dsou_surf(itg,ipg) = 0.0d0
00040 !      dsou_surf(itg,ipg) = 2.0d0*0.08*(-ycm*st*cp + xcm*st*sp)
00041 !      dsou_surf(itg,ipg) = 4.0d0*(gx-fx)*st*cp + 4.0d0*(gy-fy)*st*sp
00042     end do
00043   end do
00044 !
00045 end subroutine bh_boundary_n_Bfun