00001 subroutine printout_physq_WL_MHD_v1(iseq,flag_restmass)
00002 use phys_constant, only : long, pi
00003 use grid_parameter, only : ntfeq, npfxzp, npfyzp, ntfpolp, sw_mass_iter
00004 use def_matter, only : emd, rs
00005 use def_matter_parameter, only : ome, radi, pinx
00006 use def_quantities
00007 implicit none
00008 real(long) :: fixeddlm, fixedvir, fixedvir_asymp
00009 real(long) :: omega
00010 integer, intent(in) :: iseq, flag_restmass
00011
00012 if (iseq.eq.1) then
00013 open(100,file='rnsphyseq.dat',status='unknown')
00014 else
00015 open(100,file='rnsphyseq.dat',status='old',position='append')
00016 end if
00017
00018 if (flag_restmass.le.1) write(100,*) '#### Solution did not converge #### '
00019 write(100,*) '== Sequence number == ', iseq
00020 write(100,*) '## Coordinate and Proper Radii in G = c = Msol = 1 unit ##'
00021 write(100,'(a22,1p,2e23.15)') ' NS radius along x = ', &
00022 & coord_radius_x, proper_radius_x
00023 write(100,'(a22,1p,2e23.15)') ' NS radius along y = ', &
00024 & coord_radius_y, proper_radius_y
00025 write(100,'(a22,1p,2e23.15)') ' NS radius along z = ', &
00026 & coord_radius_z, proper_radius_z
00027 write(100,'(a22,1p,2e23.15)') ' Axis ratio y/x = ', &
00028 & coord_radius_y/coord_radius_x, proper_radius_y/proper_radius_x
00029 write(100,'(a22,1p,2e23.15)') ' Axis ratio z/x = ', &
00030 & coord_radius_z/coord_radius_x, proper_radius_z/proper_radius_x
00031 write(100,*) '## Coordinate and Proper Radii in [km] ##'
00032 write(100,'(a22,1p,2e23.15)') ' NS radius along x = ', &
00033 & coord_radius_x_km, proper_radius_x_km
00034 write(100,'(a22,1p,2e23.15)') ' NS radius along y = ', &
00035 & coord_radius_y_km, proper_radius_y_km
00036 write(100,'(a22,1p,2e23.15)') ' NS radius along z = ', &
00037 & coord_radius_z_km, proper_radius_z_km
00038
00039 write(100,*) '## Radii in Req = 1 unit ##'
00040 write(100,'(a22,1p,2e23.15)') ' NS radius along x = ', rs(ntfeq,npfxzp)
00041 write(100,'(a22,1p,2e23.15)') ' NS radius along y = ', rs(ntfeq,npfyzp)
00042 write(100,'(a22,1p,2e23.15)') ' NS radius along z = ', rs(ntfpolp,npfxzp)
00043 write(100,*) '## Eccentricity in coordinate and Proper Radii ##'
00044 write(100,'(a22,1p,2e23.15)') ' sqrt(1-(Rz/Rx)^2) = ', &
00045 & sqrt(1.0d0-(coord_radius_z/coord_radius_x)**2), &
00046 & sqrt(1.0d0-(proper_radius_z/proper_radius_x)**2)
00047
00048 write(100,*) '## M = spherical gravitational mass (G = c = Msol = 1 unit) ##'
00049 write(100,'(a22,1p,2e23.15)') ' Omega M and Omega = ', &
00050 & ome/radi*gravmass_sph, ome/radi
00051 write(100,'(a22,1p,2e23.15)') ' M_ADM = ', admmass,admmass_asymp
00052 write(100,'(a22,1p,2e23.15)') ' M_K Komar mass = ', komarmass, &
00053 & komarmass_asymp
00054 write(100,'(a22,1p,2e23.15)') ' Rest mass M_0 1star= ', restmass
00055 write(100,'(a22,1p,2e23.15)') ' Proper mass M_p = ', propermass
00056 write(100,'(a22,1p,2e23.15)') ' Angular momentum J = ', angmom, angmom_asymp
00057 write(100,'(a22,1p,2e23.15)') ' Electric charge Q = ', charge, charge_asymp
00058
00059 write(100,'(a22,1p,2e23.15)') ' Spherical rest mass= ', restmass_sph
00060 write(100,'(a22,1p,2e23.15)') ' Spherical grav mass= ', gravmass_sph
00061 write(100,'(a22,1p,2e23.15)') ' Spherical M/R = ', MoverR_sph
00062 write(100,'(a22,1p,2e23.15)') ' Schwarzschildradius= ', schwarz_radi_sph
00063 write(100,'(a22,1p,2e23.15)') ' Schwarz radi in[km]= ', schwarz_radi_sph_km
00064
00065 write(100,'(a22,1p,2e23.15)') ' E/M = (M_ADM-M)/M = ', &
00066 & admmass/gravmass_sph - 1.0d0
00067 write(100,'(a22,1p,2e23.15)') ' J/M^2 and J/M_ADM^2= ', &
00068 & angmom/gravmass_sph**2.0d0, angmom/admmass**2.0d0
00069
00070 write(100,'(a22,1p,2e23.15)') ' T kinetic energy = ', T_kinene, &
00071 & T_kinene_omeJ
00072 write(100,'(a22,1p,2e23.15)') ' P internal energy = ', P_intene
00073 write(100,'(a22,1p,2e23.15)') ' M em field energy = ', M_emfene
00074 write(100,'(a22,1p,2e23.15)') ' W grav energy = ', W_gravene, &
00075 & W_gravene_omeJ
00076 write(100,'(a22,1p,2e23.15)') ' T/|W| = ', ToverW, &
00077 & ToverW_omeJ
00078 write(100,'(a22,1p,2e23.15)') ' P/|W| = ', PoverW
00079 write(100,'(a22,1p,2e23.15)') ' M/|W| = ', MoverW
00080 write(100,'(a22,1p,2e23.15)') ' Virial relation = ', Virial
00081 write(100,'(a22,1p,2e23.15)') ' MtorB field energy = ', M_torBene
00082 write(100,'(a22,1p,2e23.15)') ' MpolB field energy = ', M_polBene
00083 write(100,'(a22,1p,2e23.15)') ' MeleE field energy = ', M_eleEene
00084 write(100,'(a22,1p,2e23.15)') ' MtorB/|W| = ', MtorBoverW
00085 write(100,'(a22,1p,2e23.15)') ' MpolB/|W| = ', MpolBoverW
00086 write(100,'(a22,1p,2e23.15)') ' MeleE/|W| = ', MeleEoverW
00087 write(100,'(a22,1p,2e23.15)') ' Moment of inertia = ', I_inertia
00088
00089 write(100,*) '## Virial and Rest mass accuracy ##'
00090 fixedvir = (admmass - komarmass)/admmass
00091 fixedvir_asymp = (admmass_asymp - komarmass_asymp)/admmass_asymp
00092 fixeddlm = (restmass - restmass_sph)/restmass_sph
00093 write(100,'(a22,1p,2e23.15)') ' 1 - M_K/M_ADM = ', fixedvir, &
00094 & fixedvir_asymp
00095 write(100,'(a22,1p,2e23.15)') ' M_0/M_0sph - 1 = ', fixeddlm
00096
00097 write(100,*) '## The density, pressure and p/rho at the NS center ##'
00098 write(100,*) '## in G = c = M = 1 unit and cgs unit ##'
00099 write(100,'(a22,1p,2e23.15)') ' rho_c = ', rho_c, rho_c_cgs
00100 write(100,'(a22,1p,2e23.15)') ' pre_c = ', pre_c, pre_c_cgs
00101 write(100,'(a22,1p,2e23.15)') ' epsilon_c = ', epsi_c, epsi_c_cgs
00102 write(100,'(a22,1p,2e23.15)') ' (p/rho)_c = ', q_c, q_c_cgs
00103
00104 write(100,*) '## Red and blue shift ##'
00105 write(100,'(a22,1p,2e23.15)') ' surface on x axis = ', &
00106 & zrb_xp_plus, zrb_xp_minus
00107 write(100,'(a22,1p,2e23.15)') ' surface on y axis = ', &
00108 & zrb_yp_plus, zrb_yp_minus
00109 write(100,'(a22,1p,2e23.15)') ' surface on z axis = ', &
00110 & zrb_zp_plus, zrb_zp_minus
00111
00112 write(100,*) '## Ratio of enthalpy gradient at x/z and y/z ##'
00113 write(100,'(a22,1p,2e23.15)') ' dhdx/dhdz,dhdy/dhdz= ', &
00114 & dhdr_x/dhdr_z, dhdr_y/dhdr_z
00115
00116 write(100,*) '## Angular velocity and frequency in cgs unit ##'
00117 write(100,'(a22,1p,2e23.15)') ' Omega[rad/s], f[Hz]= ', ome_cgs, ome_cgs/2.0d0/pi
00118
00119 write(100,*) '## Quasilocal spin estimate on the surface and at nrf+1 ##'
00120 write(100,'(a22,1p,2e23.15)') ' Quasi local spin S = ', qua_loc_spin_surf, qua_loc_spin
00121 write(100,'(a22,1p,2e23.15)') ' S/M^2 and S/M_ADM^2= ', &
00122 & qua_loc_spin/gravmass_sph**2.0d0, qua_loc_spin/admmass**2.0d0
00123
00124 write(100,*) '## Circulation line and surface integrals at 3 planes ##'
00125 write(100,'(a22,1p,2e23.15)') ' Circul at xy plane = ', circ_line_xy, circ_surf_xy
00126 write(100,'(a22,1p,2e23.15)') ' Circul at yz plane = ', circ_line_yz, circ_surf_yz
00127 write(100,'(a22,1p,2e23.15)') ' Circul at zx plane = ', circ_line_zx, circ_surf_zx
00128
00129 close(100)
00130
00131 end subroutine printout_physq_WL_MHD_v1
