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