heating.f90 2.23 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
!>@author Lena Noack
!>@date Autumn 2018
module heating
use precision
use parameters


implicit none
contains

  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
  !!  Radiogenic heating depending on time !!
  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
  function radiogenic_heating(H0,lambda,t,pI)
    real(dp),intent(in)::H0,lambda,t
    type(param_I),intent(in)::pI
    real(dp)::radiogenic_heating

    real(dp)::HtU238,HtU235,HtTh232,HtK40

    if (pI%C_U.eq.0.0_dp) then
Lena Noack's avatar
Lena Noack committed
22
      radiogenic_heating = H0*exp(-lambda*t) ! H0 is in nondimensional W/kg; hence radiogenic_heating is in ND W/kg
23
24
25
    else
      ! calculate from C_U, C_Th and C_K, lambda is 1 (cooling) or 0 (constant heat sources)

Lena Noack's avatar
Lena Noack committed
26
      HtU238 = pI%H0U238*Exp(-Log(2.0_dp)*t*lambda/pI%t12U238) ! in ND W/kg 
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
      HtU235 = pI%H0U235*Exp(-Log(2.0_dp)*t*lambda/pI%t12U235)
      HtTh232 = pI%H0Th232*Exp(-Log(2.0_dp)*t*lambda/pI%t12Th232)
      HtK40 = pI%H0K40*Exp(-Log(2.0_dp)*t*lambda/pI%t12K40)

      ! cell%h (here H0) traces enrichment/depletion factor if C_U>0, otherwise cell%h is already total heat sources
      !write(*,*) H0,pI%X0U238*HtU238 + pI%X0U235*HtU235 + pI%X0Th232*HtTh232 + pI%X0K40*HtK40
      radiogenic_heating = H0 * (pI%X0U238*HtU238 + pI%X0U235*HtU235 + pI%X0Th232*HtTh232 + pI%X0K40*HtK40) ! multiplication with rho directly in thermal solvers
    endif
  end function radiogenic_heating

  ! new version, allows for redistribution of radioactive nuclei
  ! is only called if (pX%dU238.ne.0.0_dp), else radiogenic_heating is used
  function radiogenic_nuclei(X_U238,X_U235,X_Th232,X_K40,t,pI)
    real(dp),intent(in)::X_U238,X_U235,X_Th232,X_K40,t
    type(param_I),intent(in)::pI
    real(dp)::radiogenic_nuclei

    real(dp)::HtU238,HtU235,HtTh232,HtK40

    ! calculate from C_U, C_Th and C_K, lambda is 1 (cooling) or 0 (constant heat sources)
    HtU238 = pI%H0U238*Exp(-Log(2.0_dp)*t*pI%lambda/pI%t12U238)
    HtU235 = pI%H0U235*Exp(-Log(2.0_dp)*t*pI%lambda/pI%t12U235)
    HtTh232 = pI%H0Th232*Exp(-Log(2.0_dp)*t*pI%lambda/pI%t12Th232)
    HtK40 = pI%H0K40*Exp(-Log(2.0_dp)*t*pI%lambda/pI%t12K40)

    radiogenic_nuclei = X_U238*HtU238 + X_U235*HtU235 + X_Th232*HtTh232 + X_K40*HtK40 ! multiplication with rho directly in thermal solvers
  end function radiogenic_nuclei

end module