* This program was written by Gene Preston 08/06/2015. Please distribute and use as you wish. In using this program, * you agree that Gene Preston is not responsible for any consequences that may arise from the use of this program. * This program calculates the annual LOLE and ALLR using the direct approach. * The sequential Monte Carlo is also added to calculate the annual LOLP or ALOLP. * A conversational feature is added in which D1 D2 and D3 can be adjusted to achieve a target index. * * BEGIN THE DIRECT SOLUTION CHARACTER*50 A/' '/ REAL*8 F(0:1000)/1.D0,1000*0.D0/,LOLE,ALLR,EUE,ENGY,LOLH,FOR(10) & /.15D0,.10D0,.09D0,.07D0,.06D0,.05D0,.04D0,.04D0,.04D0,.0454D0/ INTEGER PMAX(10)/5,10,13,22,43,57,101,198,248,276/ REAL HLOAD(24,365),DAILY(24) &/.56,.60,.64,.68,.72,.76,.80,.84,.88,.92,.96,1.0, & 1.0,.96,.92,.88,.84,.80,.76,.72,.68,.64,.60,.56/ * ADDITIONAL INFORMATION FOR THE SEQUENTIAL MONTECARLO SOLUTION REAL MTTR(10)/5.,5.,5.,7.,7.,14.,21.,28.,28.,42./ ! MEAN TIME TO REPAIR IN DAYS REAL RR(10),FR(10) ! PROBABILITY OF HOURLY REPAIR AND FAILURE RATES INTEGER ISTATE(10) ! 0, GEN IS OFF 1, GEN IS ON REAL*8 YRS ! COUNT THE NUMBER OF YEARS SIMULATED * CALCULATE F(X) CUMULATIVE CAPACITY DISTRIBUTION CURVE USING THE DIRECT METHOD DO 1 I=1,10 DO 1 IX=1000,1,-1 J=IX-PMAX(I) IF(J.LT.0) J=0 1 F(IX)=(1.-FOR(I))*F(J)+FOR(I)*F(IX) * INITIALIZE CONSTANTS AND VARIABLES D1=260. ! 341 MIN LD DAYS PEAK DEMANDS D2=520. ! 10 WINTER DAYS PEAK DEMANDS D3=485. ! 14 SUMMER DAYS PEAK DEMANDS WRITE(*,*) 'THESE DEMANDS WERE SELECTED TO MAKE THE LOLE=0.1 D/Y' WRITE(*,*) 'D1 MIN, D2 WIN, D3 SUM =',IFIX(D1),IFIX(D2),IFIX(D3) WRITE(*,*) 'ENTER NEW D1,D2,D3 OR RETURN TO CONTINUE' READ(*,'(A50)') A 9 IF(A.NE.' ') READ(A,*) D1,D2,D3 LOLE=0.D0 ! LOSS OF LOAD EXPECTATION LOLH=0.D0 ! LOSS OF LOAD HOURS ALLR=1.D0 ! ANNUAL LOAD LOSS RISK EUE=0.D0 ! EXPECTED UNSERVED ENERGY ENGY=0.D0 ! TOTAL LOAD ENERGY * SWEEP THROUGH THE YEAR CALCULATING LOADS AND INDICES DO 3 J=1,365 DO 3 I=1,24 HLOAD(I,J)=IFIX(D1*DAILY(I)+.01) ! LIGHT LOAD PERIOD IF(J.LE.10) HLOAD(I,J)=IFIX(D2*DAILY(I)+.01) ! 10 WINTER PEAK DEMANDS IF(J.GE.201.AND.J.LE.214) HLOAD(I,J)=IFIX(D3*DAILY(I)+.01) ! 14 SUMMER PEAK DEMANDS IF(I.EQ.12) THEN ! THIS IS THE PEAK HOUR LOLE=LOLE+(1.-F(HLOAD(I,J))) ! LOLE = RUNNING SUM OF THE DAILY LOLP ALLR=ALLR*F(HLOAD(I,J)) ! PROBABILITY OF DAILY GEN UP STATES ENDIF DO 2 K=1,HLOAD(I,J) 2 EUE=EUE+(1.-F(K)) ! EXPECTED UNSERVED ENERGY - THE AREA ABOVE F(X) AND < PR=1 LOLH=LOLH+(1.-F(HLOAD(I,J))) ! LOLH = RUNNING SUM OF THE HOURLY LOLP 3 ENGY=ENGY+HLOAD(I,J) ALLR=1.D0-ALLR * DISPLAY RESULTS WRITE(*,*) ' LOLH LOLE puEUEppm' WRITE(*,'(2F11.6,F11.3)') LOLH,LOLE,EUE/ENGY*1.D6 WRITE(*,*) 'ENTER NEW D1,D2,D3 OR RETURN TO CONTINUE' READ(*,'(A50)') A IF(A.NE.' ') GOTO 9 OPEN(3,FILE='OP') WRITE(3,*) 'D1 MIN, D2 WIN, D3 SUM =',IFIX(D1),IFIX(D2),IFIX(D3) WRITE(3,*) ' ' WRITE(3,*) ' ' WRITE(3,*) ' DIRECT SOLUTION' WRITE(3,*) ' ' WRITE(3,*) ' LOLH LOLE ALLR' WRITE(3,'(3F11.6)') LOLH,LOLE,ALLR WRITE(3,*) ' ' WRITE(3,*) ' MWHEUE MWHTOTAL puEUEppm' WRITE(3,'(F11.5,F11.0,F11.5)') EUE,ENGY,EUE/ENGY*1.D6 * SHOW THE DIRECT SOLUTION LOLP FOR WINTER, SUMMER, AND OFF PEAK HOURS WRITE(3,*) ' ' WRITE(3,*)' WINTER SUMMER OFF PEAK' WRITE(3,*)'HR MW LOLP MW LOLP MW LOLP' DO 4 I=1,24 4 WRITE(3,'(I3,3(F7.0,F12.8))') I, & HLOAD(I, 1),1.-F(HLOAD(I, 1)), & HLOAD(I,201),1.-F(HLOAD(I,201)), & HLOAD(I, 11),1.-F(HLOAD(I, 11)) * * BEGIN THE SEQUENTIAL MONTECARLO SIMULATION, WILL BE REPORTED IN INTERVALS OF 1000, 10K, 100K, ETC YEARS DO 5 I=1,10 ISTATE(I)=1 ! START OUT WITH ALL GENERATORS TURNED ON MTTR(I)=MTTR(I)*24. ! CONVERT MEAN TIME TO REPAIR FROM DAYS TO HOURS RR(I)=1./MTTR(I) ! CALCULATE THE REPAIR RATE TRANSITION PROBABILITY FOR EACH HOUR 5 FR(I)=RR(I)*FOR(I)/(1.-FOR(I)) ! CALCULATE THE FAILURE RATE FROM THE REPAIR RATE AND FORCED OUTAGE RATE LOLE=0.D0 ! NUMBER OF LOSS OF LOAD EXPECTATION EVENTS, NO MORE THAN 1 PER DAY WILL BE COUNTED. * NOTE THAT UNITS OF DAYS PER YEAR MEANS IF A DAY HAD AN OUTAGE IT HAS AN OUTAGE. * HAVING TWO OUTAGES IN THE SAME DAY DOESN'T CHANGE THE FACT THAT IT HAD AN OUTAGE * IN THAT DAY. THIS IS TO MAINTAIN CONSISTENCY WITH THE DEFINITION OF DAYS PER YEAR. LOLH=0.D0 ! LOSS OF LOAD HOURS, ADD 1 TO EVERY HOUR THERE IS AN OBSERVED LOSS OF LOAD ALLR=0.D0 ! ANNUAL LOAD LOSS RISK, WILL COUNT THE NUMBER OF YEARS THERE ARE NO LOSS OF LOAD EVENTS EUE=0.D0 ! EXPECTED UNSERVED ENERGY, WILL SUM THE DIFFERENCE BETWEEN THE LOAD AND GENERATION AVAILABLE EACH HOUR M=3 ! SEED FOR RANDOM NUMBER GENERATOR MULT=1000. ! START OUT WITH 1000 YEARS OF SIMULATION BEFORE REPORTING RESULTS YRS=0.D0 WRITE(3,*) ' ' WRITE(3,*) ' ' WRITE(3,*) ' MONTECARLO SEQUENTIAL SOLUTION' WRITE(3,*) ' ' WRITE(3,*) ' YRS LOLH LOLE ALOLP MWHEUE' WRITE(*,*) ' MONTECARLO SEQUENTIAL SOLUTION' WRITE(*,*) ' YRS LOLH LOLE ALOLP MWHEUE' * BEGIN THE LOOP SIMULATING EACH YEAR WITH SEQUENTIAL HOURS WITHIN THAT YEAR 6 YRS=YRS+1.D0 NED=0 ! COUNT THE NUMBER OF EVENT DAYS JUST FOR THIS ONE YEAR DO 8 J=1,365 NEED=0 ! COUNT THE NUMBER OF EVENTS EACH DAY DO 8 I=1,24 GEN=0. ! THE MW GENERATION ON LINE THIS HOUR DO 7 K=1,10 ! SWEEP THROUGH ALL THE GENERATORS EACH HOUR CHANGING STATES AND SUMMING THE TOTAL GENERATION AVAILABLE R=RAN(M) IF(ISTATE(K).EQ.1) THEN IF(R.LE.FR(K)) ISTATE(K)=0 ! FLIP THE GENERATOR TO A FAILURE STATE ELSE IF(R.LE.RR(K)) ISTATE(K)=1 ! FLIP THE GENERATOR BACK ON LINE, ITS FIXED ENDIF IF(ISTATE(K).EQ.1) GEN=GEN+PMAX(K) ! SUM UP THE POWER OF THE GENS ON LINE 7 CONTINUE IF(GEN+.0001.GE.HLOAD(I,J)) THEN ! THERE IS ENOUGH GENERATION, THE .0001 INSURES THE IF TEST IS ACCURATE IN CASE REAL NUMBERS HAPPEN TO END IN .999 * INSERT SOMETHING HERE FOR THE OK CONDITION IF THERE IS ANYTHING TO CALCULATE OR REMEMBER ELSE EUE=EUE+HLOAD(I,J)-GEN ! HLOAD EXCEEDS THE GENERATION MW, ADD THE MW TIMES ONE HOUR, OR UNSERVED ENERGY IN MWH LOLH=LOLH+1. ! COUNT THE HOURS THERE IS LOSS OF LOAD IF(NEED.EQ.0) THEN NED=NED+1 ! THIS IS THE FIRST OCCURENCE OF LOSS OF LOAD THIS DAY NEED=1 ! DON'T ALLOW ANY MORE TO BE COUNTED THIS DAY ENDIF ENDIF 8 CONTINUE LOLE=LOLE+NED ! SUM UP THE EVENTS FOR EACH YEAR AND LOLE STORES THE RUNNING SUM IF(NED.EQ.0) ALLR=ALLR+1.D0 ! COUNT THE NUMBER OF YEARS THERE WERE NO EVENTS IF(YRS.LT.MULT) GOTO 6 * DISPLAY RESULTS FOR CHOSEN NUMBER OF YEARS 1K 10K 100K ETC WRITE(3,'(F8.0,3F11.6,F10.3)') & YRS,LOLH/YRS,LOLE/YRS,1.-ALLR/YRS,EUE/YRS WRITE(*,'(F8.0,3F11.6,F10.3)') & YRS,LOLH/YRS,LOLE/YRS,1.-ALLR/YRS,EUE/YRS MULT=MULT*10. WRITE(*,*) 'CONTINUE MONTECARLO? 1=YES' READ(*,*) I IF(I.EQ.1) GOTO 6 WRITE(*,*) 'FILE OP CONTAINS THE OUTPUT REPORT' WRITE(*,*) 'PRESS ENTER TO END THIS PROGRAM' READ(*,'(A1)') A END D1 MIN, D2 WIN, D3 SUM = 260 520 485 DIRECT SOLUTION LOLH LOLE ALLR 0.595133 0.100009 0.095285 MWHEUE MWHTOTAL puEUEppm 27.55358 1880818. 14.64978 WINTER SUMMER OFF PEAK HR MW LOLP MW LOLP MW LOLP 1 291. 0.00007446 271. 0.00007314 145. 0.00000098 2 312. 0.00008053 291. 0.00007446 156. 0.00000312 3 332. 0.00008740 310. 0.00008047 166. 0.00000316 4 353. 0.00015289 329. 0.00008640 176. 0.00000340 5 374. 0.00016711 349. 0.00014907 187. 0.00000421 6 395. 0.00028038 368. 0.00015913 197. 0.00000705 7 416. 0.00042523 388. 0.00020789 208. 0.00000783 8 436. 0.00064142 407. 0.00040851 218. 0.00001082 9 457. 0.00213947 426. 0.00046550 228. 0.00001179 10 478. 0.00233678 446. 0.00108486 239. 0.00002078 11 499. 0.00283150 465. 0.00216419 249. 0.00003253 12 520. 0.00410284 485. 0.00243745 260. 0.00007289 13 520. 0.00410284 485. 0.00243745 260. 0.00007289 14 499. 0.00283150 465. 0.00216419 249. 0.00003253 15 478. 0.00233678 446. 0.00108486 239. 0.00002078 16 457. 0.00213947 426. 0.00046550 228. 0.00001179 17 436. 0.00064142 407. 0.00040851 218. 0.00001082 18 416. 0.00042523 388. 0.00020789 208. 0.00000783 19 395. 0.00028038 368. 0.00015913 197. 0.00000705 20 374. 0.00016711 349. 0.00014907 187. 0.00000421 21 353. 0.00015289 329. 0.00008640 176. 0.00000340 22 332. 0.00008740 310. 0.00008047 166. 0.00000316 23 312. 0.00008053 291. 0.00007446 156. 0.00000312 24 291. 0.00007446 271. 0.00007314 145. 0.00000098 MONTECARLO SEQUENTIAL SOLUTION YRS LOLH LOLE ALOLP MWHEUE 1000. 0.887000 0.115000 0.015000 35.108 10000. 0.586300 0.102900 0.013900 27.959 100000. 0.575800 0.101980 0.013740 25.367 In the simulation below the MTTR are 10% of the original values. MONTECARLO SEQUENTIAL SOLUTION YRS LOLH LOLE ALOLP MWHEUE 1000. 0.485000 0.101000 0.052000 19.174 10000. 0.596700 0.114600 0.059200 27.417 100000. 0.614950 0.117180 0.059710 29.101 In the simulation below the MTTR are 1% of the original values. MONTECARLO SEQUENTIAL SOLUTION YRS LOLH LOLE ALOLP MWHEUE 1000. 0.685000 0.269000 0.228000 32.824 10000. 0.585600 0.243300 0.212300 27.082 100000. 0.588310 0.244330 0.214150 27.116