* typed in from Appendix C of Samuel R. Pate, AEDC-TR-77-107, March 1978
*  typed in by Carl Durkes, AAE secy asst, 1 April 2003.  
*  Edited S.P. Schneider, April 2003, to debug.
*  changed all ARSIN to ASIN, 2 Apr 2003.
*  change input format to list-directed, sps 4-2-2003
*       also wrote output to filename
*       also remove looping back for repeated input data
C    **********************************************************************
C    *********PREDICTION OF BOUNDARY LAYER TRANSITION *********************
C    *********ON SHARP FLAT PLATES AND SHARP CONES USING*******************
C    *********PATE’S AERODYNAMIC-NOISE-TRANSITION CORRELATION**************
C    **********************************************************************
C    
C    **********************GENERAL COMMENTS *******************************
C            **THIS COMPUTER CODE ALLOWS REASONABLY ACCURATE
C              PREDICTIONS OF TRANSITION REYNOLDS NUMBERS AND
C              LOCATIONS TO BE MADE ON SHARP FLAT PLATES AND
C              SLENDER CONES AT ZERO INCIDENCE FOR ALL SIZE
C              WIND TUNNELS AND 3 <M <20. 
C
C            **THIS COMPUTER CODE IS VALID FOR AIR OR NITROGEN
C
C            **THE AERODYNAMICS-NOISE-TRANSITION CORRELATION AND  
C              CONSEQUENTLY, THIS COMPUTER PROGRAM IS APPLICABLE
C              ONLY TO CONVENTIONAL SUPERSONIC-HYPERSONIC WIND
C	         TUNNELS HAVING TURBULENT BOUNDARY LAYERS ON THE
C              NOZZLE WALLS AND 3 <M <20.
C
C            **THE PREDICTED LOCATION OF TRANSITION CORRESPONDS TO 
C              THE END OF TRANSITION AS DEFINED BY THE PEAK IN 
C              A SURFACE PITOT PROBE PRESSURE TRACE
C
C            **THE VALUE OF THE TUNNEL WALL MEAN TURBULENT SKIN  
C              FRICTION COEFFICIENT USED IN THE AERODYNAMICS-NOISE	   	
C              TRANSITION CORRELATION IS COMPUTED USING THE METHOD
C              OF VAN DRIEST-II, INCLUDING NON-ADIABATIC WALL EFFECTS
C
C            **THE TUNNEL TEST SECTION WALL TURBULENT BOUNDARY
C              LAYER DISPLACMENT THICKNESS USED IN THE AERODYNAMIC-
C              NOISE-TRANSITION CORRELATION IS COMPUTED USING
C              CORRELATION DEVELOPED FROM 2-D AND 3-D NOZZLE DATA
C
C            **FOUR PROGRAM OPTIONS ARE AVAILABLE  
C             *KOPT = 1 AND 2 ARE FOR IDEAL GASES,GAMMA=1,4
C             *KOPT = 3 AND 4 ARE FOR REAL GAS NOZZLE EXPANSIONS
C
C    *********************************************************************
C                   PROGRAM OPTIONS
C    *********************************************************************
C
C            KOPT = 1 
C                **TRANSITION REYNOLDS NUMBERS AND LOCATIONS ON SHARP        
C                  FLAT PLATE OR HOLLOW CYLINDER AT ZERO INCIDENCE
C                **IDEAL GAS, GAMMA = 1.4   
C                **3 < RMACH < 10 
C                **** REQUIRED INPUT DATA ****
C                   * KOPT = 1
C                   * KGEOM = 1 OR 2
C                   * C=TUNNEL TEST SECTION CIRCUMFERENCE
C                   * XL= AXIAL DISTANCE FROM TUNNEL THROAT
C                         TO MODEL LEADING EDGE LOCATIONS,INCHES
C                   * P0= TUNNEL STILLING CHAMBER PRESSURE,PSIA
C                   * T0= TUNNEL STILLING CHAMBER TEMPERATURE ,DEG R
C                   * RMACH= TUNNEL TEST SECTION MACH NUMBER
C                   * TW=TUNNEL WALL TEMPERATURE, DEG R
C
C            KOPT = 2
C                **TRANSITION REYNOLDS NUMBERS AND LOCATIONS ON
C                  SHARP SLENDER CONES AT ZERO INCIDENCE 
C                **IDEAL GAS, GAMMA =1.4
C                **3 < RMACH < 10 
C                ****REQUIRED INPUT DATA ****
C                   *KOPT = 2
C                   *KGEOM = 1 OR 2 
C                   *C=TUNNEL TEST SECTION CIRCUMFERENCE 
C                   *XL= AXIAL DISTANCE FROM TUNNEL THROAT
C                        TO MODEL LEADING EDGE LOCATIONS,INCHES
C                   *P0= TUNNEL STILLING CHAMBER PRESSURE,PSIA
C                   *T0= TUNNEL STILLING CHAMBER TEMPERATURE ,DEG R    
C                   *RMACH= TUNNEL TEST SECTION MACH NUMBER
C                   *DELTA = CONE HALF ANGLE ,DEG
C                   *TW=TUNNEL WALL TEMPERATURE, DEG R 
C
C            KOPT = 3
C                ** TRANSITION REYNOLDS NUMBERS AND LOCATIONS
C                   ON SHARP FLAT PLATES AT ZERO INCIDENCE
C                ** REAL GAS EFFECTS CONSIDERED IN NOZZLE EXPANSION
C                   PROCESS 
C                ** 7 < RMACH < 20 , P0 > 5000 PSIA
C                ** FREE-STREAM STATIC TEMP, SET EQUAL TO 90 DEG R
C                **** REQUIRED INPUT DATA ****
C                   * KOPT = 3
C                   * KGEOM = 3 
C                   * C = TUNNEL TEST SECTION CIRCUMFERENCE
C                   * XL= AXIAL DISTANCE FROM TUNNEL THROAT TO
C                        MODEL LEADING EDGE LOCATIONS, INCHES
C                   * RMACH
C                   * REINF = TUNNEL TEST SECTION UNIT REYNOLDS
C                             NUMBER ( OBTAIN FROM FIG      )
C                   * TW = TUNNEL WALL TEMPERATURE, DEG R
C
C            KOPT = 4
C                ** TRANSITION REYNOLDS NUMBERS AND LOCATIONS
C                   ON SHARP SLENDER CONES AT ZERO INCIDENCE
C                ** REAL GAS EFFECTS CONSIDERED IN NOZZLE
C                   EXPANSION PROCESS
C                ** 7 < RMACH < 20 , PP > 5000 PSIA
C                ** FREE-STREAM STATIC TEMPERATURE SET EQUAL
C                   TO 90 DEG R
C                ****REQUIRED INPUT DATA ****
C                   *KOPT = 4
C                   *KGEOM = 3
C                   *C = TUNNEL TEST SECTION CIRCUMFERENCE  
C                   *XL =AXIAL DISTANCE FROM TUNNEL TO THROAT TO 
C                        MODEL LEADING EDGE LOCATION, INCHES  
C                   *RMACH
C                   *REINF = TUNNEL TEST SECTION UNIT REYNOLDS
C                            NUMBER ( OBTAIN FROM FIG      ) 
C                   *RRESC = RATIO OF CONE SURFACE UNIT REYNOLDS
C                            NUMBER TO FREE-STREAM VALUE
C                            (USE FIG        )
C                   *TW = TUNNEL WALL TEMPERATURE
C                   *DELTA= CONE HALF ANGLE DEG
C    ********************************************************************* 
C        TUNNEL NOZZLE GEOMETRIES
C          KGEOM = 1
C                * FOR TWO-DIMENSIONAL CONTOURED NOZZLES
C                * 1.5 < RMACH <5 ADIABATIC WALLS
C                *   5 < RMACH * 10 NON-ADIABATIC WALLS
C                * METHOD OF MAXWELL IS USED TO COMPUTE
C                  BOUNDARY LAYER DISPLACEMENT THICKNESS ON 
C                  TUNNEL TEST SECTION WALL
C                * IDEAL GAS ,GAMMA = 1.4
C                * USE WITH KOPT = 1 OR 2
C          KGEOM = 2 
C                * FOR AXISYMMETRIC CONTOURED NOZZLES
C                * 1.5< RMACH<5  ADIABATIC WALLS
C                *  5 < RMACH * 10 NON-ADIABATIC WALLS
C                * NON-ADIABATIC WALL
C
C          KGEOM = 3 
C                * FOR CONCIAL AND CONTOURED AXISYMMETRIC NOZZLES
C                * NON-ADIABATIC WALL
C                * 7 < RMACH < 20
C                * USE WITH KOPT = 3 OR 4 
C                * METHOD OF WHITFIEL1,1S USED TO COMPUTE BOUNDARY
C                  LAYER DISPLACEMENT THICKNESS ON TUNNERL TEST
C                  SECTION WALL 
C     ********************************************************************
C
      character*20 outfil
