00001 subroutine calc_mass_BBH_CF
00002 use phys_constant, only : long, pi
00003 use grid_parameter, only : nrg, ntg, npg
00004 use make_array_1d
00005 use make_array_2d
00006 use def_quantities, only : admmass, komarmass, admmass_asymp, komarmass_asymp
00007 use interface_surf_source_adm_mass
00008 use interface_surf_source_komar_mass
00009
00010
00011 use interface_surf_int_grav_rg
00012 implicit none
00013 real(long) :: fac2pi , fac4pi
00014 real(long) :: surf_int_adm, surf_int_kom
00015 real(long) :: adm_em, komar_em
00016 real(long),pointer :: surf_sou_adm(:,:), surf_sou_kom(:,:)
00017 real(long),pointer :: adm_mass(:), komar_mass(:)
00018 integer :: mass_ir, ii, irg, Nmass
00019
00020 Nmass = 8
00021 call alloc_array1d(adm_mass ,1,Nmass)
00022 call alloc_array1d(komar_mass,1,Nmass)
00023 adm_mass(1:Nmass) =0.0d0
00024 komar_mass(1:Nmass)=0.0d0
00025
00026 call alloc_array2d(surf_sou_adm, 1, ntg, 1, npg)
00027 call alloc_array2d(surf_sou_kom, 1, ntg, 1, npg)
00028
00029 fac2pi = 0.5d0/pi
00030 fac4pi = 0.25d0/pi
00031 call calc_mass_ir(mass_ir)
00032 write(6,'(a10,i5)') 'mass_ir = ', mass_ir
00033
00034 do ii=1,Nmass
00035 surf_sou_adm(1:ntg,1:npg) = 0.0d0
00036 surf_sou_kom(1:ntg,1:npg) = 0.0d0
00037 surf_int_adm = 0.0d0
00038 surf_int_kom = 0.0d0
00039 adm_em = 0.0d0
00040 komar_em = 0.0d0
00041
00042 irg = mass_ir+ii-Nmass/2
00043
00044 call surf_source_adm_mass(surf_sou_adm,irg)
00045 call surf_int_grav_rg(surf_sou_adm, surf_int_adm, irg)
00046 adm_mass(ii) = -fac2pi*surf_int_adm
00047
00048 call surf_source_komar_mass(surf_sou_kom, irg)
00049 call surf_int_grav_rg(surf_sou_kom, surf_int_kom, irg)
00050 komar_mass(ii) = fac4pi*surf_int_kom
00051
00052
00053
00054 if (ii.ge.2) then
00055 adm_em = dabs(adm_mass(ii)-adm_mass(ii-1))/adm_mass(ii)
00056 komar_em = dabs(komar_mass(ii)-komar_mass(ii-1))/komar_mass(ii)
00057 end if
00058 write(6,'(a6,i5,a15,1p,e14.6,a15,1p,e14.6,a15,1p,e14.6,a15,1p,e14.6)') 'irg = ',irg, &
00059 & ' ADM mass = ',adm_mass(ii), ' Komar mass = ',komar_mass(ii), &
00060 & ' ADM error= ',adm_em, ' Komar error= ',komar_em
00061 end do
00062 admmass = adm_mass(Nmass/2)
00063 komarmass = komar_mass(Nmass/2)
00064
00065 admmass_asymp = admmass
00066 komarmass_asymp = komarmass
00067
00068 deallocate(adm_mass)
00069 deallocate(komar_mass)
00070 deallocate(surf_sou_adm)
00071 deallocate(surf_sou_kom)
00072 end subroutine calc_mass_BBH_CF
