C
C
C
      SUBROUTINE  TEVAL
C
C     ------------------------------------------------------------------
C
C       SUBROUTINE PURPOSES:
C       1) EVALUATION OF MEAN AND MEAN-SQUARE AIR TEMP.
C       2) CALCULATION OF THE INTEGRAL OF AIR DENSITY*ELEVATION
C          VARIATION FOR EACH AIRWAY.
C       3) EVALUATION OF VARIED AIRFLOW RESISTANCE DUE TO VARIED AIR
C          TEMP.
C       THIS SUBROUTINE IS CALLED FOR IN THE QUASI-EQUILIBRIUM
C     SIMULATION PART OF THE PROGRAM. ITS COUNTERPART IN THE DYNAMIC
C     SIMULATION PART IS PARAM.
C
C     ------------------------------------------------------------------
C
      INCLUDE 'CMMNW.DAT'
C
      DO 40 I=1,NB
         CP=0.2376+0.000012*(TRD(I)+TJS(I))
         E=DZRD(I)/(2.*778.26*CP)
         B=TJS(I)-TROCK(I)
         G=XNEW(I)*LA(I)/CP
         IF (G.GE.25.0) THEN
            GX=0.0
            GXX=0.0
         ELSE
            GX=EXP(-G)
            GXX=EXP(-2.*G)
         ENDIF
         F1=ABS(TJS(I)-TROCK(I))
         F2=ABS(TRD(I)-TROCK(I))
         FF=MAX(F1,F2)
         IF (FF.LT.0.1) THEN
            TMRD(I)=TJS(I)
            FRNVP(I)=2116.0*DZRD(I)/(53.352*(TJS(I)+460.0))
C         write (8,*) 'ifrnvpdztjs11',i,frnvp(i),dzrd(i),tjs(i)
         ELSE
            TMRD(I)=TROCK(I)-B*(GX-1.0)/G-E*(0.5-(GX+(GX-1.0)/G)/G)
            TTROCK=TROCK(I)
            IF (ABS(TJS(I)-TTROCK).LT.0.01) TTROCK=TJS(I)-5.0
            TT=(TRD(I)-TTROCK)/(TJS(I)-TTROCK)
            IF (TT.LT.0.01) THEN
               TTROCK=AMIN1(TJS(I),TRD(I))-5.0
               TT=(TRD(I)-TTROCK)/(TJS(I)-TTROCK)
            ENDIF
            IF (ABS(TT-1.0).LT.1.E-3) THEN
               FRNVP(I)=2116.0*DZRD(I)/(53.352*(TJS(I)+460.0))
C              write (8,*) 'ifrnvpdztjs22',i,frnvp(i),dzrd(i),tjs(i)
            ELSE
               YQ=ALOG((TRD(I)+460.0)/(TJS(I)+460.0))/ALOG(TT)
C                     write (8,*) 'yqtrdtjstt',yq,trd(i),tjs(i),tt
               FRNVP(I)=2116.0*DZRD(I)*(1.0-YQ)/(53.352*(TTROCK+460.))
C             write (8,*) 'ifrnvpdztjs33',i,frnvp(i),dzrd(i),yq,ttrock
            ENDIF
         ENDIF
         TRA1=TROCK(I)+460.0
         A1=2*E*E/(2*G)**3+2*E*B/(4*G*G)-B*B/(2*G)
         A2=4*E*E/G**3+(920.*E+2*E*TJS(I))/(G*G)-2*TRA1*B/G
         A3=2*E*E/G+4*E*E/(G*G)+(920.*E+2*E*TJS(I))/G
         A4=E*E/(2*G)+2*E*E/(4*G*G)+E*B/G
         TMSQR=TRA1*TRA1-TRA1*E-E*LA(I)/3.0+A1*(GXX-1.0)
     .            +A2*(GX-1.0)+A3*GX+A4*GXX
C
C       AIRFLOW RESISTANCE FOR A FAN BRANCH (NWTYP(I)=1) IS REGARDED
C     AS ZERO AND R FOR THAT BRANCH STANDS FOR FAN PRESSURE HEAD.
C
         IF (NWTYP(I).NE.1) R(I)=RSTD(I)*TMSQR/(460.+TR)**2
40    CONTINUE
C
      RETURN
      END

