4 double precision,
private :: Te_H_CL(1:100)
5 double precision,
private :: iz_H_CL(1:100)
6 double precision,
dimension(:),
allocatable ::
te_h_table
7 double precision,
dimension(:),
allocatable ::
iz_h_table
10 data te_h_cl /2.1000005d+03, 2.2349495d+03, 2.3785708d+03, 2.5314199d+03, 2.6940928d+03, &
11 2.8672190d+03, 3.0514692d+03, 3.2475613d+03, 3.4562544d+03, 3.6783584d+03, &
12 3.9147351d+03, 4.1662998d+03, 4.4340322d+03, 4.7189697d+03, 5.0222153d+03, &
13 5.3449502d+03, 5.6884248d+03, 6.0539712d+03, 6.4430083d+03, 6.8570420d+03, &
14 7.2976860d+03, 7.7666460d+03, 8.2657383d+03, 8.7969062d+03, 9.3622090d+03, &
15 9.9638330d+03, 1.0604130d+04, 1.1285561d+04, 1.2010781d+04, 1.2782619d+04, &
16 1.3604041d+04, 1.4478250d+04, 1.5408652d+04, 1.6398826d+04, 1.7452648d+04, &
17 1.8574172d+04, 1.9767766d+04, 2.1038082d+04, 2.2390010d+04, 2.3828838d+04, &
18 2.5360100d+04, 2.6989764d+04, 2.8724182d+04, 3.0570023d+04, 3.2534518d+04, &
19 3.4625215d+04, 3.6850262d+04, 3.9218293d+04, 4.1738543d+04, 4.4420699d+04, &
20 4.7275266d+04, 5.0313219d+04, 5.3546391d+04, 5.6987395d+04, 6.0649453d+04, &
21 6.4546914d+04, 6.8694758d+04, 7.3109148d+04, 7.7807289d+04, 8.2807258d+04, &
22 8.8128625d+04, 9.3791852d+04, 9.9819000d+04, 1.0623346d+05, 1.1306024d+05, &
23 1.2032560d+05, 1.2805797d+05, 1.3628709d+05, 1.4504503d+05, 1.5436592d+05, &
24 1.6428561d+05, 1.7484295d+05, 1.8607852d+05, 1.9803609d+05, 2.1076208d+05, &
25 2.2430608d+05, 2.3872045d+05, 2.5406084d+05, 2.7038703d+05, 2.8776234d+05, &
26 3.0625456d+05, 3.2593475d+05, 3.4688000d+05, 3.6917081d+05, 3.9289406d+05, &
27 4.1814181d+05, 4.4501247d+05, 4.7360988d+05, 5.0404450d+05, 5.3643488d+05, &
28 5.7090725d+05, 6.0759431d+05, 6.4663956d+05, 6.8819319d+05, 7.3241712d+05, &
29 7.7948294d+05, 8.2957412d+05, 8.8288331d+05, 9.3961919d+05, 1.0000000d+06 /
32 data iz_h_cl /1.0000000d+00, 1.0000000d+00, 1.0000000d+00, 1.0000000d+00, 1.0000000d+00, &
33 1.0000000d+00, 1.0000000d+00, 9.9999928d-01, 9.9999774d-01, 9.9999624d-01, &
34 9.9999428d-01, 9.9999219d-01, 9.9998862d-01, 9.9998313d-01, 9.9997944d-01, &
35 9.9995774d-01, 1.0000324d+00, 9.9822283d-01, 9.6391600d-01, 8.8995951d-01, &
36 8.0025935d-01, 7.4219799d-01, 7.1378660d-01, 6.9565988d-01, 6.7326778d-01, &
37 6.5032905d-01, 6.2487972d-01, 5.9925246d-01, 5.7028764d-01, 5.3673893d-01, &
38 4.9144468d-01, 4.4332290d-01, 3.9158964d-01, 3.4621361d-01, 2.7826238d-01, &
39 1.9633104d-01, 1.2593637d-01, 8.1383914d-02, 5.0944358d-02, 2.8770180d-02, &
40 1.4449068d-02, 7.6071853d-03, 4.3617180d-03, 2.7694483d-03, 1.8341533d-03, &
41 1.2467797d-03, 8.6159195d-04, 6.0232455d-04, 4.2862241d-04, 3.0666622d-04, &
42 2.2256430d-04, 1.6328457d-04, 1.2179965d-04, 9.2065704d-05, 7.0357783d-05, &
43 5.4168402d-05, 4.2026968d-05, 3.2807184d-05, 2.5801779d-05, 2.0414085d-05, &
44 1.6271379d-05, 1.3041737d-05, 1.0538992d-05, 8.5329248d-06, 6.9645471d-06, &
45 5.7069005d-06, 4.7155136d-06, 3.9203474d-06, 3.2734958d-06, 2.7415674d-06, &
46 2.2911979d-06, 1.9045801d-06, 1.5692771d-06, 1.2782705d-06, 1.0281080d-06, &
47 8.1567373d-07, 6.4200310d-07, 5.0668103d-07, 4.0577518d-07, 3.3245340d-07, &
48 2.8006056d-07, 2.4119515d-07, 2.1033340d-07, 1.8429903d-07, 1.6189058d-07, &
49 1.4264889d-07, 1.2616209d-07, 1.1193421d-07, 9.9408616d-08, 8.8335305d-08, &
50 7.8593857d-08, 7.0042269d-08, 6.2532941d-08, 5.5898361d-08, 5.0049657d-08, &
51 4.4888207d-08, 4.0317172d-08, 3.6192560d-08, 3.2330732d-08, 2.9387850d-08 /
58 double precision :: fact1, fact2, fact3, dL1, dL2
59 integer :: i, j, ntable, n_interpolated_table
63 n_interpolated_table=4000
74 te_table_step = (te_h_cl(ntable)-te_h_cl(1))/dble(n_interpolated_table-1)
79 te_h_table(n_interpolated_table) = te_h_cl(ntable)
80 iz_h_table(n_interpolated_table) = iz_h_cl(ntable)
82 do i=2,n_interpolated_table
88 if(j.eq. ntable-1 )
then
90 /(te_h_cl(j)-te_h_cl(j+1))
93 /(te_h_cl(j+1)-te_h_cl(j))
95 iz_h_table(i) = iz_h_cl(j)*fact1 + iz_h_cl(j+1)*fact2
98 dl1 = iz_h_cl(j+1)-iz_h_cl(j)
99 dl2 = iz_h_cl(j+2)-iz_h_cl(j+1)
100 jump =(max(dabs(dl1),dabs(dl2)) > 2.d0*min(dabs(dl1),dabs(dl2)))
105 /(te_h_cl(j)-te_h_cl(j+1))
108 /(te_h_cl(j+1)-te_h_cl(j))
110 iz_h_table(i) = iz_h_cl(j)*fact1 + iz_h_cl(j+1)*fact2
115 / ((te_h_cl(j)-te_h_cl(j+1)) &
116 * (te_h_cl(j)-te_h_cl(j+2)))
120 / ((te_h_cl(j+1)-te_h_cl(j)) &
121 * (te_h_cl(j+1)-te_h_cl(j+2)))
125 / ((te_h_cl(j+2)-te_h_cl(j)) &
126 * (te_h_cl(j+2)-te_h_cl(j+1)))
128 iz_h_table(i) = iz_h_cl(j)*fact1 + iz_h_cl(j+1)*fact2 &
149 integer,
intent(in) :: ixI^L, ixO^L
150 double precision,
intent(in) :: Te(ixI^S)
151 double precision,
intent(out) :: iz_H(ixO^S),iz_He(ixO^S)
155 {
do ix^db=ixomin^db,ixomax^db\}
167 iz_he(ix^d)=1.0084814d-4
178 double precision,
intent(in) :: te
181 double precision :: iz_h, iz_he
This module contains definitions of global parameters and variables and some generic functions/subrou...
integer, dimension(:), allocatable, parameter d
double precision unit_temperature
Physical scaling factor for temperature.
module ionization degree - get ionization degree for given temperature
double precision, dimension(:), allocatable te_h_table
double precision te_table_max
double precision, dimension(:), allocatable iz_h_table
double precision function r_ideal_gas_law_partial_ionization(Te)
subroutine ionization_degree_init()
subroutine ionization_degree_from_temperature(ixIL, ixOL, Te, iz_H, iz_He)
double precision inv_te_table_step
double precision te_table_min
double precision te_low_iz_he
double precision te_table_step