00001 subroutine printout_physq_peos(iseq,flag_restmass)
00002 use grid_parameter, only : ntfeq, npfxzp, npfyzp, ntfpolp, sw_mass_iter
00003 use def_matter, only : emd, rs
00004 use def_matter_parameter, only : ome, radi, pinx
00005 use def_quantities
00006 implicit none
00007 real(8) :: fixeddlm, fixedvir
00008 integer, intent(in) :: iseq, flag_restmass
00009
00010 if (iseq.eq.1) then
00011 open(100,file='rnsphyseq.dat',status='unknown')
00012 end if
00013 if (flag_restmass.le.1) write(100,*) '#### Solution did not converge #### '
00014 write(100,*) '== Sequence number == ', iseq
00015 write(100,*) '## Coordinate and Proper Radii in G = c = Msol = 1 unit ##'
00016 write(100,'(a22,1p,2e23.15)') ' NS radius along x = ', &
00017 & coord_radius_x, proper_radius_x
00018 write(100,'(a22,1p,2e23.15)') ' NS radius along y = ', &
00019 & coord_radius_y, proper_radius_y
00020 write(100,'(a22,1p,2e23.15)') ' NS radius along z = ', &
00021 & coord_radius_z, proper_radius_z
00022 write(100,'(a22,1p,2e23.15)') ' Axis ratio y/x = ', &
00023 & coord_radius_y/coord_radius_x, proper_radius_y/proper_radius_x
00024 write(100,'(a22,1p,2e23.15)') ' Axis ratio z/x = ', &
00025 & coord_radius_z/coord_radius_x, proper_radius_z/proper_radius_x
00026 write(100,*) '## Coordinate and Proper Radii in [km] ##'
00027 write(100,'(a22,1p,2e23.15)') ' NS radius along x = ', &
00028 & coord_radius_x_km, proper_radius_x_km
00029 write(100,'(a22,1p,2e23.15)') ' NS radius along y = ', &
00030 & coord_radius_y_km, proper_radius_y_km
00031 write(100,'(a22,1p,2e23.15)') ' NS radius along z = ', &
00032 & coord_radius_z_km, proper_radius_z_km
00033
00034 write(100,*) '## Radii in Req = 1 unit ##'
00035 write(100,'(a22,1p,2e23.15)') ' NS radius along x = ', rs(ntfeq,npfxzp)
00036 write(100,'(a22,1p,2e23.15)') ' NS radius along y = ', rs(ntfeq,npfyzp)
00037 write(100,'(a22,1p,2e23.15)') ' NS radius along z = ', rs(ntfpolp,npfxzp)
00038 write(100,*) '## Eccentricity in coordinate and Proper Radii ##'
00039 write(100,'(a22,1p,2e23.15)') ' sqrt(1-(Rz/Rx)^2) = ', &
00040 & sqrt(1.0d0-(coord_radius_z/coord_radius_x)**2), &
00041 & sqrt(1.0d0-(proper_radius_z/proper_radius_x)**2)
00042
00043 write(100,*) '## M = spherical gravitational mass (G = c = Msol = 1 unit) ##'
00044 write(100,'(a22,1p,2e23.15)') ' Omega M and Omega = ', &
00045 & ome/radi*gravmass_sph, ome/radi
00046 write(100,'(a22,1p,2e23.15)') ' M_ADM = ', admmass,admmass_asymp
00047 write(100,'(a22,1p,2e23.15)') ' M_K Komar mass = ', komarmass,komarmass_asymp
00048 write(100,'(a22,1p,2e23.15)') ' Rest mass M_0 1star= ', restmass
00049 write(100,'(a22,1p,2e23.15)') ' Proper mass M_p = ', propermass
00050 write(100,'(a22,1p,2e23.15)') ' Angular momentum J = ', angmom, angmom_asymp
00051 write(100,'(a22,1p,2e23.15)') ' Spherical rest mass= ', restmass_sph
00052 write(100,'(a22,1p,2e23.15)') ' Spherical grav mass= ', gravmass_sph
00053 write(100,'(a22,1p,2e23.15)') ' Spherical M/R = ', MoverR_sph
00054 write(100,'(a22,1p,2e23.15)') ' Schwarzschildradius= ', schwarz_radi_sph
00055 write(100,'(a22,1p,2e23.15)') ' Schwarz radi in[km]= ', schwarz_radi_sph_km
00056
00057 write(100,'(a22,1p,2e23.15)') ' E/M = (M_ADM-M)/M = ', &
00058 & admmass/gravmass_sph - 1.0d0
00059 write(100,'(a22,1p,2e23.15)') ' J/M^2 and J/M_ADM^2= ', &
00060 & angmom/gravmass_sph**2.0d0, angmom/admmass**2.0d0
00061
00062 write(100,'(a22,1p,2e23.15)') ' T kinetic energy = ', T_kinene
00063 write(100,'(a22,1p,2e23.15)') ' W grav energy = ', W_gravene
00064 write(100,'(a22,1p,2e23.15)') ' T/|W| = ', ToverW
00065 write(100,'(a22,1p,2e23.15)') ' Moment of inertia = ', I_inertia
00066
00067 write(100,*) '## Virial and Rest mass accuracy ##'
00068 fixedvir = (admmass - komarmass)/admmass
00069 fixeddlm = (restmass - restmass_sph)/restmass_sph
00070 write(100,'(a22,1p,2e23.15)') ' 1 - M_K/M_ADM = ', fixedvir
00071 write(100,'(a22,1p,2e23.15)') ' M_0/M_0sph - 1 = ', fixeddlm
00072
00073 write(100,*) '## The density, pressure and p/rho at the NS center ##'
00074 write(100,*) '## in G = c = M = 1 unit and cgs unit ##'
00075 write(100,'(a22,1p,2e23.15)') ' rho_c = ', rho_c, rho_c_cgs
00076 write(100,'(a22,1p,2e23.15)') ' pre_c = ', pre_c, pre_c_cgs
00077 write(100,'(a22,1p,2e23.15)') ' epsilon_c = ', epsi_c, epsi_c_cgs
00078 write(100,'(a22,1p,2e23.15)') ' (p/rho)_c = ', q_c, q_c_cgs
00079
00080 write(100,*) '## Red and blue shift ##'
00081 write(100,'(a22,1p,2e23.15)') ' surface on x axis = ', &
00082 & zrb_xp_plus, zrb_xp_minus
00083 write(100,'(a22,1p,2e23.15)') ' surface on y axis = ', &
00084 & zrb_yp_plus, zrb_yp_minus
00085 write(100,'(a22,1p,2e23.15)') ' surface on z axis = ', &
00086 & zrb_zp_plus, zrb_zp_minus
00087
00088 write(100,*) '## Ratio of enthalpy gradient at x/z and y/z ##'
00089 write(100,'(a22,1p,2e23.15)') ' dhdx/dhdz,dhdy/dhdz= ', &
00090 & dhdr_x/dhdr_z, dhdr_y/dhdr_z
00091
00092 end subroutine printout_physq_peos
