00001 subroutine calc_TrpBH
00002 use phys_constant, only : long
00003 use grid_parameter, only : nrg, ntg, npg
00004 use coordinate_grav_r, only : rg, hrg
00005 use def_metric, only : psi
00006 use def_metric_pBH, only : R_trpBH, hR_trpBH, psi_trpBH, hpsi_trpBH, &
00007 & alph_trpBH, halph_trpBH, bvr_trpBH, hbvr_trpBH
00008 use def_bh_parameter, only : mass_pBH
00009 implicit none
00010 integer :: irg, itg, ipg
00011 real(long) :: rgc, rini, rval, faca, facb, MoverR
00012
00013
00014
00015 faca = 27.0d0/16.0d0
00016 facb = 3.0d0*dsqrt(3.0d0)/4.0d0
00017 R_trpBH(0) = 3.0d0*mass_pBH/2.0d0
00018 psi_trpBH(0) = 1.0d+20
00019 alph_trpBH(0) = 0.0d+00
00020 bvr_trpBH(0) = 0.0d+00
00021 do irg = 1, nrg
00022 rini = R_trpBH(irg-1)
00023
00024 rgc = rg(irg)
00025 call calc_eqsolver_TrpBH(rgc,rini,rval)
00026 MoverR = mass_pBH/rval
00027 R_trpBH(irg) = rval
00028 psi_trpBH(irg) = dsqrt(rval/rgc)
00029 alph_trpBH(irg)= dsqrt(1.0d0 - 2.0d0*MoverR + faca*MoverR**4)
00030 bvr_trpBH(irg) = facb*MoverR**2*rgc/rval
00031
00032 rgc = hrg(irg)
00033 call calc_eqsolver_TrpBH(rgc,rini,rval)
00034 MoverR = mass_pBH/rval
00035 hR_trpBH(irg) = rval
00036 hpsi_trpBH(irg) = dsqrt(rval/rgc)
00037 halph_trpBH(irg)= dsqrt(1.0d0 - 2.0d0*MoverR + faca*MoverR**4)
00038 hbvr_trpBH(irg) = facb*MoverR**2*rgc/rval
00039 end do
00040
00041
00042
00043
00044
00045
00046
00047
00048
00049
00050 end subroutine calc_TrpBH