00001 subroutine iteration(iter_count)_peos
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 make_array_2d
00009 use make_array_3d
00010 use make_array_4d
00011 use def_metric
00012 use def_matter
00013 use interface_interpo_fl2gr
00014 use interface_poisson_solver
00015 use interface_sourceterm_HaC_peos
00016 use interface_sourceterm_trG_peos
00017 use interface_sourceterm_MoC_peos
00018 use interface_compute_shift
00019 use interface_compute_alps2alph
00020 use interface_hydrostatic_eq_peos
00021 use interface_calc_surface
00022 use interface_update_grfield
00023 use interface_update_matter
00024 use interface_update_parameter_peos
00025 use interface_update_surface
00026 use interface_error_matter
00027 implicit none
00028 real(long), pointer :: sou(:,:,:), pot(:,:,:), sou2(:,:,:)
00029 real(long), pointer :: potf(:,:,:), emd_bak(:,:,:)
00030 real(long), pointer :: potrs(:,:)
00031 real(long), pointer :: potx(:,:,:), poty(:,:,:), potz(:,:,:)
00032 real(long), pointer :: souvec(:,:,:,:)
00033 real(long) :: work(0:nnrg,0:nntg,0:nnpg)
00034 real(long) :: error_emd, count
00035 integer :: iter_count, flag = 0, hydro_iter = 4
00036 integer :: irf, itf, ipf, irg, itg, ipg, ihy
00037
00038 call alloc_array3d(sou,0,nrg,0,ntg,0,npg)
00039 call alloc_array3d(sou2,0,nrg,0,ntg,0,npg)
00040 call alloc_array3d(potx,0,nrg,0,ntg,0,npg)
00041 call alloc_array3d(poty,0,nrg,0,ntg,0,npg)
00042 call alloc_array3d(potz,0,nrg,0,ntg,0,npg)
00043 call alloc_array3d(pot,0,nrg,0,ntg,0,npg)
00044 call alloc_array3d(potf,0,nrf,0,ntf,0,npf)
00045 call alloc_array3d(emd_bak,0,nrf,0,ntf,0,npf)
00046 call alloc_array2d(potrs,0,ntf,0,npf)
00047 call alloc_array4d(souvec,0,nrg,0,ntg,0,npg,1,3)
00048
00049 iter_count = 0
00050 do
00051 iter_count = iter_count + 1
00052 count = dble(iter_count)
00053 conv_gra = dmin1(conv0_gra,conv_ini*count)
00054 conv_den = dmin1(conv0_den,conv_ini*count)
00055
00056 emdg = 0.0d0
00057 call interpo_fl2gr(emd, emdg)
00058 call calc_vector_x_grav(1)
00059 call calc_vector_x_matter(1)
00060 call calc_vector_phi_grav(1)
00061 call calc_vector_phi_matter(1)
00062 call excurve
00063
00064 call sourceterm_HaC_peos(sou)
00065 call poisson_solver(sou,pot)
00066 pot = pot + 1.0d0
00067 call update_grfield(pot,psi,conv_gra)
00068 call update_parameter_peos(conv_den)
00069
00070
00071 call sourceterm_trG_peos(sou)
00072 call poisson_solver(sou,pot)
00073 pot = pot + 1.0d0
00074 call compute_alps2alph(pot,psi)
00075 call update_grfield(pot,alph,conv_gra)
00076 call update_parameter_peos(conv_den)
00077
00078 call sourceterm_MoC_peos(souvec,sou)
00079 work(1:nrg, 1:ntg, 1:npg) = souvec(1:nrg, 1:ntg, 1:npg, 1)
00080 sou2(1:nrg, 1:ntg, 1:npg) = work(1:nrg, 1:ntg, 1:npg)
00081 call poisson_solver(sou2,potx)
00082 work(1:nrg, 1:ntg, 1:npg) = souvec(1:nrg, 1:ntg, 1:npg, 2)
00083 sou2(1:nrg, 1:ntg, 1:npg) = work(1:nrg, 1:ntg, 1:npg)
00084 call poisson_solver(sou2,poty)
00085 work(1:nrg, 1:ntg, 1:npg) = souvec(1:nrg, 1:ntg, 1:npg, 3)
00086 sou2(1:nrg, 1:ntg, 1:npg) = work(1:nrg, 1:ntg, 1:npg)
00087 call poisson_solver(sou2,potz)
00088 call poisson_solver(sou,pot)
00089 work( 0:nrg, 0:ntg, 0:npg) = potx(0:nrg, 0:ntg, 0:npg)
00090 souvec(0:nrg, 0:ntg, 0:npg, 1) = work(0:nrg, 0:ntg, 0:npg)
00091 work( 0:nrg, 0:ntg, 0:npg) = poty(0:nrg, 0:ntg, 0:npg)
00092 souvec(0:nrg, 0:ntg, 0:npg, 2) = work(0:nrg, 0:ntg, 0:npg)
00093 work( 0:nrg, 0:ntg, 0:npg) = potz(0:nrg, 0:ntg, 0:npg)
00094 souvec(0:nrg, 0:ntg, 0:npg, 3) = work(0:nrg, 0:ntg, 0:npg)
00095 call compute_shift(potx, poty, potz, souvec, pot)
00096 call update_grfield(potx,bvxd,conv_gra)
00097 call update_grfield(poty,bvyd,conv_gra)
00098 call update_grfield(potz,bvzd,conv_gra)
00099 call update_parameter_peos(conv_den)
00100
00101 do ihy = 1, hydro_iter
00102 call hydrostatic_eq_peos(potf)
00103 call calc_surface(potrs, potf)
00104 emd_bak = emd
00105 call update_matter(potf,emd,conv_den)
00106 call update_surface(potrs,rs,conv_den)
00107 call reset_fluid
00108 call update_parameter_peos(conv_den)
00109 call calc_vector_x_matter(1)
00110 call calc_vector_phi_matter(1)
00111
00112 if (ihy.eq.1) then
00113 call error_matter(emd,emd_bak,error_emd,flag)
00114 call printout_error(iter_count,error_emd)
00115 end if
00116 end do
00117
00118 if (flag == 0) exit
00119 if (iter_count >= iter_max) exit
00120 flag = 0
00121
00122 end do
00123
00124 deallocate(sou)
00125 deallocate(sou2)
00126 deallocate(potx)
00127 deallocate(poty)
00128 deallocate(potz)
00129 deallocate(pot)
00130 deallocate(potf)
00131 deallocate(emd_bak)
00132 deallocate(potrs)
00133 deallocate(souvec)
00134 end subroutine iteration_peos
