00001 subroutine test_analytic_bns_solution_mpt(impt)
00002   use phys_constant, only  :   long, pi
00003   use grid_parameter, only  :   nrg, ntg, npg, nrf
00004   use def_metric, only  :   bvxd
00005   use coordinate_grav_r, only : rg
00006   use grid_points_binary_excision, only : rb
00007   implicit none
00008   integer     ::   irg,itg,ipg, impt
00009   real(long)  ::   zfac, small = 1.0d-15
00010 
00011   do ipg = 0, npg
00012     do itg = 0, ntg
00013       do irg = 0, nrg
00014       if(impt.eq.1) then
00015         if(irg.lt.nrf) &
00016         & bvxd(irg,itg,ipg) = - 2.0d0*pi/3.0d0*(3.0d0-rg(irg)**2) &
00017         &                     -(8.0d0*pi/3.0d0)*1.0d0/rb(irg,itg,ipg)
00018         if(irg.ge.nrf) &
00019         & bvxd(irg,itg,ipg) = -(4.0d0*pi/3.0d0)*1.0d0/rg(irg)  &
00020         &                     -(8.0d0*pi/3.0d0)*1.0d0/rb(irg,itg,ipg)
00021       end if
00022       if(impt.eq.2) then
00023         if(irg.lt.nrf) &
00024         & bvxd(irg,itg,ipg) = - 4.0d0*pi/3.0d0*(3.0d0-rg(irg)**2) &
00025         &                     -(4.0d0*pi/3.0d0)*1.0d0/rb(irg,itg,ipg)
00026         if(irg.ge.nrf) &
00027         & bvxd(irg,itg,ipg) = -(8.0d0*pi/3.0d0)*1.0d0/rg(irg)  &
00028         &                     -(4.0d0*pi/3.0d0)*1.0d0/rb(irg,itg,ipg)
00029       end if
00030 
00031 
00032 
00033 
00034 
00035 
00036 
00037 
00038 
00039 
00040       end do
00041     end do
00042   end do
00043 
00044 end subroutine test_analytic_bns_solution_mpt