C
C
C
      SUBROUTINE  CDJN (ISQC,MARKS,NSTOP)
C
C     ------------------------------------------------------------------
C
C       SUBROUTINE PURPOSES:
C       1) EVALUATION OF THERMAL STATES IN JUNCTIONS.
C       2) DEVELOPING A NEW DATA RECORD INTO AIRWAYS (A) WHOSE BEGINNING
C          JUNCTION HAS DRASTIC CONDITION CHANGES; OR (B) WHERE AIRFLOW
C          REVERSAL HAPPENED; OR (C) WHERE A FIRE SOURCE EXISTS; OR (D)
C          WHICH GOT DETAILED RECORDS IN PREVIOUS INTERVALS ALREADY.
C
C     ------------------------------------------------------------------
C
      INCLUDE 'CMMNW.DAT'
C
C  ** UPDATING JUNCTION CONDITIONS.
C
      CP=0.24
      DO 60 I=1,NJ
         K=JLR(I)
C
C  ** DEVELOPING DATA POINTS INTO THE DOWN-STREAM AIRWAYS.
C
         JY=1
         IF (I.GT.1) JY=LOUT(I-1)+1
         JZ=LOUT(I)
         DO 50 J=JY,JZ
            NW=NGOUT(J)
            CO2=0.0
            H2O=0.0
            WR=0.0
C
C     NWTYP(I)=10 INDICATES THAT AIRWAY NO(I) CONTAINS A FIRE SOURCE.
C     NSAC(I,2): NUMBER OF CONTROL VOLUMES IN AIRWAY NO(I).
C     IBTN(I) STORES THE ADDRESS OF DATA RECORDS FOR AIRWAY NO(I) IN
C              ARRAY DPPA, DPPB ETC.
C
            IF (A(NW).GE.100000.0) GO TO 50
C                        ********
            IF (MARKS.NE.0.AND.NWTYP(NW).NE.10) THEN
               IF (JCH(K).EQ.0) GO TO 50
C
C       WHEN THE NEWLY OBTAINED DATA NEED TO BE UPDATED AGAIN IN THE
C     PRESENT TIME INTERVAL, THE DATA OBTAINED IN THE PREVIOUS SEARCH
C     IN THE INTERVAL ARE ABANDONED TO AVOID DUPLICATED DATA
C     RECORDS.
C
               IF (IBTN(NW).GT.0) THEN
                  IZ=IBTN(NW)
                  NSAC(NW,2)=NSACB(NW,2)
                  DO 5 IX=1,10
                     DO 5 IY=1,5
                        DPPA(IZ,IX,IY)=DPPB(IZ,IX,IY)
5                 CONTINUE
               ENDIF
            ELSE IF (NWTYP(NW).NE.10) THEN
C
C     NON-QUALIFIED AIRWAYS GET NO DETAILED DATA RECORDS AND ARE IGNORED
C
               IF (JDP.GT.0) THEN
                  DO 7 IX=1,JDP
                     IF (NO(NW).EQ.JDPP(IX)) GO TO 8
7                 CONTINUE
               ENDIF
               IF (JCH(K).EQ.0.AND.IBTN(NW).EQ.0) GO TO 50
            ENDIF
C
C  ** CHECK TO SEE IF A FIRE SOURCE EXISTS IN THE AIRWAY.
C
8           DO 10 L=1,INFLOW
               IF (NCENT(L).EQ.NO(NW)) GO TO 20
10          CONTINUE
C
C     NO FIRE SOURCE CONTAINED IN THE AIRWAY.
C
            CONTAM=0.0
            CONTQ=0.0
            HEATAD=0.0
            GO TO 30
C
C     HAVING A FIRE SOURCE WITH USER-SPECIFIED HEAT GENERATION RATE.
C
20          FACT=1.0
            IF (QCENT(L).GT.10.0) THEN
               Q1=QQ2(NW)
               IF (Q1.LT.QCENT(L)) Q1=Q1+QCENT(L)
               FACT=1.0+(QQ2(NW)-QCENT(L))/Q1
            ENDIF
            IF (ABS(CONT(L)).GT.1.E-3) THEN
               CONTAM=CONT(L)*CONC(L)/100.0
               CONTQ=CONT(L)
            ELSE
               CONTAM=0.0
               CONTQ=0.0
            ENDIF
            O2BEH(L)=(0.21-PROP(K)-CONTAM/QQ2(NW))*100.0
            HEATAD=HEAT(L)*FACT
