00001 subroutine source_ang_mom_qeos(souf)
00002 use phys_constant, only : long, pi
00003 use grid_parameter, only : nrg, ntg, npg, nrf, ntf, npf
00004 use def_matter, only : rhof, rs, utf, omef
00005 use def_matter_parameter
00006 use def_metric, only : tfkijkij, psi, alph, &
00007 bvxd, bvyd, bvzd
00008 use coordinate_grav_r, only : rg
00009 use trigonometry_grav_theta, only : sinthg, costhg
00010 use trigonometry_grav_phi, only : sinphig, cosphig
00011 use def_vector_phi, only : vec_phif
00012 use make_array_3d
00013 use interface_interpo_gr2fl
00014 implicit none
00015 real(long), pointer :: souf(:,:,:)
00016 integer :: ir,it,ip
00017 real(long) :: alphw, psiw, rhow, prew, hhw, utw, omew, ene
00018 real(long) :: rjjx, rjjy, rjjz, rjjphi, dummy
00019 real(long) :: zfac, small = 1.0d-15
00020 real(long) :: vphif(1:3)
00021 real(long) :: otermx, otermy, otermz, bvxufw, bvyufw, bvzufw
00022 real(long), pointer :: alphf(:,:,:), psif(:,:,:)
00023 real(long), pointer :: bvxdf(:,:,:), bvydf(:,:,:), bvzdf(:,:,:)
00024
00025 call alloc_array3d(psif, 0, nrf, 0, ntf, 0, npf)
00026 call alloc_array3d(alphf, 0, nrf, 0, ntf, 0, npf)
00027 call alloc_array3d(bvxdf, 0, nrf, 0, ntf, 0, npf)
00028 call alloc_array3d(bvydf, 0, nrf, 0, ntf, 0, npf)
00029 call alloc_array3d(bvzdf, 0, nrf, 0, ntf, 0, npf)
00030
00031 call interpo_gr2fl(alph, alphf)
00032 call interpo_gr2fl(psi , psif )
00033 call interpo_gr2fl(bvxd, bvxdf)
00034 call interpo_gr2fl(bvyd, bvydf)
00035 call interpo_gr2fl(bvzd, bvzdf)
00036
00037 do ip = 0, npf
00038 do it = 0, ntf
00039 do ir = 0, nrf
00040 rhow = rhof(ir,it,ip)
00041 psiw = psif(ir,it,ip)
00042 alphw = alphf(ir,it,ip)
00043 bvxufw = bvxdf(ir,it,ip)
00044 bvyufw = bvydf(ir,it,ip)
00045 bvzufw = bvzdf(ir,it,ip)
00046 call quark_rho2phenedpdrho(rhow, prew, hhw, ene, dummy)
00047
00048 utw = utf(ir,it,ip)
00049 omew = omef(ir,it,ip)
00050
00051 vphif(1) = vec_phif(ir,it,ip,1)
00052 vphif(2) = vec_phif(ir,it,ip,2)
00053 vphif(3) = vec_phif(ir,it,ip,3)
00054
00055 otermx = bvxufw + omew*vphif(1)
00056 otermy = bvyufw + omew*vphif(2)
00057 otermz = bvzufw + omew*vphif(3)
00058
00059 zfac = 1.0d0
00060 if (rhow <= rhos_qs) zfac = 0.0d0
00061 rjjx = hhw*rhow*alphw*utw**2*psiw**4*otermx
00062 rjjy = hhw*rhow*alphw*utw**2*psiw**4*otermy
00063 rjjz = hhw*rhow*alphw*utw**2*psiw**4*otermz
00064
00065 rjjphi = rjjx*vphif(1) + rjjy*vphif(2) + rjjz*vphif(3)
00066 souf(ir,it,ip) = rjjphi*psiw**6*zfac
00067
00068 end do
00069 end do
00070 end do
00071
00072 deallocate(alphf)
00073 deallocate(psif)
00074 deallocate(bvxdf)
00075 deallocate(bvydf)
00076 deallocate(bvzdf)
00077 end subroutine source_ang_mom_qeos
