DECLARE FUNCTION SATVP! (T!) DECLARE SUB CLEARSKY (S!) REM =================================== REM EQUILIBRIUM TEMPERATURE CALCULATION REM =================================== REM REM REF: P. SHANAHAN, "WATER TEMPERATURE MODELING: A REM PRACTICAL GUIDE", PROCEEDINGS OF STORMWATER AND REM WATER QUALITY MODEL USERS GROUP MEETING, REM APRIL 12-13, 1984, REPORT EPA-600/9-85-003. REM DIM U(3), H(3) DIM SHARED DAY(12), C1(21), C2(21), C3(21) REM REM SET UP CONSTANTS FOR CLEAR-SKY RADIATION CALCULATION FOR I = 1 TO 12 READ DAY(I) NEXT I DATA 0,31,59,90,120,151,181,212,243,273,304,334 FOR I = 1 TO 21 READ C1(I), C2(I), C3(I) NEXT I DATA 80.155, 29.207, 1.679 DATA 79.371, 30.236, 1.713 DATA 78.566, 31.219, 1.710 DATA 77.604, 32.145, 1.740 DATA 76.655, 33.156, 1.728 DATA 76.041, 34.133, 1.694 DATA 75.060, 35.194, 1.737 DATA 74.046, 35.938, 1.734 DATA 73.161, 36.834, 1.727 DATA 72.248, 37.699, 1.738 DATA 71.390, 38.598, 1.721 DATA 70.394, 39.413, 1.730 DATA 69.350, 40.188, 1.741 DATA 68.362, 40.982, 1.739 DATA 67.281, 41.706, 1.742 DATA 66.240, 42.442, 1.736 DATA 65.197, 43.128, 1.740 DATA 64.113, 43.788, 1.739 DATA 63.010, 44.471, 1.739 DATA 61.911, 45.020, 1.740 DATA 60.782, 45.639, 1.735 REM REM FORMATTED OUTPUT FILE REM OPEN "THERMAL.OUT" FOR OUTPUT AS #1 REM A = 0 B = 0 REM 10 PRINT " INPUT T, W, R, S, C, F" PRINT " T IS THE AIR TEMPERATURE, DEG. F" PRINT " W IS THE WIND SPEED, MPH AT 2 M" PRINT " R IS THE RELATIVE HUMIDITY, PERCENT" PRINT " S IS THE OBSERVED SOLAR RADIATION, BTU/SQ FT/DAY" PRINT " (IF S IS ZERO, A VALUE WILL BE COMPUTED)" PRINT " C IS THE CLOUD COVER, FRACTION" PRINT " F IS THE WIND SPEED FUNCTION:" PRINT " F=1 FOR LAKE HEFNER" PRINT " F=2 FOR MODIFIED RYAN" PRINT " F=3 FOR MODIFIED MEYER" PRINT " F=4 FOR USER SPECIFIED" 50 PRINT " " PRINT " INPUT T, W, R, S, C, F" PRINT " ENTER ZEROS TO STOP" REM REM GET USER INPUT DATA REM INPUT T, W, R, S, C, F IF (T = 0 AND W = 0 AND R = 0 AND S = 0 AND C = 0 AND F = 0) THEN STOP CLOSE #1 END IF IF (W < 0 OR R <= 0 OR R > 100 OR S < 0 OR C < 0 OR C > 1 OR (F <> 1 AND F <> 2 AND F <> 3 AND F <> 4)) THEN PRINT " " PRINT " INPUT ERROR - PLEASE RE-ENTER DATA" PRINT " " GOTO 10 END IF IF (F = 4) THEN PRINT " INPUT A,B FOR WIND FUNCTION F(W) = A+B*W" INPUT A, B END IF REM REM COMPUTE METEOROLOGICAL CONSTANTS REM REM CONVERT UNITS R = R / 100 T4 = T + 460 REM COMPUTE SOLAR RADIATION IF (S = 0) THEN CALL CLEARSKY(S2) S = S2 * (1 - .65 * C ^ 2) END IF H3 = .94 * S REM COMPUTE AIR VAPOR PRESSURE AND VIRTUAL TEMPERATURE E1 = SATVP(T) * R T3 = T4 / (1 - .378 * E1 / 760) - 460 REM COMPUTE LONG WAVE (ATMOSPHERIC) RADIATION H2 = 2.1E-08 * (1 + .149 * SQR(E1)) * T4 ^ 4 * (1 + .17 * C ^ 2) REM REM USE BISECTION METHOD TO COMPUTE EQUILIBRIUM TEMPERATURE REM I1 = 1 I2 = 3 U(1) = T - 5 U(2) = T U(3) = T + 5 Z1 = 0 100 FOR I = I1 TO I2 Z1 = Z1 + 1 IF (Z1 > 100) THEN PRINT " " PRINT " FAILURE TO CONVERGE -- RUN ENDING" CLOSE #1 STOP END IF U1 = U(I) + 460 REM COMPUTE LONGWAVE BACK RADIATION H1 = 4E-08 * U1 ^ 4 E2 = SATVP(U(I)) IF (F = 1) THEN REM LAKE HEFNER EQUATION F1 = 17 * W ELSEIF (F = 2) THEN REM MODIFIED RYAN EQUATION D2 = U1 / (1 - .378 * E2 / 760) - 460 D2 = D2 - T3 D3 = .0024 * W ^ 3 IF (D2 < D3) THEN D2 = D3 F1 = 17.7 * D2 ^ .3334 + 11.1 * W ELSEIF (F = 3) THEN REM MODIFIED MEYER EQUATION F1 = 68 + 8.5 * W ELSE REM USER SPECIFIED WIND FUNCTION F1 = A + B * W END IF REM COMPUTE EVAPORATIVE HEAT FLUX H4 = F1 * (E2 - E1) IF (H4 < 0) THEN H4 = 0 REM COMPUTE HEAT CONDUCTION H5 = F1 * .255 * (U(I) - T) REM COMPUTE TOTAL HEAT FLUX H(I) = H2 + H3 - H1 - H4 - H5 REM DETERMINE IF H HAS CONVERGED TO 1 BTU/SQ FT/DAY OR LESS IF (ABS(H(I)) < 1) THEN 200 NEXT I REM PREPARE FOR NEXT ITERATION I1 = 2 I2 = 2 IF (H(1) < 0) THEN U(1) = U(1) - 10 H(1) = 0 I1 = 1 ELSEIF (H(3) > 0) THEN U(3) = U(3) + 10 H(3) = 0 I2 = 3 ELSE J1 = 1 IF ((H(1) * H(2)) < 0) THEN J1 = 3 U(J1) = U(2) H(J1) = H(2) END IF U(2) = .5 * (U(1) + U(3)) GOTO 100 REM REM EQUILIBRIUM TEMPERATURE FOUND REM 200 PRINT " " PRINT USING " EQUILIBRIUM TEMPERATURE = ###.# DEG F"; U(2) PRINT " " PRINT " " PRINT #1, USING "###.#"; U(2) GOTO 50 END REM ============================== REM CALCULATION OF CLEAR SKY SOLAR REM RADIATION AT WATER SURFACE REM IN UNITS OF BTU/SQ FT/DAY REM ============================== REM REM REF: THACKSTON, 1974, REM REPORT EPA-660/2-74-85 REM SUB CLEARSKY (S) REM REM INPUT LATITUDE AND DATE 1010 PRINT " " PRINT " ENTER SITE LATITUDE" PRINT " (26 TO 46 DEGREES)" INPUT L L1 = INT(L) L3 = L1 - 25 L2 = INT(L + .99) L4 = L2 - 25 IF (L1 < 26 OR L1 > 46) THEN PRINT "INVALID VALUE - PLEASE RE-ENTER" GOTO 1010 END IF IF (L2 < 26 OR L2 > 46) THEN PRINT "INVALID VALUE - PLEASE RE-ENTER" GOTO 1010 END IF PRINT " " PRINT " ENTER DATE AS MONTH,DAY" INPUT D7, D8 D = DAY(D7) + D8 REM REM COMPUTE CLEAR SKY SOLAR RADIATION REM D6 = 2 * 3.14159 * D / 366 S = C1(L3) - C2(L3) * SIN(C3(L3) + D6) IF (L1 <> L2) THEN S2 = C1(L4) - C2(L4) * SIN(C3(L4) + D6) L1 = L - L1 L2 = L2 - L S = L2 * S + L1 * S2 END IF S = 24 * S REM CORRECT FOR ABSORBTION S = S / .94 REM REM DISPLAY RESULT REM PRINT " " PRINT USING " AT LATITUDE ##.# ON ##/##"; L; D7; D8 PRINT USING " COMPUTED CLEAR SKY SOLAR RADIATION IS #### BTU/SQ FT/DAY"; S PRINT " " EXIT SUB END SUB REM ================================== REM SATURATION VAPOR PRESSURE FUNCTION REM ================================== REM FUNCTION SATVP (TS) SATVP = 25.4 * EXP(17.62 - 9500.8 / (TS + 460)) END FUNCTION