Calculates Fermi energy of a white dwarf star using ruby-dimensional module

calculate_fermi_energy_of_white_dwarf.rb

#!/usr/bin/env ruby

# Calculates the Fermi energy and temperature of a typical white dwarf star.

require "ruby-dimensional-standard"
require "ruby-dimensional-physics"

#define constants
G=Physics::Constant::GravitationalConstant
Msun=Physics::Constant::Msun
m_p=1.672621898e-27*:kg
m_e=m_p/2000
k_B=Physics::Constant::BoltzmannConstant
hbar=Physics::Constant::DiracConstant

#---------------------------------------------
#white dwarf star
#---------------------------------------------
puts "============================================"
puts "Example of white dwarf star"
puts "============================================"

#define WD radius/mass
Mwd=0.6*Msun
Rwd=8667*:km

#calculate densities
mean_nucleon_density = Mwd/(4*Math::PI*Rwd**3)
mean_nucleon_number_density = (mean_nucleon_density/m_p)
mean_electron_density = mean_nucleon_number_density/2
puts "Mean density = %.1e kg/cm3" % mean_nucleon_density.in(:kg*:cm**(-3))
puts "Mean number density = %.1e nucleons/cm3" % mean_nucleon_number_density.in(:cm**(-3))
puts "Mean electron density = %.1e electrons/cm3" % mean_electron_density.in(:cm**(-3))

#calculate Fermi energy
fermi_wave_number = (3 * Math::PI**2 * mean_electron_density)**(1.0/3.0)
fermi_energy = hbar**2 / (2*m_e) * fermi_wave_number**2
puts "E_K = %.1e erg (or %.1e K) := Fermi energy (of electron)" % [fermi_energy.in(:erg), (fermi_energy/k_B).in(:K)]

puts "  U = G*Mwd/Rwd = %.1e erg    := Gravitational potential energy (of single proton)" % (G*Mwd*m_p/Rwd).in(:erg)
puts " 2*K = %.1e erg ~= U          := virial theorem" % (2*(fermi_energy.in(:erg)))

#---------------------------------------------
#sodium
#---------------------------------------------
puts "============================================"
puts "Example of sodium "
puts "============================================"
sodium_density = 0.917*:g/:cm**3
sodium_atomic_mass = 23*:g/(6.02e23)
sodium_atomic_number = 11
sodium_number_density = (sodium_density/sodium_atomic_mass)
puts "Number density = %.1e 1/cm3" % sodium_number_density.in(:cm**-3)

sodium_fermi_wave_number = (3 * Math::PI**2 * sodium_number_density)**(1.0/3.0)
sodium_fermi_energy = hbar**2 / (2*m_e) * sodium_fermi_wave_number**2
puts "Fermi energy = %.1e erg (%.1e K)" % [sodium_fermi_energy.in(:erg),(sodium_fermi_energy/k_B).in(:K)]