*
*
      write (*,*) 'enter KOPT,KGEOM: '
    1 READ (5,*,END=8000)  KOPT, KGEOM
c   2 FORMAT(I1,I1)
      write (*,*) 'enter filename for output: '
	read (*,2) outfil
    2 format(a)
      open (unit=6,file=outfil)
    3 IF(KOPT-2)  4, 100, 200
C     ********************************************************************
C     KOPT = 1
C     WHITE COLUMN HEADLINES
    4 WRITE (6,12) KOPT , KGEOM   
   12 FORMAT ( '1' , ' ESTIMATION OF BOUNDARY-LAYER TRANSITION ON SHARP 
     1FLAT PLATES' //
     2,2X,'IN CONVENTIONAL SUPERSONIC-HYPERSONIC WIND TUNNELS USING'//
     3,1X, ' PATES AERODYNAMIC NOISE TRANSITION CORRELATION ' //
     4,2X,'KOPT=',I2,3X,'KGEOM=',I2//
     5,4X,'C,IN.',2X'XL,IN.'2X,'RMACH',2X,'P0,PSIA',3X'T0,R',4X,'TW,R
     6',3X,'REINF/FT',5X,'REL',8X,'CF',5X,'DELS,IN.',3X,' RETFP',4X,
     7'XTFP,IN.',2X,'TW/TAW')
	write (*,*) 'enter C(in.),XL(in.),P0 (psia),T0 (R),RMACH,TW (R): '
   14 READ (5,*,END= 1 )   C,XL,P0,T0,RMACH,TW
