SUBROUTINE SBAMFN G (ESF , CDAT , FDAT , LDAT , MDAT , NS G ,X Y ,FX ) C C *************************************************************** C * UNCLASSIFIED * C *************************************************************** C C CALLING SEQUENCE C SUBROUTINE SBAMFN C G (ESF , CDAT , FDAT , LDAT , MDAT , NS C G ,X C Y ,FX ) C C AUTHOR C JAMES S. VANDERGRAFT (COMPUTATIONAL ENGINEERING, INC) C CODE STRUCTURING BY C PAUL J. KRAUS & CINDY A. WELLS (COMPUTATIONAL ENGINEERING, INC) C CODE READ BY C WILLIAM H. FARR (B10, NSWCDD) & OLIVER D. SMITH (EG&G WASC) C PURPOSE C TO CALCULATE THE LIKELIHOOD FUNCTION FOR X FOR EITHER THE C BROOKS AND MOTLEY BINOMIAL OR POISSON MODEL. C DESCRIPTION C IF ESF IS SET TO 1, THE LIKELIHOOD FUNCTION FOR BROOKS AND C MOTLEY BINOMIAL MODEL IS CALCULATED. IF ESF IS SET TO 2, THE C LIKELIHOOD FUNCTION FOR BROOKS AND MOTLEY POISSON MODEL IS C CALCULATED. FOR EITHER MODEL, THE RESULT IS STORED IN FX. C ASSUMPTIONS C (NONE) C RESTRICTIONS C (NONE) C PARAMETER GLOSSARY C (NONE) C GLOBAL GLOSSARY C (NONE) C ARGUMENT LIST C GIVEN C CDAT (R) = INTERVAL FAULT COUNT DATA C ESF (I) = ESTIMATION SELECTION FLAG C FDAT (R) = FRACTIONAL PART OF CODE TEST DATA VECTOR C LDAT (R) = INTERVAL TESTING LENGTH DATA VECTOR C MDAT (R) = PREVIOUS CUMULATIVE FAULT DATA VECTOR C NS (I) = SIZE OF CDAT, FDAT, LDAT AND MDAT C X (R) = PARAMETER ESTIMATE VALUES C YIELDED C FX (R) = MINIMIZED FUNCTION VALUE C LOCAL GLOSSARY C SM (R) = TEMPORARY SUMMATION C TMP1 (R) = INTERMEDIATE CALCULATION 1 C TMP2 (R) = INTERMEDIATE CALCULATION 2 C ERRORS C (NONE) C ASSOCIATED SUBPROGRAMS C SEE AUTOMATED CODE EXAMINER (ACE) OUTPUT C REFERENCES C SEE PROLOGUE OF MAIN PROGRAM C LANGUAGE C FORTRAN 77 C C *********************** START OF DECLARATIONS ************************ C C PRE-DECLARATION DEFAULT ASSIGNMENT IMPLICIT DOUBLE PRECISION (A-Z) C C PARAMETER SPECIFICATIONS C (NONE) C C GLOBAL SPECIFICATIONS C (NONE) C C ARGUMENT SPECIFICATIONS INTEGER ESF ,NS DIMENSION CDAT(NS) ,FDAT(NS) ,LDAT(NS) 1 ,MDAT(NS) ,X(3) C C LOCAL SPECIFICATIONS INTEGER I ,J C C ********************* START OF FORMAT STATEMENTS ********************* C C (NONE) C C ********************** START OF EXECUTABLE CODE ********************** C C INITIALIZE THE MINIMIZED FUNCTION VALUE. FX = 0.0 C C CALCULATE THE MINIMIZED FUNCTION VALUE FOR THE BROOKS AND MOT- C LEY BINOMIAL MODEL OR THE BROOKS AND MOTLEY POISSON MODEL. DO 1025 I = 1, NS TMP1 = 1.0 - (1.0 - X(2))**LDAT(I) TMP2 = FDAT(I) * X(1) - X(3) * MDAT(I) IF (ESF.EQ.1) THEN SM = 0.0 DO 1000 J = 1, INT(CDAT(I)) SM = SM + LOG(TMP2 - REAL(J) + 1.0) 1000 CONTINUE C ENDDO FX = FX + SM + CDAT(I) * LOG(TMP1) + (TMP2 - CDAT(I)) 1 * LDAT(I) * LOG(1.0 - X(2)) ELSE FX = FX - TMP2 * TMP1 + CDAT(I) 1 * (LOG(TMP2) + LOG(TMP1)) ENDIF 1025 CONTINUE C ENDDO FX = -FX RETURN END SUBROUTINE SBAMGD G (ESF , CDAT , FDAT , LDAT , MDAT G ,NS , X Y ,G ) C C *************************************************************** C * UNCLASSIFIED * C *************************************************************** C C CALLING SEQUENCE C SUBROUTINE SBAMGD C G (ESF , CDAT , FDAT , LDAT , MDAT C G ,NS , X C Y ,G ) C C AUTHOR C JAMES S. VANDERGRAFT (COMPUTATIONAL ENGINEERING, INC) C CODE STRUCTURING BY C PAUL J. KRAUS & CINDY A. WELLS (COMPUTATIONAL ENGINEERING, INC) C CODE READ BY C WILLIAM H. FARR (B10, NSWCDD) & OLIVER D. SMITH (EG&G WASC) C PURPOSE C TO CALCULATE THE GRADIENT OF THE LIKELIHOOD FUNCTION FOR X FOR C EITHER THE BROOKS AND MOTLEY BINOMIAL OR POISSON MODEL. C DESCRIPTION C IF ESF IS SET TO 1, THE GRADIENT OF THE LIKELIHOOD FUNCTION FOR C BROOKS AND MOTLEY BINOMIAL MODEL IS CALCULATED AND IF ESF IS C SET TO 2, THE GRADIENT OF THE LIKELIHOOD FUNCTION FOR BROOKS C AND MOTLEY POISSON MODEL IS CALCULATED. FOR EITHER MODEL, THE C RESULTING GRADIENT VALUES ARE STORED IN G. C ASSUMPTIONS C (NONE) C RESTRICTIONS C (NONE) C PARAMETER GLOSSARY C (NONE) C GLOBAL GLOSSARY C (NONE) C ARGUMENT LIST C GIVEN C CDAT (R) = INTERVAL FAULT COUNT DATA C ESF (I) = ESTIMATION SELECTION FLAG C FDAT (R) = FRACTIONAL PART OF CODE TEST DATA VECTOR C LDAT (R) = INTERVAL TESTING LENGTH DATA VECTOR C MDAT (R) = PREVIOUS CUMULATIVE FAULT DATA VECTOR C NS (I) = SIZE OF CDAT, FDAT, LDAT AND MDAT C X (R) = PARAMETER ESTIMATE VALUES C YIELDED C G (R) = GRADIENT VALUES FOR THE FUNCTION C LOCAL GLOSSARY C SM (R) = TEMPORARY SUMMATION C TMP1 (R) = INTERMEDIATE CALCULATION 1 C TMP2 (R) = INTERMEDIATE CALCULATION 2 C TMP3 (R) = INTERMEDIATE CALCULATION 3 C ERRORS C (NONE) C ASSOCIATED SUBPROGRAMS C SEE AUTOMATED CODE EXAMINER (ACE) OUTPUT C REFERENCES C SEE PROLOGUE OF MAIN PROGRAM C LANGUAGE C FORTRAN 77 C C *********************** START OF DECLARATIONS ************************ C C PRE-DECLARATION DEFAULT ASSIGNMENT IMPLICIT DOUBLE PRECISION (A-Z) C C PARAMETER SPECIFICATIONS C (NONE) C C GLOBAL SPECIFICATIONS C (NONE) C C ARGUMENT SPECIFICATIONS INTEGER ESF ,NS DIMENSION CDAT(NS) ,FDAT(NS) ,G(3) 1 ,LDAT(NS) ,MDAT(NS) ,X(3) C C LOCAL SPECIFICATIONS INTEGER I ,J C C ********************* START OF FORMAT STATEMENTS ********************* C C (NONE) C C ********************** START OF EXECUTABLE CODE ********************** C C INITIALIZE THE GRADIENT VALUES. G(1) = 0.0 G(2) = 0.0 C C CALCULATE THE GRADIENT VALUES OF THE LIKELIHOOD FUNCTION FOR C BROOKS AND MOTLEY BINOMIAL OR POISSON MODEL. DO 1025 I = 1, NS TMP1 = 1.0 - (1.0 - X(2))**LDAT(I) TMP2 = FDAT(I) * X(1) - X(3) * MDAT(I) IF (ESF.EQ.1) THEN SM = 0.0 DO 1000 J = 1, INT(CDAT(I)) SM = SM + 1.0 / (TMP2 - REAL(J) + 1.0) 1000 CONTINUE C ENDDO G(1) = G(1) + FDAT(I) 1 * (SM + LDAT(I) * LOG(1.0 - X(2))) G(2) = G(2) + CDAT(I) * LDAT(I) / TMP1 - TMP2 1 * LDAT(I) ELSE TMP3 = CDAT(I) / TMP2 - TMP1 G(1) = G(1) + FDAT(I) * TMP3 G(2) = G(2) + (TMP2 - CDAT(I) / TMP1) * LDAT(I) 1 * (1.0 - X(2))**(LDAT(I) - 1.0) ENDIF 1025 CONTINUE C ENDDO C G(1) = -G(1) IF (ESF.EQ.1) THEN G(2) = -G(2) / (1.0 - X(2)) ENDIF RETURN END SUBROUTINE SBAMHS G (ESF , CDAT , FDAT , LDAT , MDAT G ,N , NS , X Y ,H ) C C *************************************************************** C * UNCLASSIFIED * C *************************************************************** C C CALLING SEQUENCE C SUBROUTINE SBAMHS C G (ESF , CDAT , FDAT , LDAT , MDAT C G ,N , NS , X C Y ,H ) C C AUTHOR C JAMES S. VANDERGRAFT (COMPUTATIONAL ENGINEERING, INC) C CODE STRUCTURING BY C PAUL J. KRAUS & CINDY A. WELLS (COMPUTATIONAL ENGINEERING, INC) C CODE READ BY C WILLIAM H. FARR (B10, NSWCDD) & OLIVER D. SMITH (EG&G WASC) C PURPOSE C TO CALCULATE THE HESSIAN OF THE LIKELIHOOD FUNCTION FOR X FOR C EITHER THE BROOKS AND MOTLEY BINOMIAL OR POISSON MODEL. C DESCRIPTION C IF ESF IS SET TO 1, THE HESSIAN OF THE LIKELIHOOD FUNCTION FOR C BROOKS AND MOTLEY BINOMIAL MODEL IS CALCULATED AND IF ESF IS C SET TO 2, THE HESSIAN OF THE LIKELIHOOD FUNCTION FOR BROOKS C AND MOTLEY POISSON MODEL IS CALCULATED. FOR EITHER MODEL, THE C RESULTING HESSIAN VALUES ARE STORED IN H. C ASSUMPTIONS C (NONE) C RESTRICTIONS C (NONE) C PARAMETER GLOSSARY C (NONE) C GLOBAL GLOSSARY C (NONE) C ARGUMENT LIST C GIVEN C CDAT (R) = INTERVAL FAULT COUNT DATA C ESF (I) = ESTIMATION SELECTION FLAG C FDAT (R) = FRACTIONAL PART OF CODE TEST DATA VECTOR C LDAT (R) = INTERVAL TESTING LENGTH DATA VECTOR C MDAT (R) = PREVIOUS CUMULATIVE FAULT DATA VECTOR C N (I) = NUMBER OF UNKNOWN PARAMETERS C NS (I) = SIZE OF CDAT, FDAT, LDAT AND MDAT C X (R) = PARAMETER ESTIMATE VALUES C YIELDED C H (R) = HESSIAN VALUES FOR THE FUNCTION C LOCAL GLOSSARY C D (R) = EIGENVALUES OF THE HESSIAN C DUM1 (R) = DUMMY STORAGE ARRAY 1 FOR EGNVAL CALL C DUM2 (R) = DUMMY STORAGE ARRAY 2 FOR EGNVAL CALL C EMIN (R) = VALUE TO BE ADDED TO THE DIAGONAL ELEMENTS OF C THE HESSIAN WHEN THE HESSIAN IS NOT POSITIVE C DEFINITE. C HTMP (R) = INTERMEDIATE HESSIAN STORAGE FOR EGNVAL CALL C SM (R) = TEMPORARY SUMMATION C TMP1 (R) = INTERMEDIATE CALCULATION 1 C TMP2 (R) = INTERMEDIATE CALCULATION 2 C TMP3 (R) = INTERMEDIATE CALCULATION 3 C ERRORS C (NONE) C ASSOCIATED SUBPROGRAMS C SEE AUTOMATED CODE EXAMINER (ACE) OUTPUT C REFERENCES C SEE PROLOGUE OF MAIN PROGRAM C LANGUAGE C FORTRAN 77 C C *********************** START OF DECLARATIONS ************************ C C PRE-DECLARATION DEFAULT ASSIGNMENT IMPLICIT DOUBLE PRECISION (A-Z) C C PARAMETER SPECIFICATIONS C (NONE) C C GLOBAL SPECIFICATIONS C (NONE) C C ARGUMENT SPECIFICATIONS INTEGER ESF ,N ,NS DIMENSION CDAT(NS) ,FDAT(NS) ,H(3,3) 1 ,LDAT(NS) ,MDAT(NS) ,X(3) C C LOCAL SPECIFICATIONS INTEGER I ,J DIMENSION D(3) ,DUM1(3,3) ,DUM2(3,3) 1 ,HTMP(3,3) C C ********************* START OF FORMAT STATEMENTS ********************* C C (NONE) C C ********************** START OF EXECUTABLE CODE ********************** C C INITIALIZE THE HESSIAN VALUES. DO 1025 I = 1, N DO 1000 J = 1, N H(I,J) = 0.0 1000 CONTINUE C ENDDO 1025 CONTINUE C ENDDO C C CALCULATE THE HESSIAN VALUES OF THE LIKELIHOOD FUNCTION FOR C BROOKS AND MOTLEY BINOMIAL OR POISSON MODEL. DO 1075 I = 1, NS TMP1 = 1.0 - (1.0 - X(2))**LDAT(I) TMP2 = FDAT(I) * X(1) - X(3) * MDAT(I) IF (ESF.EQ.1) THEN SM = 0.0 DO 1050 J = 1, INT(CDAT(I)) SM = SM + 1.0 / (TMP2 - REAL(J) + 1.0)**2 1050 CONTINUE C ENDDO H(1,1) = H(1,1) + SM * FDAT(I)**2 H(1,2) = H(1,2) + LDAT(I) * FDAT(I) H(2,2) = H(2,2) + TMP2 * LDAT(I) + CDAT(I) 1 * LDAT(I)**2 * (1.0 - X(2))**LDAT(I) 2 / TMP1 **2 - CDAT(I) * LDAT(I) 3 / TMP1 ELSE TMP3 = LDAT(I) * (1.0 - X(2))**(LDAT(I) - 1.0) H(1,1) = H(1,1) + CDAT(I) * (FDAT(I) / TMP2)**2 H(1,2) = H(1,2) + FDAT(I) * TMP3 H(2,2) = H(2,2) + CDAT(I) * (TMP3 / TMP1)**2 1 - (TMP2 - CDAT(I) / TMP1) * (LDAT(I) - 1.0) 2 * TMP3 / (1.0 - X(2)) ENDIF 1075 CONTINUE C ENDDO C IF (ESF.EQ.1) THEN H(1,2) = H(1,2) / (1.0 - X(2)) H(2,2) = H(2,2) / (1.0 - X(2))**2 ENDIF C C STORE THE SYMMETRIC HESSIAN VALUES. H(2,1) = H(1,2) C C COMPUTE THE EIGENVALUES OF THE HESSIAN. DO 1125 I = 1, N DO 1100 J = 1, N HTMP(I,J) = H(I,J) 1100 CONTINUE C ENDDO 1125 CONTINUE C ENDDO CALL EGNVAL G (3 , N , HTMP , DUM1 , DUM2 Y ,D ) C C WHEN THE HESSIAN IS NOT POSITIVE DEFINITE, REPLACE THE DIAGONAL C ELEMENTS OF THE HESSIAN WITH THE HESSIAN VALUE PLUS EMIN. THIS C FORCES THE HESSIAN TO BE POSITIVE DEFINITE. IF (D(1).LE.0.0) THEN EMIN = ABS(D(1)) + D(N) * 0.0001 DO 1150 I = 1, N H(I,I) = H(I,I) + EMIN 1150 CONTINUE C ENDDO ENDIF RETURN END SUBROUTINE SBAMMA G (CDAT , LDAT , FDAT , MDAT , NS , NSB G ,NSE , NSR , BOPFLG, MAXIC , PEC Y ,STAT , RFLAG , INDX , V ) C C *************************************************************** C * UNCLASSIFIED * C *************************************************************** C C CALLING SEQUENCE C SUBROUTINE SBAMMA C G (CDAT , LDAT , FDAT , MDAT , NS , NSB C G ,NSE , NSR , BOPFLG, MAXIC , PEC C Y ,STAT , RFLAG , INDX , V ) C C AUTHOR C OLIVER D. SMITH (EG&G WASC) C CODE READ BY C WILLIAM H. FARR (B10, NSWCDD) C PURPOSE C TO PERFORM THE PREQUENTIAL LIKELIHOOD ANALYSIS FOR THE BROOKS C AND MOTLEY DISCRETE SOFTWARE RELIABILITY MODELS FOR A SOFTWARE C SYSTEM. C DESCRIPTION C THE ROUTINE CONTROLS THE EXECUTIONS OF THE APPLICABLE SMFLIB C ROUTINE, SBAMMD. THE PROCESSING, AS DEFINED IN THE LITTLEWOOD C PAPER (SEE REFERENCES IN THE MAIN PROGRAM) BASICALLY INVOLVES A C SERIES OF EXECUTIONS WHERE NSB TO NSE DATA POINTS ARE UTILIZED C IN THE MODEL. C C FOR EACH ITERATION, THE MODEL PARAMETERS (ESTIMATES) FROM THE C CURRENT ITERATION ARE USED IN CONJUNCTION WITH THE OBSERVED DA- C TA FROM THE NEXT ITERATION POINT TO DETERMINE AN OVER-ALL PIC- C TURE OF THE APPLICABILITY OF THE MODEL FOR THE PARTICULAR DATA C SET. C ASSUMPTIONS C (NONE) C RESTRICTIONS C (NONE) C PARAMETER GLOSSARY C (NONE) C GLOBAL GLOSSARY C (NONE) C ARGUMENT LIST C GIVEN C BOPFLG (I) = BINOMIAL OR POISSON MODEL C CDAT (R) = ARRAY TO HOLD OBSERVED FAULT COUNT DATA C FDAT (R) = ARRAY TO HOLD FRACTION OF PROGRAM TESTED C LDAT (R) = ARRAY TO HOLD INTERVAL TESTING LENGTHS C MAXIC (I) = MAXIMUM NUMBER OF ITERATIONS C MDAT (R) = ARRAY TO HOLD CUMULATIVE FAULT CORRECTIONS C NS (I) = THE NUMBER OF ENTRIES C NSB (I) = ANALYSIS ITERATION INDEX BEGINNING C NSE (I) = ANALYSIS ITERATION INDEX ENDING C NSR (I) = ANALYSIS ITERATION RANGE C PEC (R) = PROBABILITY OF CORRECTING FAULTS WITHOUT IN- C TRODUCING NEW FAULTS C YIELDED C INDX (I) = ANALYSIS ITERATION COUNTER (AT RETURN) C RFLAG (I) = MODEL RESULTS ERROR FLAG C STAT (R) = PREQUENTIAL LIKELIHOOD STATISTIC C V (R) = ANALYSIS RESULTS ARRAY OF SIZE NSR C LOCAL GLOSSARY C COUNT (I) = ITERATION COUNTER C LOGFN (R) = LOG(N(I))! C LOGFNX (R) = LOG(N(I)-X(I))! C LOGFX (R) = LOG(X(I))! C PED (R) = PROB. OF DETECTING FAULTS INITIAL ESTIMATE C STATS (R) = STATISTICS RESULTS VECTOR C TMP1 (R) = INTERMEDIATE CALCULATION 1 C TMP2 (R) = INTERMEDIATE CALCULATION 2 C TNOF (R) = TOTAL NUMBER OF FAULTS INITIAL ESTIMATE C ERRORS C (NONE) C ASSOCIATED SUBPROGRAMS C SEE AUTOMATED CODE EXAMINER (ACE) OUTPUT C REFERENCES C SEE PROLOGUE OF MAIN PROGRAM C LANGUAGE C FORTRAN 77 C C *********************** START OF DECLARATIONS ************************ C C PRE-DECLARATION DEFAULT ASSIGNMENT IMPLICIT DOUBLE PRECISION (A-Z) C C GLOBAL SPECIFICATIONS C (NONE) C C ARGUMENT SPECIFICATIONS INTEGER BOPFLG ,INDX ,NS 1 ,NSB ,NSE ,NSR 2 ,MAXIC ,RFLAG DIMENSION CDAT(NS) ,FDAT(NS) ,LDAT(NS) 1 ,MDAT(NS) ,V(NSR) C C LOCAL SPECIFICATIONS INTEGER COUNT ,I ,J 1 ,K DIMENSION STATS(3) C C ********************* START OF FORMAT STATEMENTS ********************* C C (NONE) C C ********************** START OF EXECUTABLE CODE ********************** C C INITIALIZE THE P.L. ITERATION COUNTER, STORAGE INDEX, THE MODEL C ERROR FLAG, AND THE P.L. STATISTIC. I = NSB - 1 INDX = 0 RFLAG = 0 STAT = 0.0 C C DOWHILE (PROCESSING LEFT AND NO ERRORS) 1000 I = I + 1 IF ((I.GT.NSE) .OR. (RFLAG.NE.0)) GO TO 1150 C C COMPUTE THE INITIAL ESTIMATES OF THE PROBABILITY OF DETECT- C ING FAULTS AND THE TOTAL NUMBER OF FAULTS. PED = 0.0 TNOF = 0.0 DO 1025 J = 1, I PED = PED + CDAT(J) / ((CDAT(I) + MDAT(I)) 2 * LDAT(J) * FDAT(J)) IF (TNOF.LT.(MDAT(J) + CDAT(J)) / FDAT(J)) THEN TNOF = (MDAT(J) + CDAT(J)) / FDAT(J) ENDIF 1025 CONTINUE C ENDDO PED = PED / REAL(I) TNOF = TNOF + 1.0 C C COMPUTE MODEL ESTIMATES. CALL SBAMMD G (BOPFLG, MAXIC , I , PEC , PED G ,TNOF , FDAT , MDAT , CDAT , LDAT Y ,STATS , COUNT , RFLAG) C C INCREMENT THE ITERATION COUNTER FOR STORAGE. INDX = INDX + 1 C IF (RFLAG.EQ.0) THEN TMP1 = STATS(1) - (STATS(3) * MDAT(I)) TMP2 = 1.0 - (1.0 - STATS(2))**LDAT(I+1) C IF (BOPFLG.EQ.1) THEN C CALCULATE THE SUMS FOR THE LOG(N(I))!, LOG(X(I))!, C AND LOG(N(I)-X(I))!. K = INT(TMP1) LOGFN = 0.0 DO 1050 J = 1, K LOGFN = LOGFN + LOG(REAL(J)) 1050 CONTINUE C ENDDO K = INT(CDAT(I+1)) LOGFX = 0.0 DO 1075 J = 1, K LOGFX = LOGFX + LOG(REAL(J)) 1075 CONTINUE C ENDDO K = INT(TMP1 - CDAT(I+1)) LOGFNX = 0.0 DO 1100 J = 1, K LOGFNX = LOGFNX + LOG(REAL(J)) 1100 CONTINUE C ENDDO C V(INDX) = 0.0 1 - (LOGFN - LOGFX - LOGFNX 2 + CDAT(I+1) * LOG(TMP2) 3 + (TMP1 - CDAT(I+1)) * LOG(1.0 - TMP2)) STAT = STAT + V(INDX) ELSE C CALCULATE THE SUM FOR THE LOG(X(I))!. K = INT(CDAT(I+1)) LOGFX = 0.0 DO 1125 J = 1, K LOGFX = LOGFX + LOG(REAL(J)) 1125 CONTINUE C ENDDO C V(INDX) = 0.0 1 - (CDAT(I+1) * LOG(TMP1) 2 + CDAT(I+1) * LOG(TMP2) 3 - TMP1 * TMP2 4 - LOGFX) STAT = STAT + V(INDX) ENDIF ENDIF GO TO 1000 1150 CONTINUE C ENDWHILE RETURN END SUBROUTINE SBAMMD G (ESF , MAXIC , NS , PEC , PED G ,TNOF , FDAT , MDAT , CDAT , LDAT Y ,STATS , COUNT , RFLAG ) C C *************************************************************** C * UNCLASSIFIED * C *************************************************************** C C CALLING SEQUENCE C SUBROUTINE SBAMMD C G (ESF , MAXIC , NS , PEC , PED C G ,TNOF , FDAT , MDAT , CDAT , LDAT C Y ,STATS , COUNT , RFLAG ) C C AUTHOR C OLIVER D. SMITH (EG&G WASC) C TRUST REGION IMPLEMENTATION BY C JAMES S. VANDERGRAFT (COMPUTATIONAL ENGINEERING, INC) C TRUST REGION CODE STRUCTURING BY C PAUL J. KRAUS & CINDY A. WELLS (COMPUTATIONAL ENGINEERING, INC) C CODE READ BY C WILLIAM H. FARR (B10, NSWCDD) & OLIVER D. SMITH (EG&G WASC) C PURPOSE C TO CALCULATE THE ESTIMATES USING THE BINOMIAL AND POISSON MOD- C ELS OF THE BROOKS AND MOTLEY DISCRETE SOFTWARE RELIABILITY MOD- C EL FOR A SOFTWARE SYSTEM USING THE TRUST REGION MINIMIZATION C ALGORITHM. C DESCRIPTION C THE SBAMMD ROUTINE CAN BE EXECUTED UNDER TWO DIFFERENT MODES. C THE FIRST ARGUMENT OF THE CALL LINE (ESF) WILL BE USED BY THE C ROUTINE TO DETERMINE WHETHER THE BROOKS AND MOTLEY BINOMIAL OR C POISSON MODEL IS DESIRED. C C ONCE THE ITERATIVE EXECUTION IS INITIATED, PROCESSING WILL CON- C TINUE UNTIL ONE OF THREE POSSIBLE TERMINATION STATES IS REACH- C ED. THE LAST ARGUMENT OF THE CALL LINE (RFLAG) WILL BE SET TO C INDICATE TO THE CALLING ROUTINE THE REASON FOR EXIT. THE POS- C SIBLE REASONS FOR SBAMMD EXIT ARE: C 0 - SUCCESSFUL CONVERGENCE WITHIN THE MODEL C 1 - THE MAXIMUM ITERATION COUNT WAS REACHED C 2 - THE TRUST REGION COULD NOT BE ADJUSTED PROPERLY C 4 - THE TOTAL FAULTS ESTIMATE WAS LESS THAN THE OBSERVED C FAULTS TO DATE C C THE NUMBER OF ITERATIONS NEEDED FOR THE SUCCESSFUL CONVERGENCE C WILL BE RETURNED IN THE NEXT-TO-LAST ARGUMENT (COUNT). C C THE RESULTANT STATISTICS WILL BE RETURNED IN THE STATS ARGUMENT C FOR RFLAG VALUES OF 0 AND 1 ONLY. C ASSUMPTIONS C (NONE) C RESTRICTIONS C (NONE) C PARAMETER GLOSSARY C (NONE) C GLOBAL GLOSSARY C (NONE) C ARGUMENT LIST C GIVEN C CDAT (R) = INTERVAL FAULT COUNT DATA VECTOR C ESF (I) = ESTIMATION SELECTION FLAG C FDAT (R) = FRACTIONAL PART OF CODE TEST DATA VECTOR C LDAT (R) = INTERVAL TESTING LENGTH DATA VECTOR C MAXIC (I) = MAXIMUM NUMBER OF ITERATIONS C MDAT (R) = PREVIOUS CUMULATIVE FAULT DATA VECTOR C NS (I) = SIZE OF CDAT, FDAT, LDAT, AND MDAT C PEC (R) = PROBABILITY OF CORRECTING FAULTS WITHOUT IN- C TRODUCING NEW FAULTS (SPECIFIED VALUE) C PED (R) = PROBABILITY OF DETECTING FAULTS INITIAL ES- C TIMATE C TNOF (R) = TOTAL NUMBER OF FAULTS INITIAL ESTIMATE C YIELDED C COUNT (I) = ITERATION COUNTER C RFLAG (I) = RETURN STATUS FLAG C STATS (R) = RESULTS VECTOR C LOCAL GLOSSARY C EPSD (R) = FUNCTION VALUE DIFFERENCE MINIMUM C EPSG (R) = GRADIENT VALUE DIFFERENCE MINIMUM C ERR (R) = TRUST REGION ADJUSTMENT VALUE C FACT (R) = TRUST REGION RADIUS FACTOR C FUNUP (R) = FUNCTION DECREASING VALUE C FX (R) = FUNCTION VALUE AT XP POINTS C FY (R) = FUNCTION VALUE AT Y POINTS C G (R) = GRADIENT OF THE FUNCTION C GPROD (R) = GRADIENT PRODUCT C H (R) = HESSIAN OF THE FUNCTION C HPROD (R) = HESSIAN PRODUCT C IRED (I) = NUMBER OF TRUST REGION REDUCTIONS C N (I) = NUMBER OF UNKNOWN PARAMETERS C NDIFF (R) = FUNCTION VALUE AT DIFFERENCE C NEWSTP (R) = NEWTON STEP VALUE C NGRAD (R) = GRADIENT VALUE AT DIFFERENCE C R (R) = TRUST REGION RADIUS C RNEW (R) = TRUST REGION RADIUS FACTOR FOR 1ST REDUCTION C RTMP (R) = INTERMEDIATE TRUST REGION RADIUS C S (R) = DIAGONALS OF THE HESSIAN INVERSE C TAU (R) = TRUST REGION RADIUS TAU VALUE C TMPMI (R) = INTERMEDIATE TRUST REGION MINIMUM C TMPMX (R) = INTERMEDIATE TRUST REGION MAXIMUM C TMAX (R) = TRUST REGION MAXIMUM C TMIN (R) = TRUST REGION MINIMUM C U (R) = UPPER TRIANGULAR FACTORS OF HESSIAN C XMAX (R) = UPPER BOUNDS OF THE SOLUTION C XMIN (R) = LOWER BOUNDS OF THE SOLUTION C XP (R) = PARAMETER ESTIMATE VALUES C Y (R) = UPDATED PARAMETER ESTIMATE VALUES C ERRORS C RFLAG = 1 : THE MAXIMUM ITERATION COUNT WAS REACHED C RFLAG = 2 : THE TRUST REGION COULD NOT BE ADJUSTED PROPERLY C RFLAG = 4 : THE TOTAL FAULTS ESTIMATE WAS LESS THAN THE C OBSERVED FAULTS TO DATE C ASSOCIATED SUBPROGRAMS C SEE AUTOMATED CODE EXAMINER (ACE) OUTPUT C REFERENCES C SEE PROLOGUE OF MAIN PROGRAM C LANGUAGE C FORTRAN 77 C C *********************** START OF DECLARATIONS ************************ C C PRE-DECLARATION DEFAULT ASSIGNMENT IMPLICIT DOUBLE PRECISION (A-Z) C C PARAMETER SPECIFICATIONS C (NONE) C C GLOBAL SPECIFICATIONS C (NONE) C C ARGUMENT SPECIFICATIONS INTEGER COUNT ,ESF ,MAXIC 1 ,NS ,RFLAG DIMENSION CDAT(NS) ,FDAT(NS) ,LDAT(NS) 1 ,MDAT(NS) ,STATS(3) C C LOCAL SPECIFICATIONS INTEGER I ,IRED ,J 1 ,N DIMENSION G(3) ,H(3,3) ,NEWSTP(3) 1 ,S(3) ,TMAX(3) ,TMIN(3) 2 ,U(3,3) ,XMAX(3) ,XMIN(3) 3 ,XP(3) ,Y(3) C C ********************* START OF FORMAT STATEMENTS ********************* C C (NONE) C C ********************** START OF EXECUTABLE CODE ********************** C C SET THE COUNTER, RETURN STATUS FLAG, FUNCTION VALUE DIFFERENCE, C GRADIENT VALUE DIFFERENCE, AND NUMBER OF UNKNOWN PARAMETERS. COUNT = 1 RFLAG = -1 EPSD = 1.0E-4 EPSG = 1.0E-4 N = 2 C C SET THE UPPER AND LOWER BOUNDS FOR THE SOLUTION. XMAX(1) = 1.0E10 XMAX(2) = 0.999 XMIN(1) = 0.0 DO 1000 I = 1, NS IF (XMIN(1).LT.(MDAT(I) + CDAT(I)) / FDAT(I)) THEN XMIN(1) = (MDAT(I) + CDAT(I)) / FDAT(I) ENDIF 1000 CONTINUE C ENDDO XMIN(1) = XMIN(1) + 0.001 XMIN(2) = 1.0E-10 C C SET THE INITIAL GUESSES FOR THE PARAMETER ESTIMATES. XP(1) = TNOF XP(2) = PED IF (XP(1).LT.XMIN(1)) THEN XP(1) = XMIN(1) ELSE IF (XP(1).GT.XMAX(1)) THEN XP(1) = XMAX(1) ENDIF ENDIF IF (XP(2).LT.XMIN(2)) THEN XP(2) = XMIN(2) ELSE IF (XP(2).GT.XMAX(2)) THEN XP(2) = XMAX(2) ENDIF ENDIF C C SET THE SPECIFIED VALUE FOR PROBABILITY OF CORRECTING FAULTS C WITHOUT INTRODUCING NEW FAULTS. XP(3) = PEC Y(3) = PEC C C INITIALIZE THE TRUST REGION. CALL SBAMFN G (ESF , CDAT , FDAT , LDAT , MDAT , NS , XP Y ,FX ) CALL SBAMGD G (ESF , CDAT , FDAT , LDAT , MDAT , NS , XP Y ,G ) CALL SBAMHS G (ESF , CDAT , FDAT , LDAT , MDAT , N , NS G ,XP Y ,H ) C C COMPUTE THE NEWTON STEP. CALL FACTOR G (H , N Y ,U ) CALL SOLVE G (N , U , G Y ,NEWSTP) C C COMPUTE THE DIAGONALS OF H INVERSE. DO 1050 I = 1, N DO 1025 J = 1, N Y(J) = 0.0 1025 CONTINUE C ENDDO Y(I) = 1.0 CALL SOLVE G (N , U , Y Y ,Y ) S(I) = SQRT(Y(I)) 1050 CONTINUE C ENDDO C C COMPUTE THE TRUST REGION RADIUS. R = 0.0 DO 1075 I = 1, N RTMP = ABS(NEWSTP(I) / S(I)) R = MAX(R,RTMP) 1075 CONTINUE C ENDDO R = 0.5 * R C C DOWHILE (RFLAG .EQ. -1) 1100 IRED = 0 IF (COUNT.NE.1) THEN CALL SBAMHS G (ESF , CDAT , FDAT , LDAT , MDAT G ,N , NS , XP Y ,H ) ENDIF C C COMPUTE THE MINIMUM WITHIN THE TRUST REGION OF THE QUADRAT- C IC MODEL. C C DOWHILE ((IRED .LE. 10) .AND. (FY .GT. FUNUP)) 1125 CONTINUE C C SET UP THE TRUST REGION BOUNDARIES. DO 1150 I = 1, N TMPMI = XP(I) - R * S(I) TMPMX = XP(I) + R * S(I) TMIN(I) = MAX(XMIN(I),TMPMI) TMAX(I) = MIN(XMAX(I),TMPMX) 1150 CONTINUE C ENDDO C C COMPUTE THE MINIMUM WITHIN THE TRUST REGION. CALL CONMIN G (G , H , N , TMAX , TMIN , XP Y ,Y ) CALL SBAMFN G (ESF , CDAT , FDAT , LDAT , MDAT , NS G ,Y Y ,FY ) C C CHECK WHETHER THE FUNCTION HAS DECREASED SUFFICIENTLY. GPROD = 0.0 DO 1175 I = 1, N GPROD = GPROD + G(I) * (Y(I) - XP(I)) 1175 CONTINUE C ENDDO FUNUP = FX + 0.0001 * GPROD C IF (FY.GT.FUNUP) THEN C REDUCE THE TRUST REGION AND TRY AGAIN. IRED = IRED + 1 C IF (IRED.EQ.1) THEN C MAKE SURE THE FIRST REDUCTION IS ENOUGH TO C CAUSE A CHANGE IN THE CONSTRAINED MINIMUM. RNEW = R DO 1200 I = 1, N FACT = ABS(XP(I) - Y(I)) IF (FACT.NE.0.0) THEN IF (FACT.LT.(RNEW * S(I))) THEN RNEW = FACT / S(I) ENDIF ENDIF 1200 CONTINUE C ENDDO RTMP = R IF (RNEW.LT.R) THEN R = 0.9 * RNEW ENDIF ENDIF C C REDUCE THE TRUST REGION RADIUS. IF (((IRED.EQ.1) .AND. (RNEW.GE.RTMP)) .OR. 1 ((IRED.GT.1) .AND. (IRED.LE.10))) THEN TAU = -GPROD / (2.0 * (FY - FX - GPROD)) FACT = TAU IF (TAU.LE.0.1) THEN FACT = 0.1 ENDIF IF (TAU.GT.0.5) THEN FACT = 0.5 ENDIF R = R * FACT ENDIF ENDIF IF ((IRED.LE.10) .AND. (FY.GT.FUNUP)) GO TO 1125 C ENDWHILE C C SET RFLAG IF TRUST REGION COULD NOT BE ADJUSTED PROPERLY. IF (IRED.GT.10) THEN RFLAG = 2 ELSE C AN ACCEPTABLE POINT, Y HAS BEEN FOUND. SET RFLAG IF C MAXIMUM ITERATIONS WAS REACHED. COUNT = COUNT + 1 IF (COUNT.GT.MAXIC) THEN RFLAG = 1 ELSE C EVALUATE THE GRADIENT AT THE NEW POINT. CALL SBAMGD G (ESF , CDAT , FDAT , LDAT , MDAT G ,NS , Y Y ,G ) C C TO TEST FOR CONVERGENCE, COMPUTE THE NORM OF THE C GRADIENT AND THE NORM OF THE DIFFERENCE BETWEEN C THIS POINT AND THE PREVIOUS POINT. NDIFF = 0.0 NGRAD = 0.0 DO 1225 I = 1, N NDIFF = NDIFF + (XP(I) - Y(I))**2 NGRAD = NGRAD + G(I)**2 1225 CONTINUE C ENDDO NDIFF = SQRT(NDIFF) NGRAD = SQRT(NGRAD) C C SET RFLAG IF CONVERGENCE OCCURRED. IF ((NDIFF.LT.EPSD) .OR. (NGRAD.LT.EPSG)) THEN RFLAG = 0 ELSE C ADJUST THE TRUST REGION. HPROD = 0.0 DO 1275 I = 1, N DO 1250 J = 1, N HPROD = HPROD + H(I,J) 1 * (Y(I) - XP(I)) 2 * (Y(J) - XP(J)) 1250 CONTINUE C ENDDO 1275 CONTINUE C ENDDO ERR = ABS(FX + GPROD + 0.5 * HPROD - FY) IF (ERR.LE.(0.1 * ABS(FY))) THEN R = 2.0 * R ENDIF IF (ERR.GE.(0.75 * ABS(FY))) THEN R = 0.5 * R ENDIF C C SWAP POINTS AND FUNCTION VALUES. DO 1300 I = 1, N XP(I) = Y(I) 1300 CONTINUE C ENDDO FX = FY ENDIF ENDIF ENDIF IF (RFLAG.EQ.-1) GO TO 1100 C ENDWHILE C C SET THE RESULTANT STATISTICS. STATS(1) = XP(1) STATS(2) = XP(2) STATS(3) = XP(3) C C REDEFINE THE ERROR FLAG IF THE ESTIMATE FOR THE TOTAL NUMBER OF C FAULTS IS LESS THAN THE OBSERVED NUMBER. NOTE, THE TEST WILL C BE CONSERVATIVE IF ALL FAULTS IDENTIFIED IN CDAT FOR THE LAST C PERIOD WERE NOT LOCATED IN THE FRACTIONAL PORTION OF THE PRO- C GRAM BEING TESTED. IF (STATS(1).LT.(MDAT(NS)+CDAT(NS))) THEN RFLAG = 4 ENDIF RETURN END SUBROUTINE SBAMPR G (NPV , NS , TNOF , PED , PEC , EXF G ,EXL , ERM , FDAT , LDAT , MDAT Y ,PDAT ) C C *************************************************************** C * UNCLASSIFIED * C *************************************************************** C C CALLING SEQUENCE C SUBROUTINE SBAMPR C G (NPV , NS , TNOF , PED , PEC , EXF C G ,EXL , ERM , FDAT , LDAT , MDAT C Y ,PDAT ) C C AUTHOR C OLIVER D. SMITH (EG&G WASC) C CODE READ BY C WILLIAM H. FARR (B10, NSWCDD) C PURPOSE C TO COMPUTE THE PREDICTED NUMBER OF FAULTS IN THE NEXT TESTING C PERIOD AND TO CREATE THE PREDICTED FAULT COUNT DATA VECTOR FOR C THE BROOKS AND MOTLEY DISCRETE SOFTWARE RELIABILITY MODELS FOR C A SOFTWARE SYSTEM. C DESCRIPTION C THE FIRST ARGUMENT OF THE SBAMPR CALL LINE INDICATES WHETHER C THE PREDICTED NUMBER OF FAULTS IN THE NEXT TESTING PERIOD OR C THE ENTIRE FAULT COUNT VECTOR IS TO BE GENERATED. IF NPV IS SET C TO 1, THE ROUTINE WILL COMPUTE THE ONE PREDICTED VALUE USING C THE ARGUMENTS TNOF, PED, PEC, EXL, EXF, AND ERM. OTHERWISE, THE C ROUTINE WILL COMPUTE THE ENTIRE DATA VECTOR OF FAULT COUNTS US- C ING THE ARGUMENTS TNOF, PED, PEC, AND ALL NS ELEMENTS OF THE C FDAT, LDAT, AND MDAT VECTORS. C ASSUMPTIONS C (NONE) C RESTRICTIONS C (NONE) C PARAMETER GLOSSARY C (NONE) C GLOBAL GLOSSARY C (NONE) C ARGUMENT LIST C GIVEN C ERM (R) = FAULTS FOUND TO DATE IN FRACTIONAL PART C EXF (R) = EXPECTED FRACTIONAL PART OF CODE TESTED IN C THE NEXT PERIOD C EXL (R) = EXPECTED TESTING LENGTH OF NEXT PERIOD C FDAT (R) = FRACTIONAL PART OF CODE TESTED VECTOR C LDAT (R) = INTERVAL TESTING LENGTH DATA VECTOR C MDAT (R) = PREVIOUS CUMULATIVE FAULT DATA VECTOR C NPV (I) = NUMBER OF PREDICTED VALUES DESIRED C NS (I) = SIZE OF FDAT, LDAT, AND MDAT C PEC (R) = PROB. OF CORRECTING FAULTS WITHOUT INTRODUC- C ING NEW FAULTS. C PED (R) = PROB. OF DETECTING FAULTS C TNOF (R) = ESTIMATE OF THE TOTAL NUMBER OF FAULTS C YIELDED C PDAT (R) = VECTOR OF PREDICTED DATA VALUE(S) C LOCAL GLOSSARY C (NONE) C ERRORS C (NONE) C ASSOCIATED SUBPROGRAMS C SEE AUTOMATED CODE EXAMINER (ACE) OUTPUT C REFERENCES C SEE PROLOGUE OF MAIN PROGRAM C LANGUAGE C FORTRAN 77 C C *********************** START OF DECLARATIONS ************************ C C PRE-DECLARATION DEFAULT ASSIGNMENT IMPLICIT DOUBLE PRECISION (A-Z) C C PARAMETER SPECIFICATIONS C (NONE) C C GLOBAL SPECIFICATIONS C (NONE) C C ARGUMENT SPECIFICATIONS INTEGER NPV ,NS DIMENSION FDAT(NS) ,LDAT(NS) ,MDAT(NS) 1 ,PDAT(NPV) C C LOCAL SPECIFICATIONS INTEGER I C C ********************* START OF FORMAT STATEMENTS ********************* C C (NONE) C C ********************** START OF EXECUTABLE CODE ********************** C IF (NPV.EQ.1) THEN C COMPUTE THE ESTIMATE OF THE NUMBER OF FAULTS IN THE NEXT C TESTING PERIOD. PDAT(1) = (1.0 - (1.0 - PED)**EXL) 1 * (TNOF * EXF - PEC * ERM) ELSE C COMPUTE THE ESTIMATES FOR THE ENTIRE FAULT COUNT VECTOR. DO 1000 I = 1, NS PDAT(I) = (1.0 - (1.0 - PED)**LDAT(I)) 1 * (TNOF * FDAT(I) - PEC * MDAT(I)) 1000 CONTINUE C ENDDO ENDIF RETURN END SUBROUTINE SESHMA G (CDAT , LDAT , NS , NSB , NSE , NSR Y ,STAT , RFLAG , INDX , V ) C C *************************************************************** C * UNCLASSIFIED * C *************************************************************** C C CALLING SEQUENCE C SUBROUTINE SESHMA C G (CDAT , LDAT , NS , NSB , NSE , NSR C Y ,STAT , RFLAG , INDX , V ) C C AUTHOR C OLIVER D. SMITH (EG&G WASC) C CODE READ BY C WILLIAM H. FARR (B10, NSWCDD) C PURPOSE C TO PERFORM THE PREQUENTIAL LIKELIHOOD ANALYSIS FOR THE S-SHAPED C RELIABILITY GROWTH MODEL. C DESCRIPTION C THE ROUTINE CONTROLS THE EXECUTIONS OF THE APPLICABLE SMFLIB C ROUTINE, SESHMD. THE PROCESSING, AS DEFINED IN THE LITTLEWOOD C PAPER (SEE REFERENCES IN THE MAIN PROGRAM) BASICALLY INVOLVES A C SERIES OF EXECUTIONS WHERE NSB TO NSE DATA POINTS ARE UTILIZED C IN THE MODEL. C C FOR EACH ITERATION, THE MODEL PARAMETERS (ESTIMATES) FROM THE C CURRENT ITERATION ARE USED IN CONJUNCTION WITH THE OBSERVED DA- C TA FROM THE NEXT ITERATION POINT TO DETERMINE AN OVER-ALL PIC- C TURE OF THE APPLICABILITY OF THE MODEL FOR THE PARTICULAR DATA C SET. C ASSUMPTIONS C (NONE) C RESTRICTIONS C (NONE) C PARAMETER GLOSSARY C (NONE) C GLOBAL GLOSSARY C (NONE) C ARGUMENT LIST C GIVEN C CDAT (R) = ARRAY TO HOLD OBSERVED FAULT COUNT DATA C LDAT (R) = ARRAY TO HOLD INTERVAL TESTING LENGTHS (CUM) C NS (I) = THE NUMBER OF ENTRIES C NSB (I) = ANALYSIS ITERATION INDEX BEGINNING C NSE (I) = ANALYSIS ITERATION INDEX ENDING C NSR (I) = ANALYSIS ITERATION RANGE C YIELDED C INDX (I) = ANALYSIS ITERATION COUNTER (AT RETURN) C RFLAG (I) = MODEL RESULTS ERROR FLAG C STAT (R) = PREQUENTIAL LIKELIHOOD STATISTIC C V (R) = ANALYSIS RESULTS ARRAY OF SIZE NSR C LOCAL GLOSSARY C LOGFX (R) = LOG(X(I))! C STATS (R) = STATISTICS RESULTS ARRAY C TMP (R) = INTERMEDIATE CALCULATION C ERRORS C (NONE) C ASSOCIATED SUBPROGRAMS C SEE AUTOMATED CODE EXAMINER (ACE) OUTPUT C REFERENCES C SEE PROLOGUE OF MAIN PROGRAM C LANGUAGE C FORTRAN 77 C C *********************** START OF DECLARATIONS ************************ C C PRE-DECLARATION DEFAULT ASSIGNMENT IMPLICIT DOUBLE PRECISION (A-Z) C C GLOBAL SPECIFICATIONS C (NONE) C C ARGUMENT SPECIFICATIONS INTEGER INDX ,NS ,NSB 1 ,NSE ,NSR ,RFLAG DIMENSION CDAT(NS) ,LDAT(NS) ,V(NSR) C C LOCAL SPECIFICATIONS INTEGER I ,J ,K DIMENSION STATS(2,3) C C ********************* START OF FORMAT STATEMENTS ********************* C C (NONE) C C ********************** START OF EXECUTABLE CODE ********************** C C INITIALIZE THE P.L. ITERATION COUNTER, STORAGE INDEX, THE MODEL C ERROR FLAG, AND THE P.L. STATISTIC. I = NSB - 1 INDX = 0 RFLAG = 0 STAT = 0.0 C C DOWHILE (PROCESSING LEFT AND NO ERRORS) 1000 I = I + 1 IF ((I.GT.NSE) .OR. (RFLAG.NE.0)) GO TO 1050 C C COMPUTE MODEL ESTIMATES. CALL SESHMD G (CDAT , LDAT , I Y ,STATS , RFLAG ) C C INCREMENT THE ITERATION COUNTER FOR STORAGE. INDX = INDX + 1 C IF (RFLAG.EQ.0) THEN TMP = STATS(2,1) 1 *((1.0 + STATS(1,1) * LDAT(I)) 2 * EXP(-STATS(1,1) * LDAT(I)) 3 - (1.0 + STATS(1,1) * LDAT(I+1)) 4 * EXP(-STATS(1,1) * LDAT(I+1))) C C CALCULATE THE SUM FOR THE LOG(X(I))!. K = INT(CDAT(I+1)) LOGFX = 0.0 DO 1025 J = 1, K LOGFX = LOGFX + LOG(REAL(J)) 1025 CONTINUE C ENDDO C V(INDX) = 0.0 - (CDAT(I+1) * LOG(TMP) - TMP - LOGFX) STAT = STAT + V(INDX) ENDIF GO TO 1000 1050 CONTINUE C ENDWHILE RETURN END SUBROUTINE SESHMD G (CDAT , LDAT , NS Y ,STATS , RFLAG ) C C *************************************************************** C * UNCLASSIFIED * C *************************************************************** C C CALLING SEQUENCE C SUBROUTINE SESHMD C G (CDAT , LDAT , NS C Y ,STATS , RFLAG ) C C AUTHOR C PAUL J. KRAUS (COMPUTATIONAL ENGINEERING, INC) C CODE READ BY C WILLIAM H. FARR (B10, NSWCDD) & OLIVER D. SMITH (EG&G WASC) C PURPOSE C TO CALCULATE THE ESTIMATES USING THE MAXIMUM LIKELIHOOD METHOD C OF THE S-SHAPED RELIABILITY GROWTH MODEL. C DESCRIPTION C ONCE THE ITERATIVE EXECUTION IS INITIATED, PROCESSING WILL CON- C TINUE UNTIL ONE OF TWO POSSIBLE TERMINATION STATES IS REACHED. C THE LAST ARGUMENT OF THE CALL LINE (RFLAG) WILL BE SET TO INDI- C CATE TO THE CALLING ROUTINE THE REASON FOR EXIT. THE POS- C SIBLE REASONS FOR SESHMD EXIT ARE: C 0 - SUCCESSFUL CONVERGENCE WITHIN THE MODEL C 3 - THE DATA ARE NOT APPROPRIATE FOR THE MODEL C 4 - THE TOTAL FAULTS ESTIMATE WAS LESS THAN THE OBSERVED C FAULTS TO DATE C C THE RESULTANT STATISTICS WILL BE RETURNED IN THE STATS ARGUMENT C FOR A RFLAG VALUE OF 0 ONLY. THIS WILL INCLUDE THE TWO STATIST- C ICS AND THEIR 95% CONFIDENCE INTERVALS (CI). C ASSUMPTIONS C (NONE) C RESTRICTIONS C THIS SUBPROGRAM CONTAINS CODE WHICH DOES NOT COMPLY WITH THE C PROGRAMMING STANDARDS; GROUP LEADER APPROVAL HAS BEEN OBTAINED. C -SPECIFICALLY: C THE ARGUMENT FC IN THE ACCESS TO THE ROUTINE ZERO, IS ACTUALLY C INITIALIZED IN THE ACCESSED ROUTINE. WHEN ICOUNT HAS A VALUE OF C MORE THAN ONE, THE ARGUMENT BECOMES A BOTH (AS MARKED). HENCE, C FC APPEARS UNDEFINED IN THIS ROUTINE'S COMPILATION LISTING. C PARAMETER GLOSSARY C (NONE) C GLOBAL GLOSSARY C (NONE) C ARGUMENT LIST C GIVEN C CDAT (R) = INTERVAL FAULT COUNT DATA VECTOR C LDAT (R) = INTERVAL TESTING LENGTH DATA VECTOR C NS (I) = SIZE OF CDAT AND LDAT C YIELDED C RFLAG (I) = RETURN STATUS FLAG C STATS (R) = RESULTS ARRAY C LOCAL GLOSSARY C A (R) = LOWER BOUND OF PROPORTIONALITY CONSTANT C B (R) = UPPER BOUND OF PROPORTIONALITY CONSTANT C CTOT (R) = TOTAL FAULTS FOUND TO DATE C COUNT (I) = ITERATION COUNTER C D (R) = DELTA VALUE FOR 95% CI COMPUTATIONS C EXPI (R) = INTERMEDIATE VALUE IN STATS ASSIGNMENT C EXPJ (R) = INTERMEDIATE VALUE IN STATS ASSIGNMENT C FA (R) = VALUE OF FUNCTION AT A C FB (R) = VALUE OF FUNCTION AT B C FC (R) = VALUE OF FUNCTION AT C C ICOUNT (I) = STEP COUNTER C IFLAG (I) = CONVERGENCE FLAG C SM1 (R) = TEMPORARY SUMMATION 1 C SM2 (R) = TEMPORARY SUMMATION 2 C SM3 (R) = TEMPORARY SUMMATION 3 C TMP1 (R) = INTERMEDIATE CALCULATION 1 C TMP2 (R) = INTERMEDIATE CALCULATION 2 C ERRORS C RFLAG = 3 : THE DATA ARE NOT APPROPRIATE FOR THE MODEL C RFLAG = 4 : THE TOTAL FAULTS ESTIMATE WAS LESS THAN THE C OBSERVED FAULTS TO DATE C ASSOCIATED SUBPROGRAMS C SEE AUTOMATED CODE EXAMINER (ACE) OUTPUT C REFERENCES C SEE PROLOGUE OF MAIN PROGRAM C LANGUAGE C FORTRAN 77 C C *********************** START OF DECLARATIONS ************************ C C PRE-DECLARATION DEFAULT ASSIGNMENT IMPLICIT DOUBLE PRECISION (A-Z) C C PARAMETER SPECIFICATIONS C (NONE) C C GLOBAL SPECIFICATIONS C (NONE) C C ARGUMENT SPECIFICATIONS INTEGER NS ,RFLAG DIMENSION CDAT(NS) ,LDAT(NS) ,STATS(2,3) C C LOCAL SPECIFICATIONS INTEGER COUNT ,I ,ICOUNT 1 ,IFLAG ,J C C ********************* START OF FORMAT STATEMENTS ********************* C C (NONE) C C ********************** START OF EXECUTABLE CODE ********************** C C SET THE COUNTER AND RETURN STATUS FLAG. COUNT = 0 RFLAG = -1 C C ZERO ALL LOCATIONS OF THE RESULTS VECTOR. DO 1000 I = 1, 3 STATS(1,I) = 0.0 STATS(2,I) = 0.0 1000 CONTINUE C ENDDO C C COMPUTE THE TOTAL NUMBER OF FAULTS FOUND TO DATE. CTOT = 0.0 DO 1025 I = 1, NS CTOT = CTOT + CDAT(I) 1025 CONTINUE C ENDDO C C COMPUTE THE VALUE OF THE FUNCTION AT ITS UPPER BOUND. B = REAL(NS) / LDAT(NS) IF (B.GT.1.0) THEN B = 1.0 ENDIF FB = SESHML G (B , CDAT , CTOT , LDAT , NS ) A = B * 0.1 C C COMPUTE THE VALUE OF THE FUNCTION AT A POINT WHERE ITS SIGN IS C THE OPPOSITE OF THE COMPUTED VALUE OF THE FUNCTIONS UPPER C BOUND. C C DOWHILE (RFLAG .EQ. -1) 1050 COUNT = COUNT + 1 FA = SESHML G (A , CDAT , CTOT , LDAT , NS ) C C SET RFLAG, IF CONVERGENCE POSSIBLE. IF (FA * FB.LE.0.0) THEN RFLAG = 0 ENDIF C C SET RFLAG, IF THE DATA ARE NOT APPROPRIATE FOR THE MODEL. IF (COUNT.GE.5) THEN RFLAG = 3 ENDIF IF (RFLAG.EQ.-1) THEN A = A * 0.1 ENDIF IF (RFLAG.EQ.-1) GO TO 1050 C ENDWHILE C IF (RFLAG.EQ.0) THEN C COMPUTE THE ESTIMATE OF THE PROPORTIONALITY CONSTANT AND C STORE THE VALUE IN STATS(1,1). ICOUNT = 1 IFLAG = -1 C C DOWHILE (IFLAG .EQ. -1) C COMPUTE THE ZERO OF THE FUNCTION IN THE INTERVAL A AND C B; NOTE THE ARGUMENT FC IS INITIALIZED BY THE ZERO ROU- C TINE WHEN ICOUNT IS SET TO ONE. <