FMS Example

Example 3. General Out-of-core Matrix

This example extends the application of FMS to problems which exceed physical or virtual memory. Subroutines CSUBLK and CNUBLK are used to define matrix coefficients in the lower triangle, diagonal and upper triangle. The matrix size is only limited by the available disk space.
C       E X A M P L E   3
C
C	Data type = complex:
	PARAMETER (IDTYPE = 2)
C
C	Number of vectors to reduce during factoring:
	PARAMETER (NUMRED = 0)
C
C	Skip operations during solving (no):
	PARAMETER (ISKIP  = 0)
C
C	FMS matrix and vector file attributes:
	INTEGER     LUF(25)
	INTEGER     LUX(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     I, LX, LENX, LDISK, NV
	REAL*8      EI, ERROR
C
C	Common block to communicate with CSUBLK:
	COMMON /MYDATA/N, NRHS
	DATA LOWEQ/-1/
C
C (1)   Initialize FMS:
	CALL FMSINI
	CALL FMSIST ('IPRF', 1026)
	CALL FMSIST ('MDATAU', 2)
	CALL FMSIGT ('MEMPTR', IMD_PTR)
	CALL FMSIGT ('MEMPTR', RMD_PTR)
	CALL FMSIGT ('MEMPTR', CMD_PTR)
	WRITE (6,*) 'Enter the number of equations'
	READ  (5,*) N
	WRITE (6,*) 'Enter the number of solution vectors'
	READ  (5,*) NRHS
C
C (2)   Open FMS files:
	CALL CNDI (LOWEQ, N, 'LUF', LUF)
	CALL FMSOV (N, IDTYPE, NRHS, 'LUX', LUX)
C
C       Populate test vector:
	LENX = LUX(4)/2
	CALL FMSCMG (CMD, LX, LENX)
	DO 10 I = 1,N
	   CMD(LX-1+I) = (0.0D0,0.0D0)
   10   CONTINUE
	CMD(LX) = (1.0D0,1.0D0)
C
C (3)   Write data to FMS files:
	LDISK  = 1
	DO 30 NV = 1,NRHS
	   CALL FMSWRT (LUX(1), LDISK, CMD(LX), LUX(4))
	   LDISK = LDISK + LUX(4)
   30   CONTINUE
C
C (4)   Perform matrix algebra:
	CALL CNDF (LUF, LUF, LUX, LUX,  NUMRED)
	CALL CNDS (LUF, LUX, LUX, NRHS, ISKIP)
C
C (5)   Read data from FMS files:
C       Check the answer:
	ERROR = 0.0D0
	LDISK  = 1
	DO 60 NV = 1,NRHS
	   CALL FMSRED (LUX(1), LDISK, CMD(LX), LUX(4))
	   LDISK  = LDISK  + LUX(4)
	   DO 50 I = 1,N
	      EI = ABS(CMD(LX-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 (LUF)
	CALL FMSCV (LUX)
	CALL FMSCMR (CMD, LX, LENX)
	CALL FMSEND
	END
	
	SUBROUTINE CSUBLK (A, D, LOWEQ, LOCEQ, IROW1, IROW2,
     1                     JCOL1, JCOL2, IJSTEP)
	INTEGER    IROW1, IROW2, JCOL1, JCOL2, IJSTEP
	INTEGER    LOWEQ(*), LOCEQ(*)
	COMPLEX*16 A(0:*), D(*), ONE
	PARAMETER  (ONE=(1.0D0,1.0D0))
	COMMON /MYDATA/N, NRS
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) = CMPLX(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
	
	SUBROUTINE CNUBLK (A, D, LOWEQ, LOCEQ, IROW1, IROW2,
     1                     JCOL1, JCOL2, IJSTEP)
	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