c  16 FORMAT ( 6F10.0)
   18 C1 = 48.0
      GAMMA = 1.4
      R = 1716.0
C     COMPUTE FREE-STREAM TEMPERATURE
      T1=T0/(1.0+0.2*(RMACH**2.0))
C     COMPUTE FREE-STREAM VELOCITY
      U1=RMACH*SQRT(GAMMA*(R*T1))
C     COMPUTE FREE-STREAM ABSOLUTE VISCOSITY 
C     USE SUTHERLANDS VISCOSITY LAW (T1>216 DEG R)
      IF(T1 .LE. 216) GO TO 30
      VIS1=2.270*(T1**1.5)*(10.0**(-8.0))/(198.6+T1)
      GO TO 32
C     USE LINEAR VISCOSITY LAW (T1<216 DEG R) 
   30 VIS1=(0.0805*T1)*(10.0**(-8.0))
   32 CONTINUE
C     COMPUTE FREE-STREAM STATIC PRESSURE
      P1=P0/((1.0+0.2*(RMACH**2.0))**3.5)
C     COMPUTE FREE-STREAM UNIT REYNOLDS NUMBER
      REINF=((144.0*P1)*U1)/((R*T1)*VIS1)
C     COMPUTE REYNOLDS NUMBER BASED ON TUNNEL NOZZLE LENGTH (XL)
   34 REL = REINF *(XL/12.0)
      IF ( KGEOM-3) 35,90,8000
   35 CONTINUE
C     COMPUTE TUNNEL TEST SECTION WALL TURBULENT BOUNDARY LAYER
C     DISPLACMENT THICKNESS COMPUTED USING MAXWELL'S CORRELATION
C     ABSOLUTE VISCOSITY (VIS) COMPUTED USING SUTHERLAND'S LAW
C     VALID FOR  216<T< 5000 ,DEG R
   36 VIS0= (2.270)*(T0**1.5)*(10.0**(-8.0))/(198.6+T0)
C     RH0 IS STILLING CHAMBER DENSITY
   38 RH0 =(144 * P0)/(R* T0)
C     A0 IS STILLING CHAMBER SPEED OF SOUND
   40 A0 = SQRT((GAMMA*R)* T0)
C     OK IS THE VARIABLE IN MAXWELL'S EQ. 
   42 OK = (0.0131)*((VIS0/(RH0* A0))**(0.14286))
