00001 subroutine bh_boundary_test(char_bc,sou_surf,dsou_surf)
00002   use phys_constant, only  : long
00003   use grid_parameter, only : nrg, ntg, npg, rgin
00004   use trigonometry_grav_theta, only : hsinthg
00005   use trigonometry_grav_phi,   only : hcosphig
00006   use def_binary_parameter,    only : sepa
00007   implicit none
00008   character(len=2), intent(in) :: char_bc
00009   real(long), pointer :: sou_surf(:,:), dsou_surf(:,:)
00010   real(long) :: st, cp, rad1, rad2, dr2dr1, bhmass
00011   integer    :: itg, ipg
00012 
00013   bhmass = 2.0d0*rgin
00014   do ipg = 1, npg
00015     do itg = 1, ntg
00016       st = hsinthg(itg)
00017       cp = hcosphig(ipg)
00018       rad1 = rgin
00019       rad2 = sqrt(rad1**2 - 2.0d0*rad1*sepa*st*cp + sepa**2)
00020       dr2dr1 = (rad1 - sepa*st*cp)/rad2
00021       if (char_bc.eq.'dh') then 
00022         sou_surf(itg,ipg) = 1.0d0 + 0.5d0*bhmass/rad1 + 0.5d0*bhmass/rad2
00023       end if
00024       if (char_bc.eq.'nh') then
00025         dsou_surf(itg,ipg) = - 0.5d0*bhmass/rad1**2 &
00026        &                     - 0.5d0*bhmass/rad2**2*dr2dr1
00027       end if
00028     end do
00029   end do
00030 
00031 end subroutine bh_boundary_test