C
C     OXYGIN RICH FIRE
C
            IF (O2MIN(L).GT.1.0.AND.HEAT(L).GT.(-1.E-5)) THEN
               O2SCL=O2MIN(L)+(21.0-O2MIN(L))*(FACT-1.0)
               TACCA=TACC+DELTAT
               IF (TACCA.GT.TPR(L)*60.0) TACCA=TPR(L)*60.0
               PROPA=PROP(K)
               IF (NHM.GT.0) PROPA=PROP(K)*(CN(L)*32.0+HN(L)*8.0)/
     .                             (CN(L)*32.0)
               O2CONS=((0.21-PROPA-O2SCL/100.0)*TACCA/
     .                (TPR(L)*60.0))*QQ2(NW)
               IF (O2CONS.LT.0.0) O2CONS=0.0
               O2BEH(L)=(0.21-PROPA-O2CONS/QQ2(NW))*100.0
               IF (NHM.GT.0) THEN
                  CO2=O2CONS*CN(L)*32.0/(CN(L)*32.0+HN(L)*8.0)
                  H2O=O2CONS*HN(L)*8.0/(CN(L)*32.0+HN(L)*8.0)
                  WR=2.0*H2O*18.0/(28.9*QQ2(NW))
                  HEATAD=O2CONS*437.0
                  CONTQ=0.0
                  CONTAM=CO2
               ELSE
                  WR=0.0
                  CONTQ=0.0
                  HEATAD=O2CONS*437.0
               ENDIF
               CONTAM=PROPA*QQ2(NW)+O2CONS
            ENDIF
C
C     FUEL RICH FIRE, THE MODELS OF THE ABOVE THREE TYPES OVERRIDE
C     THE ONES IN FRONT.
C
            IF (SMPO2(L).GT.1.E-3.AND.HEAT(L).GT.(-1.E-5)) THEN
               O2LT=0.21-0.21*(TACC+DELTAT)/(TPR(L)*60.0)
               IF (O2LT.LT.0.0) O2LT=0.0
               PROPA=PROP(K)
               IF (NHM.GT.0) PROPA=PROP(K)*(CN(L)*32.0+HN(L)*8.0)/
     .                             (CN(L)*32.0)
               O2CONS=(0.21-O2LT-PROPA)*QQ2(NW)
               IF (O2CONS.LT.0.0) O2CONS=0.0
               CONTAM=(0.21-O2LT)*QQ2(NW)
               IF ((0.21-O2LT).LT.PROPA) CONTAM=PROPA*QQ2(NW)
               O2BEH(L)=O2LT*100.0
               IF (FACT.LE.1.0) THEN
                  HTSCL=HTPO2(L)
               ELSE
                  HTSCL=HTPO2(L)/FACT
               ENDIF
               IF (NHM.GT.0) THEN
                  CO2=CONTAM*CN(L)*32.0/(CN(L)*32.0+HN(L)*8.0)
                  H2O=CONTAM*HN(L)*8.0/(CN(L)*32.0+HN(L)*8.0)
                  WR=2.0*H2O*18.0/(28.98*QQ2(NW))
                  HEATAD=O2CONS*HTSCL
                  CONTQ=0.0
                  CONTAM=CO2
               ELSE
                  WR=0.0
                  CONTQ=0.0
                  HEATAD=O2CONS*HTSCL
               ENDIF
            ENDIF
            GFS(L)=PRH2O(K)+WR
