00001 subroutine iteration_poisson_1bh_2pot_test(iter_count)
00002 use phys_constant, only : long, nnrg, nntg, nnpg
00003 use grid_parameter
00004 use coordinate_grav_r
00005 use coordinate_grav_phi
00006 use coordinate_grav_theta
00007 use weight_midpoint_grav
00008 use def_metric, only : psi, alph
00009 use make_array_2d
00010 use make_array_3d
00011 use interface_poisson_solver
00012 use interface_update_grfield
00013 use interface_error_metric_type0
00014 use interface_interpo_fl2gr
00015 use interface_sourceterm_poisson_solver_test
00016 use interface_sourceterm_exsurf_eqm_binary
00017 use interface_sourceterm_surface_int
00018
00019 use interface_sourceterm_1bh_2pot_test
00020 use interface_poisson_solver_1bh_homosol
00021
00022
00023 use interface_bh_boundary_1bh_nh_psi_test
00024 use interface_bh_boundary_1bh_nh_alph_test
00025
00026 implicit none
00027 real(long), pointer :: sou(:,:,:)
00028 real(long), pointer :: pot(:,:,:), psi_bak(:,:,:), alph_bak(:,:,:)
00029 real(long), pointer :: sou_exsurf(:,:), dsou_exsurf(:,:)
00030 real(long), pointer :: sou_bhsurf(:,:), dsou_bhsurf(:,:)
00031 real(long), pointer :: sou_outsurf(:,:), dsou_outsurf(:,:)
00032 real(long) :: error_psi, count, error_alph
00033 integer :: iter_count, flag_psi = 0, flag_alph = 0
00034
00035 call alloc_array3d(sou,0,nrg,0,ntg,0,npg)
00036 call alloc_array3d(pot,0,nrg,0,ntg,0,npg)
00037 call alloc_array3d(psi_bak,0,nrg,0,ntg,0,npg)
00038 call alloc_array3d(alph_bak,0,nrg,0,ntg,0,npg)
00039 call alloc_array2d(sou_exsurf,0,ntg,0,npg)
00040 call alloc_array2d(dsou_exsurf,0,ntg,0,npg)
00041 call alloc_array2d(sou_bhsurf,0,ntg,0,npg)
00042 call alloc_array2d(dsou_bhsurf,0,ntg,0,npg)
00043 call alloc_array2d(sou_outsurf,0,ntg,0,npg)
00044 call alloc_array2d(dsou_outsurf,0,ntg,0,npg)
00045
00046 iter_count = 0
00047 do
00048 iter_count = iter_count + 1
00049 count = dble(iter_count)
00050 conv_gra = dmin1(conv0_gra,conv_ini*count)
00051 conv_den = dmin1(conv0_den,conv_ini*count)
00052
00053 call calc_vector_x_grav(1)
00054
00055 call calc_vector_phi_grav(1)
00056
00057
00058
00059
00060
00061
00062 sou(0:nrg,0:ntg,0:npg) = 0.0d0
00063 call reset_bh_boundary('n')
00064 call sourceterm_surface_int(psi,0,sou_bhsurf,dsou_bhsurf)
00065 call sourceterm_surface_int(psi,nrg,sou_outsurf,dsou_outsurf)
00066
00067 call bh_boundary_1bh_nh_psi_test(dsou_bhsurf)
00068 call poisson_solver_1bh_homosol('nd',sou, &
00069 & sou_bhsurf,dsou_bhsurf, &
00070 & sou_outsurf,dsou_outsurf,pot)
00071 psi_bak = psi
00072 call update_grfield(pot,psi,conv_gra)
00073 call error_metric_type0(psi,psi_bak,error_psi,flag_psi,'bh')
00074 call printout_error_metric(iter_count,error_psi)
00075 call reset_bh_boundary('n')
00076
00077
00078
00079
00080
00081 sou(0:nrg,0:ntg,0:npg) = 0.0d0
00082
00083 call sourceterm_1bh_2pot_test(sou)
00084 call sourceterm_surface_int(alph,0,sou_bhsurf,dsou_bhsurf)
00085 call sourceterm_surface_int(alph,nrg,sou_outsurf,dsou_outsurf)
00086
00087 call bh_boundary_1bh_nh_alph_test(dsou_bhsurf)
00088 call poisson_solver_1bh_homosol('nd',sou, &
00089 & sou_bhsurf,dsou_bhsurf, &
00090 & sou_outsurf,dsou_outsurf,pot)
00091 alph_bak = alph
00092 call update_grfield(pot,alph,conv_gra)
00093 call error_metric_type0(alph,alph_bak,error_alph,flag_alph,'bh')
00094 call printout_error_metric(iter_count,error_alph)
00095 call reset_bh_boundary('n')
00096
00097 if ((flag_psi==0).and.(flag_alph==0)) exit
00098 if (iter_count >= iter_max) exit
00099 if (iter_count >= 50 .and. error_psi > 1.5d0) exit
00100 if (iter_count >= 50 .and. error_alph > 1.5d0) exit
00101 flag_psi = 0
00102 flag_alph = 0
00103 end do
00104
00105 deallocate(sou)
00106 deallocate(pot)
00107 deallocate(psi_bak)
00108 deallocate(alph_bak)
00109 deallocate(sou_bhsurf)
00110 deallocate(dsou_bhsurf)
00111 deallocate(sou_outsurf)
00112 deallocate(dsou_outsurf)
00113 end subroutine iteration_poisson_1bh_2pot_test