C     BARDEL IS MAXWELL CORRELATION PARAMETER
C    1BARDEL IS COMPUTED USING A THIRD DEGREE POLYNOMIAL CURVE FIT OF
C    1MAXWELL'S ORIGINAL CORRELATION FOR RMACH   .LE. 6.0
      IF (RMACH - 5.0 ) 46, 46, 50
   46 IF (KGEOM .GT. 1) GO TO 8000
   48 BARDEL =   1.140824 * RMACH  + 0.088132*((RMACH)**2.0)
     1+ 0.026978 *((RMACH)**3.0)
      GO TO 56
   50 IF ( RMACH.GT. 10.0 ) GO TO 8000
      IF (KGEOM -2 ) 52, 54, 8000
C     IF 5< RMACH<10 AND KGEOM = 1 THEN BARDEL = 2.0 + 1.8333*RMACH
   52 BARDEL  = 2.0 +  1.8333 * RMACH
      GO TO 56
   54 BARDEL = 0.167 + 1.833*RMACH
C     DELS  IS TUNNEL WALL DISPLACEMENT THICKNESS IN IN.
   56 DELS=(OK)*(BARDEL)*((XL/12.0)**0.857143)*12.0
C     COMPUTE MEAN TURBULENT SKIN FRICTION(CF) USING VAN DRIEST-II WITH
C     TW IS INPUT DATA , DEG. RANKINE
C     RT IS THE RECOVERY FACTOR FOR A TURBULENT BOUNDARY LAYER
   57 RT= 0.90
   58 A=SQRT((((GAMMA-1.0)/2.0)*RT*(RMACH**2.0))*(T1/TW))
      B = (T1/TW) + (A**2.0)-1.0
      B1=B**2
      B=SQRT(B1)
      ALPHA=(2.0*(A**2.0)-B)/SQRT(B**2.0+4.0*(A**2.0))
      BETA = B/SQRT ((B**2.0)+4.0*(A**2.0))
      C4 = (0.242)*(asin(ALPHA) + asin(BETA))/(A * SQRT(TW/T1))
      C2 = ALOG10(REL)
      W = 0.76
      C3= 1.5*ALOG10(T1/TW) + ALOG10((198.6+TW)/(198.6+T1))
C     STARTING VALUE OF CF IS 0.0050/RMACH
      CF= (0.0050)/RMACH
      K=0
   70 FCF = SQRT(CF)*(C2+C3  + ALOG10(CF)) - C4
      K=K+1
      IF (K .GT. 100) GO TO 8000
C     FPRIME IS THE DERIVATIVE OF FCF WITH RESPECT TO CF
      FPRIME = (0.50 / SQRT(CF))*( C2+C3 + ALOG10(CF)+0.8686)
C     THE NEWTON-RAPHSON METHOD IS 
      CF = CF - (FCF/FPRIME)
C     INTERATE UNTIL SUCCESSIVE APPORXIMATIONS OF ROOT IS LESS THAN
C     10.0000010 
      ROOT = ABS( FCF/FPRIME)
      IF (ROOT .GT. 0.0000010) GO TO 70
      C1=48.0
      IF (KOPT .EQ. 2) GO TO 152
      IF (KOPT .EQ. 4) GO TO 222
      AFP=(CF)**(-2.55)
      RC1C = C1/C
      IF (RC1C - 1.0 ) 71, 71, 72
   71 BFP = 0.56 + (0.44)* C1/C
      GO TO 73
   72 BFP = 1.00
   73 CFP = SQRT (DELS/ C ) 
      DFP=(0.0126)*AFP
      RETFP=(BFP*DFP)/CFP
      XTFP=(RETFP/REINF)*12.0
C     TAW= TUNNEL WALL RECOVERY TEMPERATURE
      TAW= T1*(1.0+ 0.18*(RMACH**2.0))
      RTWTAW= TW/TAW
      IF ( KOPT .EQ. 3) GO TO 210
   80 WRITE (6,82)  C,XL,RMACH,P0,T0,TW,REINF,REL,CF,DELS,RETFP,XTFP,
     1RTWTAW
   82 FORMAT (1x,F8.1,F7.1,F8.1,F9.2,2F8.1,1PE11.4,1PE11.4,0PF10.6,
     1F9.4,1PE11.4,0PF8.2,F10.3) 
	stop 'normal end of KOPT=1' !added sps, replace looking for new input data
c     GO TO 14
C     COMPUTE DELS USING WHITFIELDS CORRELATION FOR CONICAL NOZZLES OR
C     REAL GAS AXISYMMETRIC NOZZLES
C     WHEN 7< RMACH <20
   90 DELS=(0.22)*SQRT(RMACH)/SQRT(REL)
      GO TO 58
