IMPLICIT NONE INTEGER :: I, J, K, COUNT = 0, COUNT_HIER = 0 INTEGER :: MIN_DF_COLUMN, SUM_DF_COLUMN INTEGER :: N, M, Y, Z, ELEMENTS CHARACTER (LEN = 1), DIMENSION (:,:), ALLOCATABLE :: FUN, RAR_FUN_ROW, & RAR_FUN_COL INTEGER, DIMENSION (:), ALLOCATABLE :: DF_COLUMN, VARIABLES INTEGER, DIMENSION (:), ALLOCATABLE :: DF_ROW, HIER, EQNS INTEGER, DIMENSION (:), ALLOCATABLE :: ROW_FREQ INTEGER :: ROW_FREQ_MIN, HIGHEST_HIERARCHY, ROWS_ELIMINATED INTEGER :: ROWNUM_MIN_FREQ, TOP_OF_HIER CHARACTER (LEN = 30) :: MATRIX_FORMAT, COUNT_FORMAT, COL_DF_T_FORMAT, HEADER CHARACTER (LEN = 30) :: FINAL_FUN_FORMAT INTEGER :: COUNT_ROW_ORDER, DIAG_ELEMENT, COUNT_SUB_GROUP INTEGER, DIMENSION (:), ALLOCATABLE :: FINAL_ROW_ORDER, FINAL_EQNS, SUB_GROUP CHARACTER (LEN = 1), DIMENSION (:,:), ALLOCATABLE :: FINAL_FUN_MATRIX CHARACTER (LEN = 50) :: HIERARCHY_FORMAT, SUB_GROUP_FORMAT, READ_FORMAT CHARACTER (LEN = 1) :: LABEL INTEGER :: HIERARCHY_ROW_ALREADY_ELIMINATED = 0 READ_FORMAT = "( I2, 1X, I2, 1X, A)" OPEN (UNIT = 1, FILE = 'functionality.matrix', STATUS = 'OLD') READ(1,*) N READ(1,*) M READ(1,*) ELEMENTS ALLOCATE (FUN (N,M), RAR_FUN_ROW(N,M), RAR_FUN_COL(N,M)) ALLOCATE (DF_COLUMN(M), VARIABLES(M)) ALLOCATE (DF_ROW(N), HIER(N), EQNS(N), ROW_FREQ(N)) ALLOCATE (FINAL_ROW_ORDER(N), FINAL_EQNS(N), SUB_GROUP(N)) ALLOCATE (FINAL_FUN_MATRIX(N,M)) DO I = 1,N DO J = 1,M FUN(I,J) = '-' RAR_FUN_ROW(I,J) = '-' RAR_FUN_COL(I,J) = '-' END DO END DO DO I = 1,ELEMENTS READ(1,READ_FORMAT) Y, Z, LABEL FUN(Y,Z) = LABEL END DO CLOSE(1) DO I = 1,N EQNS(I) = I END DO DO I = 1,M VARIABLES(I) = I END DO HEADER = "(1X, /, 4X, 100(1X,I2)/)" FINAL_FUN_FORMAT = "(1X, I2, 1X, 100(2X,A))" MATRIX_FORMAT = "(1X, 4(2X,A))" SUB_GROUP_FORMAT = "(1X, 'SUB GROUP', I3, ' = ROW', I3)" HIERARCHY_FORMAT = "(1X, 'HIERARCHY', I3, ' = ROW', I3)" WRITE(*,HEADER) VARIABLES DO I = 1,N WRITE(*,FINAL_FUN_FORMAT) EQNS(I), (FUN(I,J), J = 1,M) END DO ! search for column degree of freedom of 1 CALL CALC_COLUMN_DF(N, M, DF_COLUMN, FUN, MIN_DF_COLUMN, SUM_DF_COLUMN) DO WHILE(SUM_DF_COLUMN /= 0) DO WHILE(MIN_DF_COLUMN == 1) CALL ONE_DF_COLUMN_ARRANGE(N, M, COUNT, DF_COLUMN, FUN, RAR_FUN_ROW, & EQNS) CALL CALC_COLUMN_DF(N, M, DF_COLUMN, FUN, MIN_DF_COLUMN, SUM_DF_COLUMN) END DO ! search for degree of freedom greater than 1 CALL MANY_DF_COLUMN_ARRANGE(N, M, COUNT, DF_COLUMN, MIN_DF_COLUMN, & FUN, RAR_FUN_ROW, COUNT_HIER, HIER, EQNS) CALL CALC_COLUMN_DF(N, M, DF_COLUMN, FUN, MIN_DF_COLUMN, SUM_DF_COLUMN) ! return to search for 1 degree of freedom END DO COL_DF_T_FORMAT = "(1X, A, I3)" COUNT_FORMAT = "(1X, A, I4)" ! ok, now that the rows have been eliminated, lets eliminate the columns WRITE(*,HEADER) VARIABLES DO I = 1,N WRITE(*,FINAL_FUN_FORMAT) EQNS(I), (RAR_FUN_ROW(I,J), J = 1,M) END DO COUNT = 0 COUNT_ROW_ORDER = 0 COUNT_SUB_GROUP = 0 DIAG_ELEMENT = 1 HIGHEST_HIERARCHY = COUNT_HIER TOP_OF_HIER = 1 IF (COUNT_HIER == 0) HIGHEST_HIERARCHY = 1 ! if the hier(element) is zero, then set it to the number of rows, because ! so that if no hierarchy exists, the algorithm will treat the entire ! system of equations as the current hierarchy ! if there are no hierarchy pointers, then highest_hierarchy set to 1 ! because hier(highest_hierarchy) would be hier(0) which doesn't exist DO I = 1,N IF (HIER(I) == 0) HIER(I) = N END DO IF (COUNT_HIER /= 0) THEN DO I = COUNT_HIER, 1, -1 WRITE(*,HIERARCHY_FORMAT) I, EQNS(HIER(I)) END DO ELSE WRITE(*,*) 'NO HIERARCHY POINTERS' END IF ! establish where the higest hierarchy is or if there are any hierarchies CALL CHECK_HIER_STATUS(N, M, RAR_FUN_ROW, HIER, COUNT_HIER, & TOP_OF_HIER, HIGHEST_HIERARCHY, HIERARCHY_ROW_ALREADY_ELIMINATED) ! set up a counter loop to end execution after all rows eliminated ROWS_ELIMINATED = 0 DO WHILE (ROWS_ELIMINATED < N) ! calculate the row frequencies on the highest hierarchy CALL CALC_ROW_FREQ(N, M, RAR_FUN_ROW, ROW_FREQ, ROW_FREQ_MIN, & ROWNUM_MIN_FREQ, TOP_OF_HIER, HIER, HIGHEST_HIERARCHY) CALL CHECK_HIER_STATUS(N, M, RAR_FUN_ROW, HIER, COUNT_HIER, & TOP_OF_HIER, HIGHEST_HIERARCHY, HIERARCHY_ROW_ALREADY_ELIMINATED) ! delete one row if the minimum frequency on the current hierarchy is 1 IF (ROW_FREQ_MIN == 1) THEN CALL ELIMINATE_ONE_ROW(N, M, RAR_FUN_ROW, RAR_FUN_COL, & ROWNUM_MIN_FREQ, COUNT, VARIABLES, COUNT_ROW_ORDER, & FINAL_ROW_ORDER, DIAG_ELEMENT, HIER, HIGHEST_HIERARCHY, & HIERARCHY_ROW_ALREADY_ELIMINATED, TOP_OF_HIER) END IF ! however, if the minumum frequency on the hierarchy is greater than one ! delete the top row IF (ROW_FREQ_MIN /= 1 .AND. ROW_FREQ_MIN /= 0) THEN CALL ELIMINATE_TOP_ROW(N, M, RAR_FUN_ROW, RAR_FUN_COL, ROW_FREQ, COUNT, & VARIABLES, COUNT_ROW_ORDER, FINAL_ROW_ORDER, DIAG_ELEMENT, & COUNT_SUB_GROUP, SUB_GROUP, EQNS, TOP_OF_HIER, HIER, & HIGHEST_HIERARCHY, HIERARCHY_ROW_ALREADY_ELIMINATED) END IF ROWS_ELIMINATED = ROWS_ELIMINATED + 1 END DO ! form the final rearranged matrix accounting for shuffled rows CALL RAR_FINAL_ROW_ORDER(N, M, FINAL_ROW_ORDER, RAR_FUN_COL, & FINAL_FUN_MATRIX, FINAL_EQNS, EQNS) WRITE(*,HEADER) VARIABLES DO I = 1,N WRITE(*,FINAL_FUN_FORMAT) FINAL_EQNS (I), (FINAL_FUN_MATRIX(I, J), J=1,M) END DO ! print out the relevant subgroups discoverd DO I = 1,COUNT_SUB_GROUP WRITE(*,SUB_GROUP_FORMAT) I, SUB_GROUP(I) END DO END SUBROUTINE CALC_COLUMN_DF(N, M, DF_COLUMN, FUN, MIN_DF_COLUMN, SUM_DF_COLUMN) IMPLICIT NONE INTEGER :: I, J INTEGER :: N, M INTEGER, DIMENSION (M) :: DF_COLUMN CHARACTER (LEN = 1), DIMENSION (N,M) :: FUN INTEGER :: MIN_DF_COLUMN, SUM_DF_COLUMN ! calculate the degrees of freedom for each row DF_COLUMN(1:M) = 0 DO I = 1,M DO J = 1,N IF (FUN(J,I) /= '-') DF_COLUMN(I) = DF_COLUMN(I) + 1 END DO END DO ! calculate the minimum degrees of freedom and the sum of all the column ! degrees of freedom so the routine knows when to stop eliminating rows, ! i.e. the sum is zero ! note that 1000 is an arbitrarily large number MIN_DF_COLUMN = 1000 DO I = 1,M IF (DF_COLUMN(I) < MIN_DF_COLUMN .AND. DF_COLUMN(I) /= 0) & MIN_DF_COLUMN = DF_COLUMN(I) END DO SUM_DF_COLUMN = SUM(DF_COLUMN) END SUBROUTINE ONE_DF_COLUMN_ARRANGE(N, M, COUNT, DF_COLUMN, FUN, & RAR_FUN_ROW, EQNS) IMPLICIT NONE INTEGER :: I, J, COUNT INTEGER :: N, M INTEGER, DIMENSION (M) :: DF_COLUMN CHARACTER (LEN = 1), DIMENSION (N,M) :: FUN, RAR_FUN_ROW INTEGER, DIMENSION (N) :: EQNS ! eliminate one row because the minimum column degree of freedom is 1 DO I = 1,M IF(DF_COLUMN(I) == 1) THEN COUNT = COUNT + 1 DO J = 1,N IF(FUN(J,I) /= '-') THEN RAR_FUN_ROW(N+1-COUNT,1:M) = FUN(J,1:M) FUN(J,1:M) = '-' EQNS(N+1-COUNT) = J END IF END DO EXIT END IF END DO END SUBROUTINE MANY_DF_COLUMN_ARRANGE(N, M, COUNT, DF_COLUMN, MIN_DF_COLUMN, & FUN, RAR_FUN_ROW, COUNT_HIER, HIER, EQNS) IMPLICIT NONE INTEGER :: I, J, COUNT, MIN_DF_COLUMN, COUNT_HIER INTEGER :: HIERARCHY INTEGER :: N, M CHARACTER (LEN = 1), DIMENSION (N,M) :: FUN, RAR_FUN_ROW INTEGER, DIMENSION (M) :: DF_COLUMN INTEGER, DIMENSION (N) :: HIER, EQNS ! eliminate all rows of one column and assign a hierarchy pointer because ! the minimum degree of freedom is greater than one HIERARCHY = 0 DO I = 1,M IF(DF_COLUMN(I) == MIN_DF_COLUMN) THEN DO J = 1,N IF(FUN(J,I) /= '-') THEN HIERARCHY = HIERARCHY + 1 COUNT = COUNT + 1 RAR_FUN_ROW(N+1-COUNT,1:M) = FUN(J,1:M) FUN(J,1:M) = '-' EQNS(N+1-COUNT) = J IF (HIERARCHY == 1) THEN COUNT_HIER = COUNT_HIER + 1 HIER(COUNT_HIER) = N+1-COUNT END IF END IF END DO EXIT END IF END DO END SUBROUTINE CALC_ROW_FREQ(N, M, RAR_FUN_ROW, ROW_FREQ, ROW_FREQ_MIN, & ROWNUM_MIN_FREQ, TOP_OF_HIER, HIER, HIGHEST_HIERARCHY) IMPLICIT NONE INTEGER :: I, J, ROW_FREQ_MIN, ROWNUM_MIN_FREQ INTEGER :: N, M, TOP_OF_HIER, HIGHEST_HIERARCHY, ken INTEGER, DIMENSION (N) :: ROW_FREQ, HIER CHARACTER (LEN = 1), DIMENSION (N,M) :: RAR_FUN_ROW ! calculate the row degree of freedom in the uppermost hierarchy ROW_FREQ(TOP_OF_HIER:HIER(HIGHEST_HIERARCHY)) = 0 DO I = TOP_OF_HIER,HIER(HIGHEST_HIERARCHY) DO J = 1,M IF (RAR_FUN_ROW(I,J) /= '-') ROW_FREQ(I) = ROW_FREQ(I) + 1 END DO END DO ! find out what row has the minimum degrees of freedom ! note that 1000 is simply an arbitrarily large number ROW_FREQ_MIN = 1000 DO I = TOP_OF_HIER,HIER(HIGHEST_HIERARCHY) IF (ROW_FREQ(I) /= 0 .AND. ROW_FREQ(I) < ROW_FREQ_MIN) THEN ROW_FREQ_MIN = ROW_FREQ(I) ROWNUM_MIN_FREQ = I END IF END DO END SUBROUTINE CALC_ROW_FREQ SUBROUTINE ELIMINATE_TOP_ROW(N, M, RAR_FUN_ROW, RAR_FUN_COL, ROW_FREQ, & COUNT, VARIABLES, COUNT_ROW_ORDER, FINAL_ROW_ORDER, DIAG_ELEMENT, & COUNT_SUB_GROUP, SUB_GROUP, EQNS, TOP_OF_HIER, HIER, HIGHEST_HIERARCHY,& HIERARCHY_ROW_ALREADY_ELIMINATED) IMPLICIT NONE INTEGER :: N, M INTEGER :: I, J, COUNT_ELIM_ROW, COUNT, L, TOP_OF_HIER, HIGHEST_HIERARCHY INTEGER, DIMENSION (N) :: ROW_FREQ CHARACTER (LEN = 1), DIMENSION (N,M) :: RAR_FUN_ROW, RAR_FUN_COL INTEGER, DIMENSION (M) :: VARIABLES INTEGER :: COUNT_ROW_ORDER, DIAG_ELEMENT, COUNT_SUB_GROUP, ROW_SHUFFLED INTEGER, DIMENSION (N) :: FINAL_ROW_ORDER, SUB_GROUP, EQNS, HIER INTEGER :: FREQ_HIERARCHY_ROW, HIERARCHY_ROW_ALREADY_ELIMINATED ! eliminate all the columns with entries in the top row because the row ! because the minimum row degree of freedom is greater than 1 COUNT_ELIM_ROW = 0 DO I = TOP_OF_HIER, HIER(HIGHEST_HIERARCHY) ROW_SHUFFLED = 0 DO L = TOP_OF_HIER, HIER(HIGHEST_HIERARCHY) IF (FINAL_ROW_ORDER(L) == I) ROW_SHUFFLED = 1 END DO IF (ROW_SHUFFLED /= 1) THEN IF (ROW_FREQ(I) /= 0) COUNT_ELIM_ROW = COUNT_ELIM_ROW + 1 IF (COUNT_ELIM_ROW == 1) THEN DO J = 1,M IF (RAR_FUN_ROW(I,J) /= '-') THEN COUNT = COUNT + 1 RAR_FUN_COL(1:N, COUNT) = RAR_FUN_ROW(1:N, J) RAR_FUN_ROW(1:N, J) = '-' VARIABLES(COUNT) = J END IF END DO COUNT_ROW_ORDER = COUNT_ROW_ORDER + 1 FINAL_ROW_ORDER(COUNT_ROW_ORDER) = I DIAG_ELEMENT = I COUNT_SUB_GROUP = COUNT_SUB_GROUP + 1 SUB_GROUP(COUNT_SUB_GROUP) = EQNS(I) END IF END IF END DO ! account for the possibility that eliminating one row above the ! hierarchy pointer may simultaneously eliminate the hierarchy row ! but not if it has already been eliminated because it was the minimum ! frequency row with a frequency of one IF (HIERARCHY_ROW_ALREADY_ELIMINATED == 0) THEN FREQ_HIERARCHY_ROW = 0 DO I = 1,M IF (RAR_FUN_ROW(HIER(HIGHEST_HIERARCHY), I) /= '-') THEN FREQ_HIERARCHY_ROW = FREQ_HIERARCHY_ROW + 1 END IF END DO IF (FREQ_HIERARCHY_ROW == 0) THEN COUNT_ROW_ORDER = COUNT_ROW_ORDER + 1 FINAL_ROW_ORDER(COUNT_ROW_ORDER) = HIER(HIGHEST_HIERARCHY) END IF END IF END SUBROUTINE ELIMINATE_TOP_ROW SUBROUTINE ELIMINATE_ONE_ROW(N, M, RAR_FUN_ROW, RAR_FUN_COL, & ROWNUM_MIN_FREQ, COUNT, VARIABLES, COUNT_ROW_ORDER, FINAL_ROW_ORDER, & DIAG_ELEMENT, HIER, HIGHEST_HIERARCHY, HIERARCHY_ROW_ALREADY_ELIMINATED, & TOP_OF_HIER) IMPLICIT NONE INTEGER :: N, M INTEGER :: I, J, COUNT, COUNT_ELIM_ROW, ROWNUM_MIN_FREQ CHARACTER (LEN = 1), DIMENSION (N,M) :: RAR_FUN_ROW, RAR_FUN_COL INTEGER, DIMENSION (M) :: VARIABLES INTEGER :: COUNT_ROW_ORDER, DIAG_ELEMENT, HIGHEST_HIERARCHY INTEGER :: FREQ_HIERARCHY_ROW INTEGER, DIMENSION (N) :: FINAL_ROW_ORDER, HIER INTEGER :: HIERARCHY_ROW_ALREADY_ELIMINATED INTEGER, DIMENSION (N) :: PRE_ROW_DF, POST_ROW_DF, POST_MINUS_PRE INTEGER :: TOP_OF_HIER ! check for the condition that eliminating one row may eliminate ! another. The possiblilty of eliminating the hierarchy pointer ! is checked for later PRE_ROW_DF(TOP_OF_HIER:HIER(HIGHEST_HIERARCHY)) = 0 DO I = TOP_OF_HIER, (HIER(HIGHEST_HIERARCHY)-1) DO J = 1,M IF (RAR_FUN_ROW(I,J) /= '-') PRE_ROW_DF(I) = PRE_ROW_DF(I) + 1 END DO END DO COUNT_ELIM_ROW = 0 ! eliminate one row because the minimum row degree of freedom is 1 DO I = 1, M IF (RAR_FUN_ROW(ROWNUM_MIN_FREQ, I) /= '-') THEN COUNT_ELIM_ROW = COUNT_ELIM_ROW + 1 IF (COUNT_ELIM_ROW == 1) THEN COUNT = COUNT + 1 RAR_FUN_COL(1:N, COUNT) = RAR_FUN_ROW(1:N, I) RAR_FUN_ROW(1:N, I) = '-' VARIABLES(COUNT) = I COUNT_ROW_ORDER = COUNT_ROW_ORDER + 1 ! shuffle a row if it is not on the current diagonal, yet above the current ! hierarchy pointer IF (ROWNUM_MIN_FREQ >= DIAG_ELEMENT) THEN FINAL_ROW_ORDER(COUNT_ROW_ORDER) = ROWNUM_MIN_FREQ DIAG_ELEMENT = DIAG_ELEMENT + 1 ELSE FINAL_ROW_ORDER(COUNT_ROW_ORDER) = I ENDIF END IF END IF END DO POST_ROW_DF(TOP_OF_HIER:HIER(HIGHEST_HIERARCHY)) = 0 DO I = TOP_OF_HIER, (HIER(HIGHEST_HIERARCHY)-1) DO J = 1,M IF (RAR_FUN_ROW(I,J) /= '-') POST_ROW_DF(I) = POST_ROW_DF(I) + 1 END DO END DO POST_MINUS_PRE = POST_ROW_DF - PRE_ROW_DF DO I = TOP_OF_HIER, (HIER(HIGHEST_HIERARCHY)-1) IF (POST_MINUS_PRE(I) < 0 .AND. I /= ROWNUM_MIN_FREQ) THEN IF (POST_ROW_DF(I) == 0) THEN COUNT_ROW_ORDER = COUNT_ROW_ORDER + 1 FINAL_ROW_ORDER(COUNT_ROW_ORDER) = I END IF END IF END DO ! account for the possibility that eliminating one row above the ! hierarchy pointer may simultaneously eliminate the hierarchy row ! but not if it has already been eliminated because it was the minimum ! frequency row with a frequency of one IF (ROWNUM_MIN_FREQ == HIER(HIGHEST_HIERARCHY)) THEN HIERARCHY_ROW_ALREADY_ELIMINATED = 1 END IF IF (HIERARCHY_ROW_ALREADY_ELIMINATED == 0) THEN FREQ_HIERARCHY_ROW = 0 DO I = 1,M IF (RAR_FUN_ROW(HIER(HIGHEST_HIERARCHY), I) /= '-') THEN FREQ_HIERARCHY_ROW = FREQ_HIERARCHY_ROW + 1 END IF END DO IF (FREQ_HIERARCHY_ROW == 0) THEN COUNT_ROW_ORDER = COUNT_ROW_ORDER + 1 FINAL_ROW_ORDER(COUNT_ROW_ORDER) = HIER(HIGHEST_HIERARCHY) END IF END IF END SUBROUTINE ELIMINATE_ONE_ROW SUBROUTINE CHECK_HIER_STATUS(N, M, RAR_FUN_ROW, HIER, COUNT_HIER, & TOP_OF_HIER, HIGHEST_HIERARCHY, HIERARCHY_ROW_ALREADY_ELIMINATED) IMPLICIT NONE INTEGER :: N, M INTEGER :: I, J, HIGHEST_HIERARCHY CHARACTER (LEN = 1), DIMENSION (N,M) :: RAR_FUN_ROW INTEGER, DIMENSION (N) :: HIER INTEGER :: COUNT_HIER, TOP_OF_HIER, HIER_STATUS INTEGER :: HIERARCHY_ROW_ALREADY_ELIMINATED ! determine which hierarchy pointer is currently the highest and ! if there are no hierarchy pointers set the current hierarchy ! to be the last row (entire matrix) HIER_STATUS = 0 DO I = TOP_OF_HIER, HIER(HIGHEST_HIERARCHY) DO J = 1, M IF (RAR_FUN_ROW(I, J) /= '-') HIER_STATUS = 1 END DO END DO IF (HIER_STATUS == 0) THEN TOP_OF_HIER = HIER(HIGHEST_HIERARCHY) + 1 COUNT_HIER = COUNT_HIER - 1 HIGHEST_HIERARCHY = COUNT_HIER IF (COUNT_HIER == 0) HIGHEST_HIERARCHY = N HIERARCHY_ROW_ALREADY_ELIMINATED = 0 END IF END SUBROUTINE CHECK_HIER_STATUS SUBROUTINE RAR_FINAL_ROW_ORDER(N, M, FINAL_ROW_ORDER, RAR_FUN_COL, & FINAL_FUN_MATRIX, FINAL_EQNS, EQNS) IMPLICIT NONE INTEGER :: N, M, I, NEXT_ROW CHARACTER (LEN = 1), DIMENSION (N,M) :: RAR_FUN_COL, FINAL_FUN_MATRIX INTEGER, DIMENSION (N) :: FINAL_ROW_ORDER, EQNS, FINAL_EQNS ! properly order the rows in the output matrix since they may have been ! shuffled during the variable ordering algorithm DO I = 1,N NEXT_ROW = FINAL_ROW_ORDER(I) FINAL_FUN_MATRIX(I, 1:M) = RAR_FUN_COL(NEXT_ROW, 1:M) FINAL_EQNS(I) = EQNS(NEXT_ROW) END DO END SUBROUTINE RAR_FINAL_ROW_ORDER