00001 subroutine calc_vector_bh(sb)
00002 use phys_constant, only : long
00003 use grid_parameter, only : nrf, ntf, npf, nrg, ntg, npg
00004 use coordinate_grav_r, only : rg, hrg
00005 use trigonometry_grav_theta, only : sinthg, hsinthg, costhg, hcosthg
00006 use trigonometry_grav_phi, only : sinphig, cosphig, hsinphig, hcosphig
00007 use def_binary_parameter, only : dis
00008 use def_bh_parameter, only : th_spin_bh, phi_spin_bh
00009 use def_vector_bh
00010 use interface_interpo_linear_type0_2Dsurf
00011 implicit none
00012 real(long) :: rrgg, sth, cth, sphi, cphi, dis_bin
00013 real(long) :: spin_matrix(3,3), spin_axis(3), spin_z(3)
00014 integer :: irg, itg, ipg, sb, ii
00015
00016 dis_bin = dis
00017 if(sb.eq.1) dis_bin = 0.0d0
00018
00019 spin_matrix(1,1) = dcos(th_spin_bh)*dcos(phi_spin_bh)
00020 spin_matrix(1,2) = -dsin(phi_spin_bh)
00021 spin_matrix(1,3) = dsin(th_spin_bh)*dcos(phi_spin_bh)
00022 spin_matrix(2,1) = dcos(th_spin_bh)*dsin(phi_spin_bh)
00023 spin_matrix(2,2) = dcos(phi_spin_bh)
00024 spin_matrix(2,3) = dsin(th_spin_bh)*dsin(phi_spin_bh)
00025 spin_matrix(3,1) = -dsin(th_spin_bh)
00026 spin_matrix(3,2) = 0.0d0
00027 spin_matrix(3,3) = dcos(th_spin_bh)
00028 spin_z(1) = 0.0d0
00029 spin_z(2) = 0.0d0
00030 spin_z(3) = 1.0d0
00031 do ii = 1, 3
00032 spin_axis(ii) = spin_matrix(ii,1)*spin_z(1) &
00033 & + spin_matrix(ii,2)*spin_z(2) &
00034 & + spin_matrix(ii,3)*spin_z(3)
00035 end do
00036
00037 do ipg = 0, npg
00038 do itg = 0, ntg
00039 irg = 0
00040 rrgg = rg(irg)
00041 sth = sinthg(itg)
00042 cth = costhg(itg)
00043 sphi = sinphig(ipg)
00044 cphi = cosphig(ipg)
00045 vec_bh_cm_xg(itg,ipg,1) = - dis_bin + rrgg*sth*cphi
00046 vec_bh_cm_xg(itg,ipg,2) = rrgg*sth*sphi
00047 vec_bh_cm_xg(itg,ipg,3) = rrgg*cth
00048 vec_bh_cm_phig(itg,ipg,1) = - rrgg*sth*sphi
00049 vec_bh_cm_phig(itg,ipg,2) = - dis_bin + rrgg*sth*cphi
00050 vec_bh_cm_phig(itg,ipg,3) = 0.0d0
00051 vec_bh_cbh_xg(itg,ipg,1) = rrgg*sth*cphi
00052 vec_bh_cbh_xg(itg,ipg,2) = rrgg*sth*sphi
00053 vec_bh_cbh_xg(itg,ipg,3) = rrgg*cth
00054 vec_bh_cbh_phig(itg,ipg,1) = - rrgg*sth*sphi
00055 vec_bh_cbh_phig(itg,ipg,2) = rrgg*sth*cphi
00056 vec_bh_cbh_phig(itg,ipg,3) = 0.0d0
00057 vec_bh_cbh_spin(itg,ipg,1) = spin_axis(2)*vec_bh_cbh_xg(itg,ipg,3)&
00058 & - spin_axis(3)*vec_bh_cbh_xg(itg,ipg,2)
00059 vec_bh_cbh_spin(itg,ipg,2) = spin_axis(3)*vec_bh_cbh_xg(itg,ipg,1)&
00060 & - spin_axis(1)*vec_bh_cbh_xg(itg,ipg,3)
00061 vec_bh_cbh_spin(itg,ipg,3) = spin_axis(1)*vec_bh_cbh_xg(itg,ipg,2)&
00062 & - spin_axis(2)*vec_bh_cbh_xg(itg,ipg,1)
00063 if (itg.eq.0.or.ipg.eq.0) cycle
00064 rrgg = rg(irg)
00065 sth = hsinthg(itg)
00066 cth = hcosthg(itg)
00067 sphi = hsinphig(ipg)
00068 cphi = hcosphig(ipg)
00069 hvec_bh_cm_xg(itg,ipg,1) = - dis_bin + rrgg*sth*cphi
00070 hvec_bh_cm_xg(itg,ipg,2) = rrgg*sth*sphi
00071 hvec_bh_cm_xg(itg,ipg,3) = rrgg*cth
00072 hvec_bh_cm_phig(itg,ipg,1) = - rrgg*sth*sphi
00073 hvec_bh_cm_phig(itg,ipg,2) = - dis_bin + rrgg*sth*cphi
00074 hvec_bh_cm_phig(itg,ipg,3) = 0.0d0
00075 hvec_bh_cbh_xg(itg,ipg,1) = rrgg*sth*cphi
00076 hvec_bh_cbh_xg(itg,ipg,2) = rrgg*sth*sphi
00077 hvec_bh_cbh_xg(itg,ipg,3) = rrgg*cth
00078 hvec_bh_cbh_phig(itg,ipg,1) = - rrgg*sth*sphi
00079 hvec_bh_cbh_phig(itg,ipg,2) = rrgg*sth*cphi
00080 hvec_bh_cbh_phig(itg,ipg,3) = 0.0d0
00081 hvec_bh_cbh_spin(itg,ipg,1) = spin_axis(2)*hvec_bh_cbh_xg(itg,ipg,3)&
00082 & - spin_axis(3)*hvec_bh_cbh_xg(itg,ipg,2)
00083 hvec_bh_cbh_spin(itg,ipg,2) = spin_axis(3)*hvec_bh_cbh_xg(itg,ipg,1)&
00084 & - spin_axis(1)*hvec_bh_cbh_xg(itg,ipg,3)
00085 hvec_bh_cbh_spin(itg,ipg,3) = spin_axis(1)*hvec_bh_cbh_xg(itg,ipg,2)&
00086 & - spin_axis(2)*hvec_bh_cbh_xg(itg,ipg,1)
00087 end do
00088 end do
00089 end subroutine calc_vector_bh