C     ********************************************************************
C     KOPT = 2 
  100 WRITE (6,101) KOPT,KGEOM
  101 FORMAT ( '1', 'ESTIMATION OF BOUNDARY LAYER TRANSITION ON SHARP SL
     1ENDER CONES IN CONVENTIONAL SUPERSONIC-HYPERSONIC',//
     1,1X,'WIND TUNNELS USING THE AERODYNAMIC NOISE CORRELATIONS BY PATE
     1',//
     1,1X,'KOPT =',I2,3X,'KGEOM =',I2,//
     1,2X,'C,IN.',2X,'XL,IN.',2X,'RMACH',3X,'P0,PSIA',1X,'T0,R',4X
     2'TW,R',3X,'REINF/FT',3X,'DELT,DEG',1X,'RTS',4X,'RPS',3X,'RMACHS'
     3,1X,'RRESC',4X,'RETSC',3X,' XTSC,IN',2X,'TW/TAW')
	write (*,*) 'enter C,XL,P0,T0,RMACH,DELTA,TW: '
  102 READ(5,*,END= 1 )  C,XL, P0, T0, RMACH, DELTA, TW
c 104 FORMAT ( 7F10.0)
C     COMPUTE CONE SURFACE STATIC PRESSURE RATIO (RPSI = PS/P1) USING
      DELTA = DELTA/ 57.296
C    1RASMUSSEN'S .EQ.
      PS1 = (SIN(DELTA)**2.0)*(1.4*(RMACH**2.0))/2.0
      PS2 = 2.4*(RMACH**2.0)*(SIN(DELTA)**2.0) + 2.0
      PS3 = 0.4*(RMACH**2.0)*(SIN(DELTA)**2.0) + 2.0
      PS4 = ALOG(1.20 + 1.0/((RMACH**2.0)*(SIN(DELTA)**2.0)))
      PS5 = (PS2/PS3)*PS4
      RPS1 = 1.0 +PS1*(1.0+PS5)
      GAMMA = 1.40
      R = 1716.0
C     COMPUTE CONE BOW SHOCK ANGLE(PHI)
  110 PHI = asin( SIN(DELTA)*((1.20 + 1.0/((RMACH*SIN(DELTA))**2.0))**
     10.50))
c     write (*,*) 'phi= ',phi  !debug
C     COMPUTE STATIC PRESSURE AND TEMPERATURE BEHIND BOW SHOCK USING
C     2-D OBLIQUE SHOCK WAVE THEORY
C     RT21 = T2/T1
  112 RT21 = (7.0*((RMACH * SIN(PHI))**2.0)-1.0)*(((RMACH*SIN(PHI))**2.0
     1)+ 5.0)/(36.0*((RMACH*SIN(PHI))**2.0))
C     RP21 = P2/P1 
  114 RP21 = (7.0*((RMACH * SIN(PHI))**2.0)-1.0 ) /6.0
C     COMPUTE PRESSURE AND TEMPERATURE ON CONE SURFACE BY USING CONDITION
C     BEHIND SHOCK WAVE AND ISENTROPHIC COMPRESSION PROCESS
C     RTS1 = TS/T1
  116 RTS1 = RT21 *( (RPS1 / RP21)**0.285714)
      if (rts1 .le. 0.0) then !sps add error trapping
	  write (*,*) 'rts1= ',rts1,'rt21= ',rt21
	  stop 'fatal error'
	end if
C     COMPUTE FREE-STREAM TEMPERATURE (T1)
  118 T1= T0/(1.0 + 0.2*(RMACH)**2.0)
  120 TS= T1*RTS1
C     COMPUTE VELOCITY AT CONE SURFACE (US)
C     RATIO OF SPECIFIC HEAT IS CP=6006 FT.SQ./ SEC. SQ.-DEGR
      CP = 6006
  122 US =(SQRT(2.0 * CP))*SQRT( T0 - (T1*RT21)*((RPS1/RP21)**0.285714))
C     COMPUTE MACH NUMBER AT CONE SURFACE (RMACHS)
  124 RMACHS = US /(49.0 * SQRT(TS))
C     COMPUTE DENSITY AT CONE SURFACE (RHOS)
      R=1716.0
C     COMPUTE FREE-STREEM STATIC PRESSURE (P1)
  126 P1 = P0/ ((1.0 + 0.2*(RMACH**2.0))**3.5)
