C E X A M P L E 7
C
C Program name:
CHARACTER*9 MYNAME
PARAMETER (MYNAME='EXAMPLE_7')
C
C Input data functions:
INTEGER ASK_I
C
C Data type = complex:
PARAMETER (IDTYPE = 2)
C
C Number of vectors to reduce during factoring:
INTEGER NUMRED
C
C Skip operations during solving (no):
INTEGER ISKIP
C
C FMS matrix and vector file attributes:
INTEGER LUA(25), LUASUB(25)
INTEGER LUX(25), LUXSUB(25)
C
C FMS memory management requires the following arrays:
POINTER (CMD_PTR, CMD)
POINTER (RMD_PTR, RMD)
POINTER (IMD_PTR, IMD)
COMPLEX*16 CMD(0:1)
REAL*8 RMD(0:1)
INTEGER IMD(0:1)
C
C Profile vector for a full matrix;
INTEGER LOWEQ(1)
C
C Local variables:
INTEGER N, NRHS, I, LDISK, NV
INTEGER L_X, L_XSUB, L_ASUB, LENX
REAL*8 EI, ERROR
C
C Common block to communicate with CSUBLK:
COMMON /MYDATA/N, NEQSUB, NRHS
DATA LOWEQ/-1/
C
C (1) Initialize FMS:
CALL FMSINI
CALL FMSPSH (MYNAME)
CALL FMSIGT ('MEMPTR', IMD_PTR)
CALL FMSIGT ('MEMPTR', RMD_PTR)
CALL FMSIGT ('MEMPTR', CMD_PTR)
N = ASK_I ('Enter the number of equations')
NEQSUB = ASK_I ('Enter the substructure equation NEQSUB')
NRHS = ASK_I ('Enter the number of solution vectors')
C
C (2) Open FMS files:
CALL FMSIST ('NEQSUB', NEQSUB)
CALL CNDI (LOWEQ, N, 'LUA', LUA)
CALL FMSOV (N, IDTYPE, NRHS, 'LUX', LUX)
C
C Populate test vector:
LENX = LUX(4)/2
CALL FMSCMG (CMD, L_X, LENX)
DO I = 1,LENX
CMD(L_X-1+I) = (0.0D0,0.0D0)
END DO
CMD(L_X) = (1.0D0,1.0D0)
C
C (3) Write data to FMS files:
LDISK = 1
DO NV = 1,NRHS
CALL FMSWRT (LUX(1), LDISK, CMD(L_X), LUX(4))
LDISK = LDISK + LUX(4)
END DO
C
C (4) Perform matrix algebra:
C
C Factor the matrix down to the substructure:
CALL FMSIST ('MDATAU', 2)
NUMRED = 0
CALL CNDF (LUA, LUA, LUX, LUX, NUMRED)
C
C Reduce the RHS vectors down to the substructure.
ISKIP = 2
CALL CNDS (LUA, LUX, LUX, NRHS, ISKIP)
C
C Complete substructure solution:
NSUB = N - NEQSUB + 1
IF(NSUB .GT. 0) THEN
C Substructuring selected:
C Turn off further substructuring:
CALL FMSIST ('NEQSUB', N+1)
C Turn off calling user subroutine:
CALL FMSIST ('MDATAU', 0)
C
C Allocate memory to store the submatrix:
CALL FMSCMG (CMD, L_ASUB, NSUB*NSUB)
C
C Extract the substructure matrix:
CALL CNDEX (LUA, CMD(L_ASUB), NSUB, NSUB)
C
C Form an incore FMS matrix for [ASUB]:
CALL CNDANN (CMD(L_ASUB), NSUB, NSUB, LUASUB)
C
C Allocate memory to store the submatrix vector:
CALL FMSCMG (CMD, L_XSUB, NSUB*NRHS)
C
C Read the reduced vectors {BSUB}:
LBUF = L_X - 1 + (NEQSUB-1)
LX = L_XSUB - 1
LDISK = 1
DO IRHS = 1,NRHS
CALL FMSRED (LUX(1), LDISK, CMD(L_X), LUX(4))
DO I=1,NSUB
CMD(LX+I) = CMD(LBUF+I)
END DO
LDISK = LDISK + LUX(4)
LX = LX + NSUB
END DO
C
C Form an incore FMS vector for {XSUB}:
CALL FMSOV2 (NSUB, IDTYPE, NRHS, CMD(L_XSUB), NSUB, LUXSUB)
C
C Factor the submatrix:
NUMRED = 0
CALL CNDF (LUASUB, LUASUB, LUXSUB, LUXSUB, NUMRED)
C
C Perform forward reduction and back substitution on {XSUB}:
ISKIP = 0
CALL CNDS (LUASUB, LUXSUB, LUXSUB, NRHS, ISKIP)
C
C Close substructure matrix and return storage:
CALL FMSCM (LUASUB)
CALL FMSCMR (CMD, L_ASUB, NSUB*NSUB)
C
C Place the solution {XSUB} back in {X}:
C
C Read the reduced vectors {BSUB}:
LBUF = L_X - 1 + (NEQSUB-1)
LX = L_XSUB - 1
LDISK = 1
DO IRHS = 1,NRHS
CALL FMSRED (LUX(1), LDISK, CMD(L_X), LUX(4))
DO I=1,NSUB
CMD(LBUF+I) = CMD(LX+I)
END DO
CALL FMSWRT (LUX(1), LDISK, CMD(L_X), LUX(4))
LDISK = LDISK + LUX(4)
LX = LX + NSUB
END DO
C
C Close substructure vector file and return storage:
CALL FMSCV (LUXSUB)
CALL FMSCMR (CMD, L_XSUB, NSUB*NRHS)
C
C Restore the value of NEQSUB:
CALL FMSIST ('NEQSUB', NEQSUB)
END IF
C
C Complete the back-substitution for {X1}:
ISKIP = 1
CALL CNDS (LUA, LUX, LUX, NRHS, ISKIP)
C
C (5) Read data from FMS files:
C Check the answer:
ERROR = 0.0D0
LDINC = LUX(4)
LDISK = 1 - LDINC
DO 60 NV = 1,NRHS
LDISK = LDISK + LDINC
CALL FMSRED (LUX(1), LDISK, CMD(L_X), LUX(4))
DO 50 I = 1,N
EI = ABS(CMD(L_X-1+I) - 1.0D0)
IF(EI .GT. ERROR) ERROR = EI
50 CONTINUE
60 CONTINUE
WRITE(6,*) 'MAXIMUM ERROR =', ERROR
C
C (6) Close FMS files:
CALL FMSCM (LUA)
CALL FMSCV (LUX)
CALL FMSCMR (CMD, L_X, LENX)
CALL FMSPOP (MYNAME)
CALL FMSEND
END
C=======================================================================
SUBROUTINE CSUBLK (A, D, LOWEQ, LOCEQ, IROW1, IROW2,
1 JCOL1, JCOL2, IJSTEP)
C=======================================================================
INTEGER IROW1, IROW2, JCOL1, JCOL2, IJSTEP
INTEGER LOWEQ(*), LOCEQ(*)
COMPLEX*16 A(0:*), D(*), ONE
PARAMETER (ONE=(1.0D0,1.0D0))
COMMON /MYDATA/N, NEQSUB, NRHS
C
C Populate the diagonal with test data:
IF(IROW2 .EQ. JCOL2) THEN
C This is a diagonal block:
DO 10 I = IROW1,IROW2
D(I) = ONE
10 CONTINUE
IF(IROW1 .EQ. 1) D(1) = DCMPLX(N,N)
END IF
C
C Populate profile of [AL] with test data:
C The term A(I,J) is addressed as A(LOCEQ(I)+IJSTEP*J)
DO 20 I = IROW1,IROW2
J = LOWEQ(I)
IF( (J .GE. JCOL1) .AND.
1 (J .LE. JCOL2) .AND.
2 (J .LT. I) ) A(LOCEQ(I) + IJSTEP*J) = -ONE
20 CONTINUE
RETURN
END
C=======================================================================
SUBROUTINE CNUBLK (A, D, LOWEQ, LOCEQ, IROW1, IROW2,
1 JCOL1, JCOL2, IJSTEP)
C=======================================================================
INTEGER IROW1, IROW2, JCOL1, JCOL2, IJSTEP
INTEGER LOWEQ(*), LOCEQ(*)
COMPLEX*16 A(0:*), D(*), ONE
PARAMETER (ONE=(1.0D0,1.0D0))
C
C Populate profile of [AU] with test data:
C The term A(I,J) is addressed as A(LOCEQ(J)+IJSTEP*I)
DO 10 J = JCOL1,JCOL2
I = LOWEQ(J)
IF( (I .GE. IROW1) .AND.
1 (I .LE. IROW2) .AND.
2 (I .LT. J) ) A(LOCEQ(J) + IJSTEP*I) = -ONE
10 CONTINUE
RETURN
END
C=======================================================================
INTEGER FUNCTION ASK_I(STRING)
C=======================================================================
CHARACTER* (*) STRING
WRITE(6,2000) STRING
READ (5,*) ASK_I
RETURN
2000 FORMAT (1X,A,'>')
END