00001 subroutine sourceterm_MWtemp_CF(sou)
00002 use phys_constant, only : long
00003 use grid_parameter, only : nrg, ntg, npg
00004 use def_metric, only : alph, psi
00005 use def_faraday_tensor, only : fxd, fyd, fzd
00006 use def_emfield_derivatives, only : Lie_bAxd_grid, Lie_bAyd_grid, &
00007 & Lie_bAzd_grid
00008 use interface_interpo_linear_type0
00009 use interface_grgrad_midpoint
00010 use make_array_3d
00011 implicit none
00012 real(long), pointer :: sou(:,:,:)
00013 real(long), pointer :: fnc2(:,:,:)
00014 real(long), pointer :: grad2x(:,:,:), grad2y(:,:,:), grad2z(:,:,:)
00015 real(long), pointer :: dxLbAx(:,:,:), dyLbAx(:,:,:), dzLbAx(:,:,:)
00016 real(long), pointer :: dxLbAy(:,:,:), dyLbAy(:,:,:), dzLbAy(:,:,:)
00017 real(long), pointer :: dxLbAz(:,:,:), dyLbAz(:,:,:), dzLbAz(:,:,:)
00018 real(long) :: psigc, alphgc, a2divp2, dxLbAxgc, dyLbAygc, dzLbAzgc
00019 real(long) :: fxdgc, fydgc, fzdgc
00020 real(long) :: dxps2al, dyps2al, dzps2al, dps2alfa
00021 integer :: irg, itg, ipg
00022
00023
00024
00025
00026 call alloc_array3d(fnc2,0,nrg,0,ntg,0,npg)
00027 call alloc_array3d(grad2x,1,nrg,1,ntg,1,npg)
00028 call alloc_array3d(grad2y,1,nrg,1,ntg,1,npg)
00029 call alloc_array3d(grad2z,1,nrg,1,ntg,1,npg)
00030 call alloc_array3d(dxLbAx,1,nrg,1,ntg,1,npg)
00031 call alloc_array3d(dyLbAx,1,nrg,1,ntg,1,npg)
00032 call alloc_array3d(dzLbAx,1,nrg,1,ntg,1,npg)
00033 call alloc_array3d(dxLbAy,1,nrg,1,ntg,1,npg)
00034 call alloc_array3d(dyLbAy,1,nrg,1,ntg,1,npg)
00035 call alloc_array3d(dzLbAy,1,nrg,1,ntg,1,npg)
00036 call alloc_array3d(dxLbAz,1,nrg,1,ntg,1,npg)
00037 call alloc_array3d(dyLbAz,1,nrg,1,ntg,1,npg)
00038 call alloc_array3d(dzLbAz,1,nrg,1,ntg,1,npg)
00039 fnc2(0:nrg,0:ntg,0:npg) = psi(0:nrg,0:ntg,0:npg)**2/alph(0:nrg,0:ntg,0:npg)
00040 call grgrad_midpoint(fnc2,grad2x,grad2y,grad2z)
00041 call grgrad_midpoint(Lie_bAxd_grid,dxLbAx,dyLbAx,dzLbAx)
00042 call grgrad_midpoint(Lie_bAyd_grid,dxLbAy,dyLbAy,dzLbAy)
00043 call grgrad_midpoint(Lie_bAzd_grid,dxLbAz,dyLbAz,dzLbAz)
00044
00045 do ipg = 1, npg
00046 do itg = 1, ntg
00047 do irg = 1, nrg
00048 call interpo_linear_type0(psigc,psi,irg,itg,ipg)
00049 call interpo_linear_type0(alphgc,alph,irg,itg,ipg)
00050 a2divp2 = alphgc**2/psigc**2
00051 dxps2al = grad2x(irg,itg,ipg)
00052 dyps2al = grad2y(irg,itg,ipg)
00053 dzps2al = grad2z(irg,itg,ipg)
00054 fxdgc = fxd(irg,itg,ipg)
00055 fydgc = fyd(irg,itg,ipg)
00056 fzdgc = fzd(irg,itg,ipg)
00057 dxLbAxgc = dxLbAx(irg,itg,ipg)
00058 dyLbAygc = dyLbAy(irg,itg,ipg)
00059 dzLbAzgc = dzLbAz(irg,itg,ipg)
00060
00061 dps2alfa = a2divp2*(dxps2al*fxdgc+dyps2al*fydgc+dzps2al*fzdgc)
00062
00063 sou(irg,itg,ipg) = dps2alfa + dxLbAxgc + dyLbAygc + dzLbAzgc
00064
00065 end do
00066 end do
00067 end do
00068
00069 deallocate(fnc2)
00070 deallocate(grad2x)
00071 deallocate(grad2y)
00072 deallocate(grad2z)
00073 deallocate(dxLbAx)
00074 deallocate(dyLbAx)
00075 deallocate(dzLbAx)
00076 deallocate(dxLbAy)
00077 deallocate(dyLbAy)
00078 deallocate(dzLbAy)
00079 deallocate(dxLbAz)
00080 deallocate(dyLbAz)
00081 deallocate(dzLbAz)
00082
00083 end subroutine sourceterm_MWtemp_CF
00084