C     COMPUTE DENSITY RATIO AT CONE SURFACE (RRHOS1 = RHOS/RHO1)
  128 RRHOS1=(RPS1/RTS1)
C     COMPUTE REYNOLDS NUMBER AT CONE SURFACE (RRES1 = RES/RE1)
C     COMPUTE FREE-STREAM VELOCITY (U1)
  130 U1 = (RMACH * 49.00)* SQRT(T1)
  132 IF ( T1 .GE. 216) GO TO 138
C     COMPUTE FREE-STREAM ABSOLUTE VISCOSITY (VIS1) USING LINEAR LAW
      VIS1 = (0.0805 * T1)*(10.0)**(-8.0)
  136 GO TO 140
C     COMPUTE FREE-STREAM ABSOLUTE VISCOSITY USING SUTHERLANDS LAW
  138 VIS1 = (2.270*(T1**1.50))*(10.0**(-8.0))/(198.6 + T1)
C     RE1 = REINF = RHO1 * U1/VIS1
C     COMPUTE FREE-STREAM-UNIT REYNOLDS NUMBER( RE1=REINF) 
  140 REINF=(U1*P1*144.0) / (R*T1 *VIS1)
  142 IF (TS .GE. 216) GO TO 146
C     COMPUTE CONE SURFACE VISCOSITY VISS USING LINEAR LAW (TS) 216 DEGR
  144 VISS = (0.0805 * TS)*(10.0**(-8.0))
      GO TO 148
C     COMPUTE CONE SURFACE VISCOSITY USING SOTHERLANDS LAW (TS>216 DEGR)
  146 vISS = (2.270*(TS**1.50))*(10.0**(-8.0))/(198.6 + TS)
C     COMPUTE CONE SURFACE REYNOLDS NUMBER (RES)
  148 RES=(RPS1*P1*144.0*US)/(R*TS*VISS)
C     COMPUTE REYNOLDS NUMBER RATIO(RRES1 = RES/REINF)
  149 RRES1 = RES/REINF
  150 GO TO 34
  152 CONTINUE
      DELTA = DELTA * 57.296
C     COMPUTE TRANSITION REYNOLDS NUMBER ON CONE SURFACE (RETSC)
  153 AC = (CF)**(-1.40)
      RC1C = C1/ C
      IF (RC1C -1.0) 154,154,155
  154 BC = 0.80 + 0.20*(C1/C)
      GO TO 156
  155 BC = 1.00
  156 CC = SQRT (DELS/C)
      RETSC=(48.5)*(AC * BC ) / CC
      XTSC = (RETSC / RES) * 12.0
C     TAW= TUNNEL WALL RECOVERY TEMPERATURE
      TAW= T1*(1.0+ 0.18*(RMACH**2.0))
      RTWTAW= TW/TAW
      WRITE (6,158) C,XL,RMACH,P0 ,T0,TW  ,REINF,DELTA,RTS1,RPS1,
     2RMACHS,RRES1,RETSC,XTSC,RTWTAW
  158 FORMAT ( 1x,F6.1,F7.1,F8.1,F9.1,F7.1,F7.1,1PE12.4,0PF8.2
     2,F8.3, F7.3,F7.2,F7.3,1PE12.4,0PF7.2,F8.3)
	stop !replace statement below where look for new data
c     GO TO 102
C     ****************************************************************
C     KOPT = 3
C
  200 IF(KOPT-4) 201,214,8000
  201 IF (KGEOM -3 ) 8000 ,202 , 8000
  202 WRITE ( 6, 204) KOPT, KGEOM
  204 FORMAT ( '1', 'PREDICTION OF BOUNDARY LAYER TRANSITION ON SHARP
     1FLAT PLATES'//
     1,2X,'FOR NON IDEAL GAS NOZZLE EXPANSION PROCES'//
     1,2X,'RMACH > 10 OR RMACH > 8 , REINF > 15*10**6 ,SEE FIG    '//
     1,2X,'USING PATES   AERODYNAMICS NOISE TRANSITION CORRELATION'//
     1,2X, 'KOPT = ',I2,4X,'KGEOM= ',I2,//
     1,4X, 'C,IN.',5X, 'XL,IN.',4X, 'RMACH',5X, 'TW',5X,'REINF',8X, 'CF'
     1,3X, 'DELS,IN.'8X, 'RETFP',5X, 'XTFP,IN.'5X,'TW/TAW')
	write (*,*) 'enter C,XL,RMACH,REINF,TW: '
  206 READ(5, *, END=1 )  C, XL, RMACH, REINF, TW
