00001 subroutine sourceterm_HaC_peos(sou)
00002 use phys_constant, only : long, pi
00003 use grid_parameter, only : nrg, ntg, npg
00004 use def_metric, only : psi, alph, tfkijkij
00005 use def_matter, only : emdg, jomeg_int
00006 use def_matter_parameter, only : ber, radi
00007 use interface_interpo_linear_type0
00008 implicit none
00009 real(long), pointer :: sou(:,:,:)
00010 real(long) :: emdgc, rhogc, pregc, hhgc, utgc, rhoHc, zfac
00011 real(long) :: psigc, alpgc, alpsigc, aijaij, ene, jomeg_intgc
00012 integer :: irg, itg, ipg
00013
00014
00015
00016
00017 do ipg = 1, npg
00018 do itg = 1, ntg
00019 do irg = 1, nrg
00020 call interpo_linear_type0(emdgc,emdg,irg,itg,ipg)
00021 call interpo_linear_type0(jomeg_intgc,jomeg_int,irg,itg,ipg)
00022 call interpo_linear_type0(psigc,psi,irg,itg,ipg)
00023 call interpo_linear_type0(alpgc,alph,irg,itg,ipg)
00024 aijaij = tfkijkij(irg,itg,ipg)
00025 zfac = 1.0d0
00026 if (emdgc <= 1.0d-15) then
00027 emdgc = 1.0d-15
00028 zfac = 0.0d0
00029 end if
00030 call peos_q2hprho(emdgc, hhgc, pregc, rhogc, ene)
00031 utgc = hhgc/ber*exp(jomeg_intgc)
00032 rhoHc = hhgc*rhogc*(alpgc*utgc)**2 - pregc
00033
00034 sou(irg,itg,ipg) = - 0.125d0*psigc**5*aijaij &
00035 & - radi**2*2.0d0*pi*psigc**5*rhoHc*zfac
00036 end do
00037 end do
00038 end do
00039 end subroutine sourceterm_HaC_peos
00040