00001 subroutine source_ang_mom_exc(sous)
00002   use phys_constant, only : long
00003   use grid_parameter, only : nrg, ntg, npg
00004   use make_array_2d
00005   use def_metric, only :  psi, tfkij
00006   use def_metric_excurve_grid, only : tfkij_grid
00007   use coordinate_grav_r, only : rg
00008   use def_vector_irg, only : hvec_irg_cm_phig, hvec_irg_cbh_xg
00009   use interface_interpo_linear_type0_2Dsurf
00010   use grid_parameter_binary_excision, only: ex_nrg
00011   implicit none
00012   real(long), external :: lagint_2nd
00013   real(long), pointer :: sous(:,:), psi_exc(:,:)
00014   real(long) :: ni, vphi_bh, vphi_cm, Aij_surf, val, psi6, work(2,2)
00015   real(long) :: x(2),y(2), v
00016   integer    :: irg, itg, ipg, ia, ib
00017 
00018   call alloc_array2d(psi_exc,0,ntg,0,npg)
00019   call calc_vector_irg(2,ex_nrg)
00020 
00021 
00022 
00023 
00024 
00025 
00026 
00027 
00028 
00029 
00030 
00031 
00032 
00033 
00034 
00035 
00036   psi_exc(0:ntg,0:npg) = psi(ex_nrg,0:ntg,0:npg)
00037   do ipg = 1, npg
00038     do itg = 1, ntg
00039       call interpo_linear_type0_2Dsurf(val,psi_exc,itg,ipg)
00040       psi6 = val**6
00041       sous(itg,ipg)=0.0d0
00042       do ib = 1, 3
00043         do ia = 1, 3
00044           work(1:2,1:2) = tfkij_grid(0,itg-1:itg,ipg-1:ipg,ia,ib)
00045           call interpo_linear1p_type0_2Dsurf(Aij_surf,work)          
00046 
00047 
00048           ni = -hvec_irg_cbh_xg(itg,ipg,ib)/rg(ex_nrg)
00049           vphi_cm = hvec_irg_cm_phig(itg,ipg,ia)
00050 
00051           sous(itg,ipg) = sous(itg,ipg) + Aij_surf*vphi_cm*ni
00052         end do
00053       end do
00054       sous(itg,ipg) = sous(itg,ipg)*psi6
00055     end do
00056   end do
00057 
00058   deallocate(psi_exc)
00059 end subroutine source_ang_mom_exc
00060