Line data Source code
1 : ! ***********************************************************************
2 : !
3 : ! Copyright (C) 2010-2021 The MESA Team
4 : !
5 : ! This program is free software: you can redistribute it and/or modify
6 : ! it under the terms of the GNU Lesser General Public License
7 : ! as published by the Free Software Foundation,
8 : ! either version 3 of the License, or (at your option) any later version.
9 : !
10 : ! This program is distributed in the hope that it will be useful,
11 : ! but WITHOUT ANY WARRANTY; without even the implied warranty of
12 : ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
13 : ! See the GNU Lesser General Public License for more details.
14 : !
15 : ! You should have received a copy of the GNU Lesser General Public License
16 : ! along with this program. If not, see <https://www.gnu.org/licenses/>.
17 : !
18 : ! ***********************************************************************
19 :
20 : module magnetic_diffusion
21 :
22 : use const_def, only: dp, pi, pi4, clight, boltz_sigma, one_third
23 : use num_lib
24 : use math_lib
25 : use utils_lib
26 :
27 : implicit none
28 :
29 : private
30 : public :: calc_sige
31 : public :: calc_eta
32 :
33 : contains
34 :
35 : !> Compute the magnetic diffusivity from the electric conductivity.
36 : !! @param sig The electrical conductivity (1/s).
37 : !! @param eta The magnetic diffusivity (output, cm^2/s).
38 0 : real(dp) function calc_eta(sig) result(eta)
39 : real(dp), intent(in) :: sig
40 :
41 0 : eta = (clight * clight / (4d0 * pi)) /sig
42 0 : end function calc_eta
43 :
44 : !> Computes the electrical conductivity following
45 : !! S.-C. YOON Oct. 10, 2003.
46 : !!
47 : !! @param abar The mean atomic mass number.
48 : !! @param zbar The mean atomic charge.
49 : !! @param rho The density (g/cm^3).
50 : !! @param T The temperature (K).
51 : !! @param Cp The specific heat at constant pressure (erg/g/K).
52 : !! @param kap_cond The electronic thermal opacity (cm^2/g).
53 : !! @param opacity The opacity (cm^2/g).
54 : !! @param sig The electrical conductivity (output, 1/s).
55 0 : real(dp) function calc_sige(abar, zbar, rho, T, Cp, kap_cond, opacity) result(sig)
56 : real(dp), intent(in) :: abar, zbar, rho, T, Cp, kap_cond, opacity
57 0 : real(dp) :: gamma, xlg, alpha, ffff, xxx, xsig1, xsig2, xsig3
58 :
59 0 : gamma = 0.2275d0*pow2(zbar) * pow(rho * 1.d-6 / abar, one_third)*1.d8/T
60 0 : xlg = log10(gamma)
61 :
62 0 : alpha = 16d0 * boltz_sigma * pow3(T) / (3d0 * opacity * pow2(rho) * Cp)
63 :
64 0 : if (xlg < -1.5d0) then
65 0 : sig = sige1(zbar,T,gamma)
66 0 : else if (xlg >= -1.5d0 .and. xlg <= 0d0) then
67 0 : xxx = (xlg + 0.75d0)*4d0/3d0
68 0 : ffff = 0.25d0*(2d0-3d0*xxx + pow3(xxx))
69 0 : xsig1 = sige1(zbar,T,gamma)
70 0 : xsig2 = sige2(T,rho,kap_cond)
71 0 : sig = (1d0-ffff)*xsig2 + ffff*xsig1
72 0 : else if (xlg > 0d0 .and. xlg < 0.5d0) then
73 0 : xsig2 = sige2(T,rho,kap_cond)
74 0 : sig = xsig2
75 0 : else if (xlg >= 0.5d0 .and. xlg < 1d0) then
76 0 : xxx = (xlg-0.75d0)*4d0
77 0 : ffff = 0.25d0*(2d0-3d0*xxx + pow3(xxx))
78 0 : xsig2 = sige2(T,rho,kap_cond)
79 0 : xsig3 = sige3(zbar,T,gamma)
80 0 : sig = (1d0-ffff)*xsig3 + ffff*xsig2
81 : else
82 0 : sig = sige3(zbar,T,gamma)
83 : end if
84 :
85 0 : end function calc_sige
86 :
87 : !> Computes one regime of the electrical conductivity.
88 : !! Written by S.-C. Yoon, Oct. 10, 2003
89 : !! See also Spitzer 1962 and Wendell et al. 1987, ApJ 313:284
90 : !! @param Z species charge
91 : !! @param T Temperature (K)
92 : !! @param xgamma The ion coupling strength (dimensionless).
93 : !! @param sige1 The electrical conductivity (1/s).
94 0 : real(dp) function sige1(z,t,xgamma)
95 : real(dp), intent(in) :: z, t, xgamma
96 0 : real(dp) :: etan, xlambda,f
97 0 : if (t >= 4.2d5) then
98 0 : f = sqrt(4.2d5/t)
99 : else
100 : f = 1.d0
101 : end if
102 0 : xlambda = sqrt(3d0*z*z*z)*pow(xgamma,-1.5d0)*f + 1d0
103 0 : etan = 3.d11*z*log(xlambda)*pow(t,-1.5d0) ! magnetic diffusivity
104 0 : etan = etan/(1.d0-1.20487d0*exp(-1.0576d0*pow(z,0.347044d0))) ! correction: gammae
105 0 : sige1 = clight*clight/(pi4*etan) ! sigma = c^2/(4pi*eta)
106 0 : end function sige1
107 :
108 : !> Computes one regime of the electrical conductivity using the conductive opacity.
109 : !! Written by S.-C. Yoon, Oct. 10, 2003
110 : !! See Wendell et al. 1987, ApJ 313:284
111 : !! @param T Temperature (K)
112 : !! @param rho Temperature (g/cm^3)
113 : !! @param kap_cond The electronic thermal opacity (cm^2/g).
114 : !! @param sige2 The electrical conductivity (1/s).
115 0 : real(dp) function sige2(T,rho,kap_cond)
116 : real(dp), intent(in) :: t,rho,kap_cond
117 0 : sige2 = 1.11d9*T*T/(rho*kap_cond)
118 : end function sige2
119 :
120 : !> Computes the electrical conductivity in degenerate matter.
121 : !! Written by S.-C. Yoon, Oct. 10, 2003
122 : !! See Nandkumar & Pethick (1984)
123 : !! @param Z species charge
124 : !! @param T Temperature (K)
125 : !! @param xgamma The ion coupling strength (dimensionless).
126 : !! @param sige3 The electrical conductivity (1/s).
127 0 : real(dp) function sige3(z,t,xgamma)
128 : real(dp), intent(in) :: z, t, xgamma
129 0 : real(dp) :: rme, rm23, ctmp, xi
130 0 : rme = 8.5646d-23*t*t*t*xgamma*xgamma*xgamma/pow5(z) ! rme = rho6/mue
131 0 : rm23 = pow(rme,2d0/3d0)
132 0 : ctmp = 1d0 + 1.018d0*rm23
133 0 : xi= sqrt(pi/3.0d0)*log(z)/3.d0 + 2.d0*log(1.32d0+2.33d0/sqrt(xgamma))/3.d0-0.484d0*rm23/ctmp
134 0 : sige3 = 8.630d21*rme/(z*ctmp*xi)
135 0 : end function sige3
136 :
137 : end module magnetic_diffusion
|
