Skip to content

Commit

Permalink
Renaming masse_mol_gaz to mu_mH
Browse files Browse the repository at this point in the history
  • Loading branch information
cpinte committed Dec 11, 2024
1 parent a0892a2 commit 9d2a986
Show file tree
Hide file tree
Showing 17 changed files with 40 additions and 41 deletions.
4 changes: 2 additions & 2 deletions src/ML_prodimo.f90
Original file line number Diff line number Diff line change
Expand Up @@ -195,15 +195,15 @@ subroutine xgb_compute_features()
feature_Tgas(2,:) = z_grid(:)
endif
feature_Tgas(3,:) = Tdust(:)
feature_Tgas(4,:) = densite_gaz(:) * masse_mol_gaz / m3_to_cm3 ! g.cm^3
feature_Tgas(4,:) = densite_gaz(:) * mu_mH / m3_to_cm3 ! g.cm^3
feature_Tgas(5:43,:) = J_ML(:,:)
feature_Tgas(44:47,:) = N_grains(:,:)
do i=1,n_directions
feature_Tgas(48+i-1,:) = CD(:,i) ! CD(n_cells, n_directions)
enddo
else if (n_features == 45) then
feature_Tgas(1,:) = Tdust(:)
feature_Tgas(2,:) = densite_gaz(:) * masse_mol_gaz / m3_to_cm3 ! g.cm^3
feature_Tgas(2,:) = densite_gaz(:) * mu_mH / m3_to_cm3 ! g.cm^3
feature_Tgas(3:41,:) = J_ML(:,:)
feature_Tgas(42:45,:) = N_grains(:,:)
endif
Expand Down
6 changes: 3 additions & 3 deletions src/SPH2mcfost.f90
Original file line number Diff line number Diff line change
Expand Up @@ -304,8 +304,8 @@ subroutine SPH_to_Voronoi(n_SPH, ndusttypes, particle_id, x,y,z,h, vx,vy,vz, T_g
call allocate_densities(n_cells_max = n_SPH + n_etoiles) ! we allocate all the SPH particule for libmcfost
! Tableau de densite et masse de gaz
!do icell=1,n_cells
! densite_gaz(icell) = rho(icell) / masse_mol_gaz * m3_to_cm3 ! rho is in g/cm^3 --> part.m^3
! masse_gaz(icell) = densite_gaz(icell) * masse_mol_gaz * volume(icell)
! densite_gaz(icell) = rho(icell) / mu_mH * m3_to_cm3 ! rho is in g/cm^3 --> part.m^3
! masse_gaz(icell) = densite_gaz(icell) * mu_mH * volume(icell)
!enddo
!masse_gaz(:) = masse_gaz(:) * AU3_to_cm3

Expand All @@ -314,7 +314,7 @@ subroutine SPH_to_Voronoi(n_SPH, ndusttypes, particle_id, x,y,z,h, vx,vy,vz, T_g
iSPH = Voronoi(icell)%id
if (iSPH > 0) then
masse_gaz(icell) = massgas(iSPH) * Msun_to_g ! en g
densite_gaz(icell) = masse_gaz(icell) / (masse_mol_gaz * volume(icell) * AU3_to_m3)
densite_gaz(icell) = masse_gaz(icell) / (mu_mH * volume(icell) * AU3_to_m3)
else ! star
masse_gaz(icell) = 0.
densite_gaz(icell) = 0.
Expand Down
6 changes: 3 additions & 3 deletions src/constants.f90
Original file line number Diff line number Diff line change
Expand Up @@ -31,9 +31,9 @@ module constantes

real(kind=dp), parameter :: Na = 6.022140857e23_dp ! Nombre d'Avogadro CODATA 2014
real(kind=dp), parameter :: amu = 1.0_dp/Na ! atomic mass unit [g]
real(kind=dp), parameter :: masseH = 1.007825032231_dp * amu ! H atom mass [g] CODATA 2014
real, parameter :: mu = 2.3 ! [g] 2.3 following Walker 2004
real, parameter :: masse_mol_gaz = mu * masseH
real(kind=dp), parameter :: mH = 1.007825032231_dp * amu ! H atom mass [g] CODATA 2014
real, parameter :: mu = 2.3 ! [mH] 2.3 following Walker 2004
real, parameter :: mu_mH = mu * mH ! mean molecular weight [g]
real, parameter :: T_Cmb = 2.7260 ! K

real(kind=dp), parameter :: e_ion_hmin = 0.754 * electron_charge !Ionization energy (affinity) Hmin in [J]
Expand Down
18 changes: 9 additions & 9 deletions src/density.f90
Original file line number Diff line number Diff line change
Expand Up @@ -88,7 +88,7 @@ subroutine define_gas_density()

! Facteur multiplicatif pour passer en masse de gaz
! puis en nH2/AU**3 puis en nH2/m**3
cst_gaz(izone) = C * dz%diskmass * dz%gas_to_dust / masse_mol_gaz * Msun_to_g / AU3_to_m3
cst_gaz(izone) = C * dz%diskmass * dz%gas_to_dust / mu_mH * Msun_to_g / AU3_to_m3
enddo

do izone=1, n_zones
Expand Down Expand Up @@ -282,7 +282,7 @@ subroutine define_gas_density()
! Calcul de la masse de gaz de la zone
mass = 0.
do icell = 1, n_cells
mass = mass + densite_gaz_tmp(icell) * masse_mol_gaz * volume(icell)
mass = mass + densite_gaz_tmp(icell) * mu_mH * volume(icell)
enddo
mass = mass * AU3_to_m3 * g_to_Msun

Expand Down Expand Up @@ -318,7 +318,7 @@ subroutine define_gas_density()

! Tableau de masse de gaz
do icell=1,n_cells
masse_gaz(icell) = densite_gaz(icell) * masse_mol_gaz * volume(icell) * AU3_to_m3
masse_gaz(icell) = densite_gaz(icell) * mu_mH * volume(icell) * AU3_to_m3
enddo
write(*,*) 'Total gas mass in model:', real(sum(masse_gaz) * g_to_Msun),' Msun'

Expand Down Expand Up @@ -442,7 +442,7 @@ subroutine define_dust_density()

do i=1, n_rad

!write(*,*) " ", rcyl, rho0*masse_mol_gaz*cm_to_m**2, dust_pop(pop)%rho1g_avg
!write(*,*) " ", rcyl, rho0*mu_mH*cm_to_m**2, dust_pop(pop)%rho1g_avg
!write(*,*) "s_opt", rcyl, s_opt/1000.

bz : do j=j_start,nz
Expand Down Expand Up @@ -678,11 +678,11 @@ subroutine define_dust_density()
do k=1, n_az
rho0 = densite_gaz(cell_map(i,1,k)) ! pour dependance en R : pb en coord sperique
!s_opt = rho_g * cs / (rho * Omega) ! cs = H * Omega ! on doit trouver 1mm vers 50AU
!omega_tau= dust_pop(ipop)%rho1g_avg*(r_grain(l)*mum_to_cm) / (rho * masse_mol_gaz/m_to_cm**3 * H*AU_to_cm)
!omega_tau= dust_pop(ipop)%rho1g_avg*(r_grain(l)*mum_to_cm) / (rho * mu_mH/m_to_cm**3 * H*AU_to_cm)
icell = cell_map(i,1,k)
rcyl = r_grid(icell)
H = dz%sclht * (rcyl/dz%rref)**dz%exp_beta
s_opt = (rho0*masse_mol_gaz*cm_to_m**3 /dust_pop(pop)%rho1g_avg) * H * AU_to_m * m_to_mum
s_opt = (rho0*mu_mH*cm_to_m**3 /dust_pop(pop)%rho1g_avg) * H * AU_to_m * m_to_mum

!write(*,*) "r=", rcyl, "a_migration =", s_opt

Expand Down Expand Up @@ -1551,7 +1551,7 @@ subroutine read_density_file()
! Calcul de la masse de gaz de la zone
mass = 0.
do icell=1,n_cells
mass = mass + densite_gaz(icell) * masse_mol_gaz * volume(icell)
mass = mass + densite_gaz(icell) * mu_mH * volume(icell)
enddo !icell
mass = mass * AU3_to_m3 * g_to_Msun

Expand All @@ -1569,7 +1569,7 @@ subroutine read_density_file()

! Tableau de masse de gaz
do icell=1,n_cells
masse_gaz(icell) = densite_gaz(icell) * masse_mol_gaz * volume(icell) * AU3_to_m3
masse_gaz(icell) = densite_gaz(icell) * mu_mH * volume(icell) * AU3_to_m3
enddo ! icell
write(*,*) 'Total gas mass in model:', real(sum(masse_gaz) * g_to_Msun),' Msun'

Expand Down Expand Up @@ -2010,7 +2010,7 @@ real(kind=dp) function omega_tau(rho,H,l)
ipop = grain(l)%pop
!write(*,*) ipop, dust_pop(ipop)%rho1g_avg, rho
if (rho > tiny_dp) then
omega_tau = dust_pop(ipop)%rho1g_avg*(r_grain(l)*mum_to_cm) / (rho * masse_mol_gaz/m_to_cm**3 * H*AU_to_cm)
omega_tau = dust_pop(ipop)%rho1g_avg*(r_grain(l)*mum_to_cm) / (rho * mu_mH/m_to_cm**3 * H*AU_to_cm)
else
omega_tau = huge_dp
endif
Expand Down
2 changes: 1 addition & 1 deletion src/disk_physics.f90
Original file line number Diff line number Diff line change
Expand Up @@ -191,7 +191,7 @@ subroutine equilibre_hydrostatique()
real, parameter :: gas_dust = 100

M_etoiles = sum(etoile(:)%M) * Msun_to_kg
M_mol = masse_mol_gaz * g_to_kg
M_mol = mu_mH * g_to_kg

cst = Ggrav * M_etoiles * M_mol / (kb * AU_to_m**2)

Expand Down
2 changes: 1 addition & 1 deletion src/for_astrochem.f90
Original file line number Diff line number Diff line change
Expand Up @@ -7,7 +7,7 @@
G(:) = compute_UV_field() ! Habing

! unit conversions for astrochem :
gas_density(:) = densite_gaz(1:n_cells) * masse_mol_gaz / m3_to_cm3 ! nH2/m**3 --> g/cm**3
gas_density(:) = densite_gaz(1:n_cells) * mu_mH / m3_to_cm3 ! nH2/m**3 --> g/cm**3

do icell=1,n_cells
dust_mass_density(icell) = sum(densite_pouss(:,icell) * M_grain(:)) ! M_grain en g
Expand Down
2 changes: 1 addition & 1 deletion src/io_prodimo.f90
Original file line number Diff line number Diff line change
Expand Up @@ -789,7 +789,7 @@ subroutine mcfost2ProDiMo()
do ri=1, n_rad
do zj=1,nz
icell = cell_map(ri,zj,1)
dens(ri,zj) = densite_gaz(icell) * masse_mol_gaz / m3_to_cm3 ! g.cm^-3
dens(ri,zj) = densite_gaz(icell) * mu_mH / m3_to_cm3 ! g.cm^-3
enddo
enddo

Expand Down
2 changes: 1 addition & 1 deletion src/mcfost2phantom.f90
Original file line number Diff line number Diff line change
Expand Up @@ -569,7 +569,7 @@ subroutine diffusion_opacity(temp,icell,kappa_diffusion)
somme2 = somme2 + B*delta_wl
enddo
kappa_diffusion = somme2/somme &
*cm_to_AU/(densite_gaz(icell)*masse_mol_gaz*(cm_to_m)**3) ! cm^2/g
*cm_to_AU/(densite_gaz(icell)*mu_mH*(cm_to_m)**3) ! cm^2/g
! check : somme2/somme * cm_to_AU /(masse_gaz(icell)/(volume(icell)*AU_to_cm**3))
else
kappa_diffusion = 0.
Expand Down
2 changes: 1 addition & 1 deletion src/mol_transfer.f90
Original file line number Diff line number Diff line change
Expand Up @@ -929,7 +929,7 @@ subroutine init_dust_mol(imol)
! AU_to_cm**2 car on veut kappa_abs_LTE en AU-1
write(*,*) "TODO : the water benchmark 3 needs to be updated for cell pointer in opacity table"
do icell=1,n_cells
kappa_abs_LTE(icell,iTrans) = kap * (densite_gaz(icell) * cm_to_m**3) * masse_mol_gaz / &
kappa_abs_LTE(icell,iTrans) = kap * (densite_gaz(icell) * cm_to_m**3) * mu_mH / &
gas_dust / cm_to_AU
enddo

Expand Down
9 changes: 4 additions & 5 deletions src/molecular_emission.f90
Original file line number Diff line number Diff line change
Expand Up @@ -33,8 +33,7 @@ module molecular_emission
real :: correct_Tgas
logical :: lcorrect_Tgas

real :: nH2, masse_mol
! masse_mol_gaz sert uniquement pour convertir masse disque en desnite de particule
real :: nH2, masse_mol ! masse_mol is the mass of the current molecule
real(kind=dp), dimension(:,:), allocatable :: kappa_mol_o_freq, kappa_mol_o_freq2 ! n_cells, nTrans
real(kind=dp), dimension(:,:), allocatable :: emissivite_mol_o_freq, emissivite_mol_o_freq2 ! n_cells, nTrans
real, dimension(:,:), allocatable :: tab_nLevel, tab_nLevel2 ! n_cells, nLevels
Expand Down Expand Up @@ -158,7 +157,7 @@ subroutine init_Doppler_profiles(imol)
! WARNING : c'est pas un sigma mais un delta, cf Cours de Boisse p47
! Consistent avec benchmark
if (trim(v_turb_unit) == "cs") then
v_turb(icell) = sqrt((kb*Tcin(icell) / (masse_mol_gaz * g_to_kg))) * v_turb(icell)
v_turb(icell) = sqrt((kb*Tcin(icell) / (mu_mH * g_to_kg))) * v_turb(icell)
endif
sigma2 = 2.0_dp * (kb*Tcin(icell) / (masse_mol * g_to_kg)) + v_turb(icell)**2
v_line(icell) = sqrt(sigma2)
Expand Down Expand Up @@ -421,8 +420,8 @@ subroutine equilibre_LTE_mol(imol)
!$omp end do
!$omp end parallel

! write(*,*) real( (sum(masse) * g_to_kg * gas_dust / masse_mol_gaz) / (4.*pi/3. * (rout * AU_to_cm)**3 ) )
! write(*,*) (sum(masse) * g_to_kg * gas_dust / masse_mol_gaz) * abundance * fAul(1) * hp * transfreq(1)
! write(*,*) real( (sum(masse) * g_to_kg * gas_dust / mu_mH) / (4.*pi/3. * (rout * AU_to_cm)**3 ) )
! write(*,*) (sum(masse) * g_to_kg * gas_dust / mu_mH) * abundance * fAul(1) * hp * transfreq(1)

return

Expand Down
2 changes: 1 addition & 1 deletion src/optical_depth.f90
Original file line number Diff line number Diff line change
Expand Up @@ -341,7 +341,7 @@ subroutine compute_column(type, column, lambda)
real(kind=dp) :: x0,y0,z0, x1,y1,z1, norme, l, u,v,w, l_contrib, l_void_before, CD_units, factor, sum

if (type==1) then
CD_units = AU_to_m * masse_mol_gaz / (m_to_cm)**2 ! g/cm^-2 and AU_to_m factor as l_contrib is in AU
CD_units = AU_to_m * mu_mH / (m_to_cm)**2 ! g/cm^-2 and AU_to_m factor as l_contrib is in AU
else if (type==3) then
CD_units = AU_to_m / (m_to_cm)**2 ! particle/cm^-2 and AU_to_m factor as l_contrib is in AU
endif
Expand Down
2 changes: 1 addition & 1 deletion src/output.f90
Original file line number Diff line number Diff line change
Expand Up @@ -1680,7 +1680,7 @@ subroutine write_disk_struct(lparticle_density,lcolumn_density,lvelocity)

dens = 0.0

dens(:) = densite_gaz(1:n_cells) * masse_mol_gaz / m3_to_cm3 ! nH2/m**3 --> g/cm**3
dens(:) = densite_gaz(1:n_cells) * mu_mH / m3_to_cm3 ! nH2/m**3 --> g/cm**3

! le e signifie real*4
call ftppre(unit,group,fpixel,nelements,dens,status)
Expand Down
4 changes: 2 additions & 2 deletions src/read_athena++.f90
Original file line number Diff line number Diff line change
Expand Up @@ -333,7 +333,7 @@ subroutine read_athena_model()
! Calcul de la masse de gaz de la zone
mass = 0.
do icell=1,n_cells
mass = mass + densite_gaz(icell) * masse_mol_gaz * volume(icell)
mass = mass + densite_gaz(icell) * mu_mH * volume(icell)
enddo !icell
mass = mass * AU3_to_m3 * g_to_Msun

Expand All @@ -344,7 +344,7 @@ subroutine read_athena_model()
! Somme sur les zones pour densite finale
do icell=1,n_cells
densite_gaz(icell) = densite_gaz(icell) * facteur
masse_gaz(icell) = densite_gaz(icell) * masse_mol_gaz * volume(icell) * AU3_to_m3
masse_gaz(icell) = densite_gaz(icell) * mu_mH * volume(icell) * AU3_to_m3
enddo ! icell
else
call error('Gas mass is 0')
Expand Down
4 changes: 2 additions & 2 deletions src/read_fargo3d.f90
Original file line number Diff line number Diff line change
Expand Up @@ -273,7 +273,7 @@ subroutine read_fargo3d_files()
! Calcul de la masse de gaz de la zone
mass = 0.
do icell=1,n_cells
mass = mass + densite_gaz(icell) * masse_mol_gaz * volume(icell)
mass = mass + densite_gaz(icell) * mu_mH * volume(icell)
enddo !icell
mass = mass * AU3_to_m3 * g_to_Msun

Expand All @@ -284,7 +284,7 @@ subroutine read_fargo3d_files()
! Somme sur les zones pour densite finale
do icell=1,n_cells
densite_gaz(icell) = densite_gaz(icell) * facteur
masse_gaz(icell) = densite_gaz(icell) * masse_mol_gaz * volume(icell) * AU3_to_m3
masse_gaz(icell) = densite_gaz(icell) * mu_mH * volume(icell) * AU3_to_m3
enddo ! icell
else
call error('Gas mass is 0')
Expand Down
4 changes: 2 additions & 2 deletions src/read_idefix.f90
Original file line number Diff line number Diff line change
Expand Up @@ -210,7 +210,7 @@ subroutine read_idefix_model()
! Calcul de la masse de gaz de la zone
mass = 0.
do icell=1,n_cells
mass = mass + densite_gaz(icell) * masse_mol_gaz * volume(icell)
mass = mass + densite_gaz(icell) * mu_mH * volume(icell)
enddo !icell
mass = mass * AU3_to_m3 * g_to_Msun

Expand All @@ -221,7 +221,7 @@ subroutine read_idefix_model()
! Somme sur les zones pour densite finale
do icell=1,n_cells
densite_gaz(icell) = densite_gaz(icell) * facteur
masse_gaz(icell) = densite_gaz(icell) * masse_mol_gaz * volume(icell) * AU3_to_m3
masse_gaz(icell) = densite_gaz(icell) * mu_mH * volume(icell) * AU3_to_m3
enddo ! icell
else
call error('Gas mass is 0')
Expand Down
4 changes: 2 additions & 2 deletions src/read_pluto.f90
Original file line number Diff line number Diff line change
Expand Up @@ -258,7 +258,7 @@ subroutine read_pluto_files()
! Calcul de la masse de gaz de la zone
mass = 0.
do icell=1,n_cells
mass = mass + densite_gaz(icell) * masse_mol_gaz * volume(icell)
mass = mass + densite_gaz(icell) * mu_mH * volume(icell)
enddo !icell
mass = mass * AU3_to_m3 * g_to_Msun

Expand All @@ -269,7 +269,7 @@ subroutine read_pluto_files()
! Somme sur les zones pour densite finale
do icell=1,n_cells
densite_gaz(icell) = densite_gaz(icell) * facteur
masse_gaz(icell) = densite_gaz(icell) * masse_mol_gaz * volume(icell) * AU3_to_m3
masse_gaz(icell) = densite_gaz(icell) * mu_mH * volume(icell) * AU3_to_m3
enddo ! icell
else
call error('Gas mass is 0')
Expand Down
8 changes: 4 additions & 4 deletions src/read_spherical_grid.f90
Original file line number Diff line number Diff line change
Expand Up @@ -191,8 +191,8 @@ subroutine read_spherical_model(filename)
! !-> wrapper for dust RT.
! !-> taking into account proper weights assuming only molecular gas in dusty regions
! if (rho_dust(i,j,k) > 0.0) then ! dusty region
! densite_gaz(icell) = rho(i,j,k) * 1d3 / masse_mol_gaz !total molecular gas density in H2/m^3
! densite_pouss(:,icell) = rho_dust(i,j,k) * 1d3 / masse_mol_gaz ! [m^-3]
! densite_gaz(icell) = rho(i,j,k) * 1d3 / mu_mH !total molecular gas density in H2/m^3
! densite_pouss(:,icell) = rho_dust(i,j,k) * 1d3 / mu_mH ! [m^-3]
! disk_zone(1)%diskmass = disk_zone(1)%diskmass + rho_dust(i,j,k) * volume(icell)
! else !No dust.
! densite_gaz(icell) = nHtot(icell) * wght_per_H !total atomic gas density in [m^-3]
Expand Down Expand Up @@ -226,8 +226,8 @@ subroutine read_spherical_model(filename)
!-> wrapper for dust RT.
!-> taking into account proper weights assuming only molecular gas in dusty regions
if (rho_dust(i,j,k) > 0.0) then ! dusty region
densite_gaz(icell) = rho(i,j,k) * 1d3 / masse_mol_gaz !total molecular gas density in H2/m^3
densite_pouss(:,icell) = rho_dust(i,j,k) * 1d3 / masse_mol_gaz ! [m^-3]
densite_gaz(icell) = rho(i,j,k) * 1d3 / mu_mH !total molecular gas density in H2/m^3
densite_pouss(:,icell) = rho_dust(i,j,k) * 1d3 / mu_mH ! [m^-3]
disk_zone(1)%diskmass = disk_zone(1)%diskmass + rho_dust(i,j,k) * volume(icell)
else !No dust.
densite_gaz(icell) = nHtot(icell) * wght_per_H !total atomic gas density in [m^-3]
Expand Down

0 comments on commit 9d2a986

Please sign in to comment.