139 lines
3.9 KiB
Fortran
139 lines
3.9 KiB
Fortran
subroutine mpartf(jatom,ion,indmol,t,u,dulog)
|
|
c =============================================
|
|
c
|
|
c yields partition functions with polynomial data from
|
|
c ref. Irwin, A.W., 1981, ApJ Suppl. 45, 621.
|
|
c ln u(temp)=sum(a(i)*(ln(temp))**(i-1)) 1<=a<=6
|
|
c
|
|
c Input:
|
|
c jatom = element number in periodic table
|
|
c ion = 1 for neutral, 2 for once ionized and 3 for twice ionized
|
|
c indmol= index of a molecular specie (Tsuji index)
|
|
c temp = temperature
|
|
c Output:
|
|
c u = partf.(linear scale) for iat,ion, or indmol, and temperature t
|
|
c dulog = d ln(u)/d ln(T)
|
|
c
|
|
c
|
|
implicit real*8 (a-h,o-z)
|
|
common/moldat/moltab,irwtab
|
|
real*8 a(6,3,92),aa(6),am(6,600)
|
|
dimension indtsu(324),irw(600),igle(28)
|
|
save iread,a,am
|
|
DATA IGLE/2,1,2,1,6,9,4,9,6,1,2,1,6,9,4,9,6,1,
|
|
* 10,21,28,25,6,25,28,21,10,21/
|
|
C
|
|
data indtsu / 2, 5, 12, 4, 8, 7, 6,
|
|
* 9, 11, 10, 29, 50, 59, 46, 133, 52, 19,
|
|
* 13, 42, 38, 39, 37, 44, 36, 14, 118, 33,
|
|
* 3, 16, 57, 32, 49, 60, 54, 41, 107,304,
|
|
* 148, 152, 153, 155, 303, 17, 24, 25, 28, 51,
|
|
* 112, 119, 102, 0, 21, 15, 43, 22,478, 64,
|
|
* 47, 65, 414, 61, 191, 62 ,109, 40, 66,214,
|
|
* 120*0, 30, 136*0/
|
|
c
|
|
data iread /0/
|
|
c
|
|
c read data if first call:
|
|
c
|
|
if(iread.ne.1) then
|
|
if(irwtab.eq.0) then
|
|
open(67,file= './data/irwin_orig.dat',status='old')
|
|
nmoli=66
|
|
else
|
|
open(67,file= './data/irwin_bc.dat',status='old')
|
|
nmoli=324
|
|
end if
|
|
read(67,*)
|
|
read(67,*)
|
|
do j=1,92
|
|
do i=1,3
|
|
if(j.eq.1.and.i.eq.3) go to 10
|
|
sp=float(j)+float(i-1)/100.
|
|
read(67,*) spec,aa
|
|
do k=1,6
|
|
a(k,i,j)=aa(k)
|
|
end do
|
|
10 continue
|
|
end do
|
|
end do
|
|
c
|
|
read(67,*)
|
|
read(67,*)
|
|
read(67,*)
|
|
do i=1,600
|
|
irw(i)=0
|
|
end do
|
|
do i=1,nmoli
|
|
read(67,*) spec,aa
|
|
indm=indtsu(i)
|
|
if(indm.gt.0) then
|
|
irw(indm)=i
|
|
do j=1,6
|
|
am(j,indm)=aa(j)
|
|
end do
|
|
end if
|
|
end do
|
|
close(67)
|
|
iread=1
|
|
endif
|
|
c
|
|
c evaluation of the partition function
|
|
c stop if T is out of limits of Irwin's tables
|
|
c
|
|
u=1.
|
|
dulog=0.
|
|
if(t.lt.1000.) then
|
|
stop 'partf; temp<1000 K'
|
|
else if(t.gt.16000.) then
|
|
c stop 'partf; temp>16000 K'
|
|
c write(6,601) t
|
|
c 601 format(' warning! T = ',f12.1, 'larger than 16000.'/)
|
|
if(indmol.eq.0) then
|
|
if(jatom.le.28.and.ion.le.jatom) u=igle(jatom-ion+1)
|
|
end if
|
|
return
|
|
endif
|
|
tl=log(t)
|
|
u=0.
|
|
c
|
|
c atomic species
|
|
c
|
|
if(jatom.gt.0.and.ion.gt.0) then
|
|
ulog= a(1,ion,jatom)+
|
|
* tl*(a(2,ion,jatom)+
|
|
* tl*(a(3,ion,jatom)+
|
|
* tl*(a(4,ion,jatom)+
|
|
* tl*(a(5,ion,jatom)+
|
|
* tl*(a(6,ion,jatom))))))
|
|
if(jatom.eq.5.and.ion.eq.3) ulog=1.
|
|
u=exp(ulog)
|
|
dulog= a(2,ion,jatom)+
|
|
* tl*(a(3,ion,jatom)*2.+
|
|
* tl*(a(4,ion,jatom)*3.+
|
|
* tl*(a(5,ion,jatom)*4.+
|
|
* tl*(a(6,ion,jatom)*5.))))
|
|
end if
|
|
c
|
|
c molecular species
|
|
c
|
|
if(indmol.gt.0) then
|
|
indm=indmol
|
|
if(irw(indm).gt.0) then
|
|
ulog= am(1,indm)+
|
|
* tl*(am(2,indm)+
|
|
* tl*(am(3,indm)+
|
|
* tl*(am(4,indm)+
|
|
* tl*(am(5,indm)+
|
|
* tl*(am(6,indm))))))
|
|
u=exp(ulog)
|
|
dulog= am(2,indm)+
|
|
* tl*(am(3,indm)*2.+
|
|
* tl*(am(4,indm)*3.+
|
|
* tl*(am(5,indm)*4.+
|
|
* tl*(am(6,indm)*5.))))
|
|
end if
|
|
end if
|
|
return
|
|
end
|