C
C  ** EVALUATION OF AIR TEMP. TJS(I) AT THE BEGINNING OF AIRWAY NO(I).
C
30          TK=T(K)
            IF (NWTYP(NW).EQ.10) TK=TAS(K)
            VART=(9900.0+TK)**2+2.0*HEATAD/(QQ2(NW)*2.4E-5*DR)
            IF (VART.LT.9.4E7) THEN
               WRITE (8,110) NO(I),Q(I)
               NSTOP=1
               RETURN
            ENDIF
            TJS(NW)=-9900.0+SQRT(VART)
            IF (TJS(I).LT.(-200.0).OR.TJS(I).GT.5000.0) THEN
               WRITE (8,100) NO(I),Q(I),TJS(I)
               NSTOP=1
               RETURN
            ENDIF
C            HTAD(L)=(TJS(NW)-T(K))*0.24*QQ2(NW)*DR
C            ************************************** 
            IF (NWTYP(NW).EQ.10) THEN
               TRD(NW)=TJS(NW)
               RDPROP(NW)=PROP(K)*(QQ2(NW)-CONTQ)/QQ2(NW)+CONTAM/QQ2(NW)
               RDCH4(NW)=PRCH4(K)
               HTAD(L)=(TJS(NW)-T(K))*0.24*QQ2(NW)*DR
C              *************************************
               GO TO 50
            ENDIF
C
C  ** CHECK AND PREPARATION OF DATA STORAGE ARRAYS.
C
            NTAUXC=0
            IF (NSAC(NW,1).EQ.0) THEN
C
C     NSAC(I,1)=1 INDICATES THAT THE FUME FRONT IS WITHIN AIRWAY NO(I).
C
               NTAUXC=1
               NSAC(NW,1)=1
            ENDIF
            NSAC(NW,2)=NSAC(NW,2)+1
            NS=NSAC(NW,2)
            IF (IBTN(NW).EQ.0) THEN
               ISQC=ISQC+1
               IF (ISQC.LT.NMZ) THEN
                  IBTN(NW)=ISQC
C
C       NSEQE(I): A BRIDGE ARRAY INDICATING THE LOCATION OF DATA FOR
C     WET CASE TEMP. CORR. CALCULATION FOR AIRWAY NO(I).
C
                  IF (ISQC.EQ.1) NSEQE(ISQC)=0
                  NSEQE(ISQC+1)=NSEQE(ISQC)+NSP(NW)
               ELSE
C
C       WHEN THE CAPACITIES OF DATA RECORDING ARRAYS ARE SATURATED, LATE
C     REQUEST WILL BE TURNED DOWN AND A MESSAGE WILL BE GIVEN.
C
                  NSAC(NW,1)=0
                  NSAC(NW,2)=0
                  IF (IAC.GT.0) THEN
                     DO 35 I1=1,IAC
                        IF (IOMIT(I1).EQ.NW) GO TO 36
35                   CONTINUE
                  ENDIF
                  IF (IAC.GE.IMX) GO TO 50
                  IAC=IAC+1
                  I1=IAC
36                IOMIT(I1)=NW
                  ROMIT(I1,1)=T(K)
                  ROMIT(I1,2)=PROP(K)
                  ROMIT(I1,3)=PRCH4(K)
                  GO TO 50
               ENDIF
            ENDIF
            NQ=IBTN(NW)
            NN=0
            IF (DCOAGE(NW).LT.0.0) THEN
               TMST=TACC+DELTAT
               CALL KALPHA (TMST,TM,FX,CP,NW,1)
               COAGE=BI(NW)-FX*BI(NW)*BI(NW)/(0.375+BI(NW))
               DCOAGE(NW)=COAGE
               XNEW(NW)=HK(NW)*O(NW)*O(NW)*COAGE/(120.*DR*QQ2(NW)*A(NW))
            ENDIF
C
C  ** ADVANCING DATA POINTS.
C
            IF (NHM.EQ.0) THEN
