00001 subroutine calc_ergo
00002 use phys_constant, only : long
00003 use def_metric
00004 use grid_parameter, only : nrg, ntg, npg
00005 use def_matter, only : ergoin, ergoout
00006 implicit none
00007 integer :: itg, ipg, irg
00008 real(long) :: gtt
00009 do ipg = 0, npg
00010 do itg = 0, ntg
00011 ergoin(itg,ipg)=0
00012 ergoout(itg,ipg)=0
00013 do irg = 0, nrg
00014 gtt = psi(irg,itg,ipg)**4*(bvxd(irg,itg,ipg)*bvxd(irg,itg,ipg)+ &
00015 & bvyd(irg,itg,ipg)*bvyd(irg,itg,ipg) + &
00016 & bvzd(irg,itg,ipg)*bvzd(irg,itg,ipg)) - &
00017 & alph(irg,itg,ipg)**2
00018 if (gtt.gt.0.and.ergoin(itg,ipg).eq.0) then
00019 ergoin(itg,ipg)=irg
00020 ergoout(itg,ipg)=irg
00021 end if
00022 if (gtt.gt.0.and.ergoin(itg,ipg).gt.0) then
00023 ergoout(itg,ipg)=irg
00024 end if
00025 end do
00026 end do
00027 end do
00028 end subroutine calc_ergo