SpectraRust/tlusty/extracted/incldy.f
2026-03-19 14:05:33 +08:00

75 lines
1.9 KiB
Fortran

SUBROUTINE INCLDY(NDPTH)
C ========================
C
C Read an initial model atmosphere from unit 8
C in Cloudy's format as provided by Katya Verner
C
INCLUDE 'IMPLIC.FOR'
INCLUDE 'BASICS.FOR'
INCLUDE 'ATOMIC.FOR'
INCLUDE 'MODELQ.FOR'
PARAMETER (MINPUT=6)
DIMENSION RS(MDEPTH),KKFIX0(MLEVEL)
COMMON POPUL0(MLEVEL,MDEPTH),ESEMAT(MLEVEL,MLEVEL),BESE(MLEVEL),
* TEMP0(MDEPTH),ELEC0(MDEPTH),DENS0(MDEPTH),PPL0(MDEPTH),
* DEPTH(MDEPTH),PPL(MDEPTH),POPLTE(MLEVEL),X(MINPUT)
dimension a(mlevel,mlevel),b(mlevel),iifor0(mlevel)
C
do iat=1,natom
kkfix0(iat)=iifix(iat)
iifix(iat)=0
end do
C
READ(8,*) NDPTH
NUMPAR=3
IF(NDPTH.gt.mdepth.and.intrpl.eq.0)
* CALL QUIT('ndpth.gt.mdepth in INCLDY',ndpth,mdepth)
NLEV0=NLEVEL
READ(8,*) TSTARY,RSTARY
IF(RSTARY.LT.1.E6) RSTARY=RSTARY*6.96E10
C
DO ID=NDPTH,1,-1
READ(8,*) (X(I),I=1,MINPUT)
RS(ID)=X(1)
TEMP0(ID)=X(2)
ELEC0(ID)=X(4)
AN=X(3)/BOLK/TEMP0(ID)
DENS0(ID)=WMM(ID)*X(3)
TEMP(ID)=TEMP0(ID)
ELEC(ID)=ELEC0(ID)
DENS(ID)=DENS0(ID)
ANMA(ID)=DENS(ID)/WMM(ID)
ANTO(ID)=ANMA(ID)+ELEC(ID)
END DO
DM(1)=0.
DEPTH(1)=0.
DO ID=2,NDPTH
DDDM=(DENS(ID-1)+DENS(ID))*(RS(ID-1)-RS(ID))
DM(ID)=DM(ID-1)+0.5*DDDM
DEPTH(ID)=DM(ID)
END DO
C
RRDIL=(RSTARY/RS(NDPTH))*(RSTARY/RS(NDPTH))
TEMPBD=TSTARY
C
DO ID=1,NDPTH
CALL WNSTOR(ID)
CALL SABOLF(ID)
DO I=1,NLEV0
IIFOR0(I)=I
END DO
CALL RATMAT(ID,IIFOR0,-1,A,B)
CALL LEVSOL(A,B,POPLTE,IIFOR0,NLEV0,1)
DO I=1,NLEV0
POPUL0(I,ID)=POPLTE(I)
END DO
END DO
C
INTRPL=0
do iat=1,natom
iifix(iat)=kkfix0(iat)
end do
C
RETURN
END