C
C     PURELY DRY CASE IS ASSUMED WHEN NHM=0.
C
               TS=DELTAT/60.0
               IF (TLEFT(K).GE.0.0) TS=TLEFT(K)+0.00001
               DPPA(NQ,NS,5)=DR*QQ2(NW)*TS
               IF (NS.GT.1.AND.DELTAT/60.0.GT.TS) THEN
                  DPPA(NQ,NS-1,5)=DPPB(NQ,NS-1,5)+DR*QQ2(NW)*
     .                            (DELTAT/60.0-TS)
               ENDIF
               TAUXA(NQ,NS)=TJS(NW)
               CALL COLUMN (TJS(NW),SEG,NW,NQ,NS)
               FE=SEG/LA(NW)
               FA=XNEW(NW)*SEG/CP
               IF (FA.LT.25.0) THEN
                  FAS=EXP(-FA)
               ELSE
                  FAS=0.0
               ENDIF
C
C     AIR TEMP. OF THE C.V. AT ITS NEWLY ARRIVED POSITION IS ESTIMATED.
C
               DPPA(NQ,NS,2)=TROCK(NW)+(TJS(NW)-TROCK(NW))*FAS-
     .                       (DZRD(NW)*FE/(2.0*778.26*CP))*(1.0+FAS)
            ELSE
C
C  ** WET CASE ASSUMED.
C
               TS=DELTAT/60.0
               IF (TLEFT(K).GE.0.0) TS=TLEFT(K)+0.00001
               DPPA(NQ,NS,5)=DR*QQ2(NW)*TS
               IF ((NS-1).GT.0.AND.(DELTAT/60.0-TS).GT.0.0) THEN
                  DPPA(NQ,NS-1,5)=DPPB(NQ,NS-1,5)+DR*QQ2(NW)*
     .                            (DELTAT/60.0-TS)
               ENDIF
               HC=BI(NW)*HK(NW)*O(NW)/(2.0*A(NW))
               HLW=1050.0
               DRR=DR
C
C     TB: AIR TEMP. BEFORE A DATA POINT ADVANCED. F.
C     TQ: AIR TEMP. AFTER THE DATA POINT ADVANCED. F.
C     TTSS: MEAN AIR TEMP. IN THE SEGMENT COVERED BY THE DATA POINT
C           WITHIN THE PRESENT TIME INTERVAL.
C     TXM: MEAN WALL TEMP. IN THE SEGMENT. F.
C
               JLE=NSP(NW)
               TAUXA(NQ,NS)=TJS(NW)
               TB=TJS(NW)
               TQ=TJS(NW)-5.0
               TTSS=0.7*TB
               IF (TB.LT.200.0) TTSS=0.9*TB
               SEG=QQ2(NW)*TS*(TTSS+460)/(A(NW)*(TR+460))
25             TXM=DCOAGE(NW)*(TROCK(NW)-TTSS)/BI(NW)+TTSS
               CALL MTCOEF (TXM,TTSS,NW,BETA)
               FE=SEG/LA(NW)
C
C     WMTR: MASS TRANSFER BETWEEN WALL FILM AND BULK FLOW AT
C           STATIONS IN THE PRESENT TIME INTERVAL, LBM/FT2.
C     GFF: VARIATION IN HUMIDITY RATIO CAUSED BY VAPOR CONDENSATION
C          AS MIST IN THE SEGMENT IN THE PRESENT INTERVAL.
C
               WC=0.0
               WD=0.0
               NC=1
               DO 45 JL=JLE,1,-1
                  SDL=LA(NW)*(JL-1)/JLE
                  IF (SEG.GE.SDL) THEN
                     WC=WC+WMTR(NQ,JL)
                     WD=WD+GFF(NQ,JL)
                     NC=NC+1
                  ENDIF
45             CONTINUE
               GADD=WC/NC
               GF=WD/NC
               DRR=(DR+DR*GADD*18/28.98)*(TR+460)/((TTSS+460)*(1+GADD))
               GDP=GADD*DRR
               IF (GDP.LT.0.0) GDP=0.0
               CNT(NW)=HC*WTCORR(NQ,1)+HK(NW)*DCOAGE(NW)*TROCK(NW)
     .                 *O(NW)/(2.0*A(NW))+BETA*GDP*CP*(TXM-TQ)-
     .                 HLW*QQ2(NW)*DR*GF/(O(NW)*SEG)
               FAW(NW)=HK(NW)*DCOAGE(NW)*O(NW)*O(NW)/(2*QQ2(NW)*A(NW)*
     .                 60.0*CP*DR)
               SEG1=SEG
               CALL COLUMN (TB,SEG,NW,NQ,NS)
               EL=FAW(NW)*SEG
               IF (EL.LT.25.0) THEN
                  EXPN(NW)=EXP(-EL)
               ELSE
                  EXPN(NW)=0.0
               ENDIF