c 208 FORMAT (5F10.0)
      REL = REINF*(XL/12.0)
C     COMPUTE DELS USING WHITFIELDS FORMULAS 
      DELS = XL*(0.22)*SQRT(RMACH)/(REL**0.25)
c     COMPUTE TUNNEL WALL SKIN FRICTION USING VAN DRIEST-II
      GAMMA=1.40
C     SET FREE-STREAM STATIN TEMPERATURE EQUAL TO 90 DEG R
      T1 = 90.0
C     TAW= TUNNEL WALL RECOVERY TEMPERATURE 
      TAW= T1*(1.0+ 0.18*(RMACH**2.0))
      RTWTAW=TW/TAW
      GO TO 57
  210 CONTINUE
      WRITE (6,211) C,XL,RMACH,TW,REINF,CF,DELS,RETFP,XTFP,RTWTAW
  211 FORMAT (1x,1F8.1,1F10.1,1F9.1,1F8.1,1PE12.4,0PF10.6,1F8.4,1PE16.4,
     10PF10.2,F8.2)
	stop !replace statement below where look for new input data
c 212 GO TO 206
C    ******************************************************************
C     KOPT = 4 
C     
  214 IF (KGEOM - 3 ) 8000 , 215 , 8000
  215 WRITE ( 6, 216) KOPT, KGEOM
  216 FORMAT ('1', 'ESTIMATION OF BOUNDARY LAYER TRANSITION ON SHARP
     1SLENDER CONES'//
     1,2X, 'FOR MACH NUMBERS GREATER THAN 10 ON VERY HIGH REYNOLDS'//
     1,2X,'AT MACH NUMBER EQUAL APPROX 8(REINF>10**7)'//
     1,2X,'USING PATES  AERODYNAMICS NOISE TRANSITION CORRELATION'//
     1,2X, 'KOPT = ',I2,4X, 'KGEOM= ',I2,//
     1,1X,'C,IN.',5X,'XL,IN.',4X,'RMACH',3X,'TW,DEGR',4X,'REINF/FT'
     2,3X,'RESC/REINF',2X,'DELS,IN.'7X,'RETSC',5X,'XTSC,IN.',3X,
     3'DELTA,DEG',5X,'TW/TAW' )
	write (*,*) 'enter C,XL,RMACH,REINF,RRESC,TW,DELTA: '
  218 READ(5,*,END= 1 ) C,XL,RMACH,REINF,RRESC,TW,DELTA
c 220 FORMAT (7F10.0)
      REL = REINF*(XL/12.0)
C     COMPUTE DELS USING WHITFIELD'S FORMULAS
      DELS = XL*(0.22)*SQRT(RMACH)/(REL**0.25)
C     COMPUTE TUNNEL WALL SKIN FRICTION USING VAN DRIEST-II
      GAMMA=1.40
C     SET FREE-STREAM STATIC TEMPERATURE EQUAL TO 90 DEG R
      T1 = 90.0 
C     TAW= TUNNEL WALL RECOVERY TEMPERATURE
      TAW= T1*(1.0+ 0.18*(RMACH**2.0))
      RTWTAW= TW/TAW
      GO TO 57
  222 CONTINUE
C     COMPUTE TRANSITION REYNOLDS NUMBER ON CONE
      AC = (CF)**(-1.40)
      RC1C = C1 /C
      IF(RC1C - 1.0 ) 224, 224, 225
  224 BC = 0.80 + 0.20 *(C1/C ) 
      GO TO 226
  225 BC = 1.00
  226 CC = SQRT (DELS/ C )
      RETSC =(48.5)*(AC * BC ) /CC
      RESC= RRESC*REINF
      XTSC= (RETSC/RESC)*12.0
      WRITE(6,230)C,XL,RMACH,TW,REINF,RRESC,DELS,RETSC,XTSC,DELTA,RTWTAW
  230 FORMAT(1x,F5.1,2F10.2,F10.2,1PE13.4,0PF9.3,F12.3,1PE16.4,0PF10.3,
     1F10.3,F13.2 )
c     GO TO 218 
      stop 'normal end'  !replace statement above where get new input data
* following is an error end
 8000 CONTINUE
      STOP 'abnormal termination 8000, check input' !text added sps 
      END 

























                 