00001 subroutine sourceterm_trG_peos(sou)
00002   use phys_constant, only : long, pi
00003   use grid_parameter, only : nrg, ntg, npg
00004   use def_metric, only : psi, alph, alps, 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, rp2s, utgc, zfac
00011   real(long) :: psigc, alpgc, 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         rp2s = 3.0d0*hhgc*rhogc*(alpgc*utgc)**2   &
00033       &      - 2.0d0*hhgc*rhogc + 5.0d0*pregc
00034 
00035         sou(irg,itg,ipg) = + 0.875d0*alpgc*psigc**5*aijaij  &
00036       &                    + radi**2*2.0d0*pi*alpgc*psigc**5*rp2s*zfac
00037       end do
00038     end do
00039   end do
00040 end subroutine sourceterm_trG_peos