C
C     AIR TEMP. OF THE DATA POINT AT ITS NEWLY ARRIVED POSITION.
C
               DPPA(NQ,NS,2)=2*A(NW)*(CNT(NW)-(CNT(NW)-HK(NW)*DCOAGE(NW)
     .                 *O(NW)*TB/(2*A(NW)))*EXPN(NW))/(HK(NW)*DCOAGE(NW)
     .                 *O(NW))-DZRD(NW)*SEG/(LA(NW)*778.26*CP)

C
C     ITERATION PROCESS IS CALLED FOR (REFER TO NOTES FOR DRY CASE).
C
               NN=NN+1
               IF (ABS(SEG1-SEG).GT.0.1) THEN
                  B=2*A(NW)/(HK(NW)*DCOAGE(NW)*O(NW)*EL)
                  E=B*(CNT(NW)-HK(NW)*DCOAGE(NW)*TB*O(NW)/(2*A(NW)))
     .              *(1.0-EXPN(NW))
                  TTSS=2*A(NW)*CNT(NW)/(HK(NW)*DCOAGE(NW)*O(NW))-E-
     .                 SEG*DZRD(NW)/(2*778.26*LA(NW)*CP)
                  GADD=WMTR(NQ,1)*SEG*O(NW)*60.0/(QQ2(NW)*DR*TS)
                  TQ=DPPA(NQ,NS,2)
C
C     GADD: CHANGE IN HUMIDITY RATIO DUE TO MASS TRANSFER BETWEEN THE
C           WALL FILM AND BULK FLOW.
C     GF: CHANGE OF HUMIDITY RATIO OF C.V. BY CONDENSATION AS FOG.
C     ** ALL IN LB OF VAPOR/LB OF DRY AIR.
C
                  IF (NN.LE.20) GO TO 25
                  WRITE (8,*) 'SUB. CDJN NN EXCEEDED.',SEG,SEG1,NW,NS
               ENDIF
            ENDIF
C
C     DATA RECORDED IN DPPA: (1) LOCATION, (2) TEMP. (3) CH4,
C                            (4) FUME, (5) AIR MASS IN C.V..
C
            FA=XNEW(NW)*SEG/CP
            FAS=0.0
            IF (FA.LT.25.0) FAS=EXP(-FA)
            FE=SEG/LA(NW)
            IF (NTAUXC.EQ.1) TAUXC(NQ)=TROCK(NW)+(TTJS(NW)-TROCK(NW))*
     .                       FAS-(DZRD(NW)*FE/(2.0*778.26*CP))*(1.+FAS)
            DPPA(NQ,NS,1)=SEG
            DPPA(NQ,NS,3)=PROP(K)*(QQ2(NW)-CONTQ)/QQ2(NW)+CONTAM/QQ2(NW)
            DPPA(NQ,NS,4)=(PRCH4(K)*QQ2(NW)+CH4V(NW)*FE)/
     .                    (QQ2(NW)+CH4V(NW)*FE)
50       CONTINUE
60    CONTINUE
C
      RETURN
100   FORMAT (/,' * ERR * ABNORMAL TEMPERATURE IN HEAT SOURCE ',I5,/,
     .       ' CURRENT AIRFLOW:',F8.0,'  TEMPERATURE: ',F8.1,/,
     .       ' "HEAT" AND/OR "QCENT" MUST BE ADJUSTED.')
110   FORMAT (/,' * ERR * ABNORMAL TEMPERATURE IN HEAT SOURCE ',I5,/,1X,
     .       'CURRENT AIRFLOW:',F8.0,'. ADJUST "HEAT" AND/OR "QCENT".')
      END
