crp/OpenMC
1
0
Fork 0
mirror of https://github.com/openmc-dev/openmc.git synced 2026-07-29 06:35:48 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions 5

Completed free gas elastic scattering implementation. Changed K_BOLTZMANN to MeV/K.

Browse source
This commit is contained in:
Paul Romano 2011-10-25 11:24:40 -04:00
parent 7a5681c4c1
commit a1e34461cb
3 changed files with 27 additions and 19 deletions
Download patch file Download diff file Expand all files Collapse all files

5
src/constants.f90
View file

@ -14,12 +14,15 @@ module constants
MASS_PROTON = 1.00727646677, & ! mass of a proton
AMU = 1.66053873e-27, & ! 1 amu in kg
N_AVOGADRO = 0.602214179, & ! Avogadro's number in 10^24/mol
K_BOLTZMANN = 8.617342e-5, & ! Boltzmann constant in eV/K
K_BOLTZMANN = 8.617342e-11, & ! Boltzmann constant in MeV/K
INFINITY = huge(0.0_8), & ! positive infinity
ZERO = 0.0_8, &
ONE = 1.0_8, &
TWO = 2.0_8
! Monoatomic ideal-gas scattering treatment threshold
real(8), parameter :: FREE_GAS_THRESHOLD = 400.0
! Boundary conditions
integer, parameter :: &
& BC_TRANSMIT = 0, & ! Transmission boundary condition (default)

2
src/doppler.f90
View file

@ -42,7 +42,7 @@ contains
real(8) :: sigma ! broadened cross section at one point
! Determine alpha parameter -- have to convert k to MeV/K
alpha = A_target/(1.0e-6 * K_BOLTZMANN * T)
alpha = A_target/(K_BOLTZMANN * T)
! Allocate memory for x and assign values
n = size(energy)

39
src/physics.f90
View file

@ -533,6 +533,7 @@ contains
real(8) :: vel ! magnitude of velocity
real(8) :: v_n(3) ! velocity of neutron
real(8) :: v_cm(3) ! velocity of center-of-mass
real(8) :: v_t(3) ! velocity of target nucleus
real(8) :: u ! x-direction
real(8) :: v ! y-direction
real(8) :: w ! z-direction
@ -544,14 +545,17 @@ contains
! Neutron velocity in LAB
v_n = vel * p % uvw
! Sample velocity of target nucleus
call sample_target_velocity(p, nuc, v_t)
! Velocity of center-of-mass
v_cm = v_n/(awr + ONE)
v_cm = (v_n + awr*v_t)/(awr + ONE)
! Transform to CM frame
v_n = v_n - v_cm
! Find speed of neutron in CM
vel = norm2(v_n)
vel = sqrt(dot_product(v_n, v_n))
! Sample scattering angle
mu = sample_angle(rxn, p % E)
@ -590,14 +594,12 @@ contains
! for this method can be found in FRA-TM-123.
!===============================================================================
subroutine sample_target_velocity(nuc, E, uvw, v_target)
subroutine sample_target_velocity(p, nuc, v_target)
type(Nuclide), pointer :: nuc
real(8), intent(in) :: E
real(8), intent(in) :: uvw(3)
real(8), intent(out) :: v_target(3)
type(Particle), pointer :: p
type(Nuclide), pointer :: nuc
real(8), intent(out) :: v_target(3)
real(8) :: awr ! atomic weight ratio of target nucleus
real(8) :: u, v, w ! direction of target
real(8) :: kT ! equilibrium temperature of target in MeV
real(8) :: alpha ! probability of sampling f2 over f1
@ -610,14 +612,17 @@ contains
real(8) :: beta_vt_sq ! (beta * speed of target)^2
real(8) :: vt ! speed of target
! Copy atomic weight ratio for this nuclide
awr = nuc % awr
! Determine equilibrium temperature in MeV
kT = K_BOLTZMANN * nuc % temp * 1e-6
kT = K_BOLTZMANN * nuc % temp
! Check if energy is above threshold
if (p % E >= FREE_GAS_THRESHOLD * kT) then
v_target = ZERO
return
end if
! calculate beta
beta_vn = sqrt(awr*E/kT)
beta_vn = sqrt(nuc%awr * p%E / kT)
alpha = ONE/(ONE + sqrt(pi)*beta_vn/2.0)
@ -658,13 +663,13 @@ contains
! determine direction of target velocity based on the neutron's velocity
! vector and the sampled angle between them
u = uvw(1)
v = uvw(2)
w = uvw(3)
u = p % uvw(1)
v = p % uvw(2)
w = p % uvw(3)
call rotate_angle(u, v, w, mu)
! determine speed of target nucleus
vt = sqrt(beta_vt_sq*kT/awr)
vt = sqrt(beta_vt_sq*kT/nuc % awr)
! determine velocity vector of target nucleus
v_target(1) = u*vt