00001 subroutine sourceterm_volume_int_bbh_2pot_test(sou)
00002 use phys_constant, only : long, pi
00003 use grid_parameter, only : nrg, ntg, npg
00004 use def_metric, only : psi, alph
00005
00006
00007 use interface_interpo_linear_type0
00008 use interface_grgrad_midpoint
00009 use make_array_3d
00010 implicit none
00011 real(long), pointer :: sou(:,:,:)
00012 real(8), pointer :: psi_grad2x(:,:,:), psi_grad2y(:,:,:), psi_grad2z(:,:,:)
00013 real(8), pointer :: alph_grad2x(:,:,:),alph_grad2y(:,:,:),alph_grad2z(:,:,:)
00014
00015 real(long) :: psigc, alpgc, alpsigc, psiinv
00016 real(long) :: dxpsi,dypsi,dzpsi,dxalph,dyalph,dzalph
00017 integer :: irg, itg, ipg
00018
00019 call alloc_array3d(psi_grad2x,1,nrg,1,ntg,1,npg)
00020 call alloc_array3d(psi_grad2y,1,nrg,1,ntg,1,npg)
00021 call alloc_array3d(psi_grad2z,1,nrg,1,ntg,1,npg)
00022 call alloc_array3d(alph_grad2x,1,nrg,1,ntg,1,npg)
00023 call alloc_array3d(alph_grad2y,1,nrg,1,ntg,1,npg)
00024 call alloc_array3d(alph_grad2z,1,nrg,1,ntg,1,npg)
00025
00026 call grgrad_midpoint(psi,psi_grad2x,psi_grad2y,psi_grad2z)
00027 call grgrad_midpoint(alph,alph_grad2x,alph_grad2y,alph_grad2z)
00028
00029 do ipg = 1, npg
00030 do itg = 1, ntg
00031 do irg = 1, nrg
00032 call interpo_linear_type0(psigc,psi,irg,itg,ipg)
00033 psiinv = 1.0d0/psigc
00034 dxpsi = psi_grad2x(irg,itg,ipg)
00035 dypsi = psi_grad2y(irg,itg,ipg)
00036 dzpsi = psi_grad2z(irg,itg,ipg)
00037 dxalph = alph_grad2x(irg,itg,ipg)
00038 dyalph = alph_grad2y(irg,itg,ipg)
00039 dzalph = alph_grad2z(irg,itg,ipg)
00040
00041 sou(irg,itg,ipg) = -2.0d0*psiinv &
00042 & *(dxpsi*dxalph + dypsi*dyalph + dzpsi*dzalph)
00043 end do
00044 end do
00045 end do
00046
00047 deallocate(psi_grad2x)
00048 deallocate(psi_grad2y)
00049 deallocate(psi_grad2z)
00050 deallocate(alph_grad2x)
00051 deallocate(alph_grad2y)
00052 deallocate(alph_grad2z)
00053 end subroutine sourceterm_volume_int_bbh_2pot_test
00054