C E X A M P L E 16
C
C Program name:
CHARACTER*10 MYNAME
PARAMETER (MYNAME='EXAMPLE_16')
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:
C Unfactored matrix:
INTEGER LUA(25)
C Factored matrix:
INTEGER LUF(25)
C R.H.S. Vectors:
INTEGER LUB(25)
C Solution Vectors:
INTEGER LUX(25)
INTEGER LUA0(25)
DATA LUA0(1)/0/
INTEGER LUX0(25)
INTEGER LUS0(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 Local variables:
C Data type:
INTEGER IDTYPE
C Loop counter:
INTEGER I
C Array subscript:
INTEGER LX
C Array length:
INTEGER LENX
C Disk location:
INTEGER LDISK
C Vector number:
INTEGER NV
C First nonfactored equation:
INTEGER LOWASM
C Current error:
REAL*8 EI
C Maximum overall error:
REAL*8 ERROR
C File delete status (Keep or Delete)
CHARACTER*6 K_OR_D
C Restart flag:
LOGICAL RESTART
C Delete files:
LOGICAL DELETE
C Input functions:
LOGICAL ASK
INTEGER ASK_I
C Dummy argument:
INTEGER IDUM
C FORTRAN unit number for restart data:
INTEGER LUREST
DATA LUREST/10/
C Values of ALPHA for xxDAF:
REAL*8 RZERO(1)
DATA RZERO/0.0D0/
COMPLEX*16 CZERO(1)
DATA CZERO/(0.0D0,0.0D0)/
C Profile vector for a full matrix:
INTEGER LOWEQ(1)
DATA LOWEQ/-1/
C
C Common block to communicate with fill routines:
COMMON /MYDATA/MOD, N, NRHS
C
C (1) Initialize FMS:
100 CONTINUE
CALL FMSINI
CALL FMSPSH (MYNAME)
IDUM = 0
CALL FMSIGT ('LUPR', LUPR)
CALL FMSIST ('IREST', 1)
RESTART = ASK('Is this job being restarted?')
IF(.NOT.RESTART) THEN
C This is a new job:
1 CONTINUE
WRITE (6,*) 'The FMS modules are numbered as follows:'
WRITE (6,*) ' 1 = Real Symmetric'
WRITE (6,*) ' 2 = Real Nonsymmetric'
WRITE (6,*) ' 3 = Complex Hermitian'
WRITE (6,*) ' 4 = Complex Symmetric'
WRITE (6,*) ' 5 = Complex Nonsymmetric'
MOD = ASK_I('Enter the FMS module number (1 to 5)')
IF( (MOD.LT.1) .OR. (MOD.GT.5) ) GO TO 1
N = ASK_I('Enter the number of equations')
NRHS = ASK_I('Enter the number of solution vectors')
WRITE (6,*) 'You may now alter any FMS parameter.'
WRITE (6,*) 'When you are finished, type the letters RETURN'
CALL FMSSET
CALL FMSIGT ('MEMPTR', IMD_PTR)
CALL FMSIGT ('MEMPTR', RMD_PTR)
CALL FMSIGT ('MEMPTR', CMD_PTR)
IF(MOD .LE. 2) THEN
IDTYPE = 1
ELSE
IDTYPE = 2
END IF
NUMEQ = N
LOWASM = 1
C Start at the beginning of the matrix:
CALL FMSIST ('LOWASM', 1)
C Create new files:
CALL FMSIST ('IEXIST', 0)
C Keep all files at the end of this run.
DELETE = .FALSE.
CALL FMSIST ('IFKEEP', 1)
ELSE
C This job is being restarted.
LOWASM = ASK_I('Enter the highest factored equation')
C Start at the first unfactored equation:
C NOTE: FMS will automatically round this down to the first
C unfactored segment.
LOWASM = LOWASM + 1
CALL FMSIST ('LOWASM', LOWASM)
WRITE(LUPR,2000) MOD, N, NRHS, LOWASM
CALL FMSIST ('IEXIST', 1)
DELETE = ASK('Should all files be deleted after this run')
IF(DELETE) THEN
CALL FMSIST ('IFKEEP', 0)
K_OR_D = 'DELETE'
ELSE
CALL FMSIST ('IFKEEP', 1)
K_OR_D = 'KEEP'
END IF
C
C Open the restart file:
OPEN (
1 UNIT = LUREST,
2 FILE = 'FMS_RESTART',
3 STATUS = 'OLD',
4 FORM = 'FORMATTED',
5 IOSTAT = ISTAT)
PRINT *,'Restart file ',LUREST,' was reopened'
IF(ISTAT .NE. 0) THEN
PRINT *,'Error opening restart file FMS_RESTART'
PRINT *,'Status of opening file=',ISTAT
PRINT *,'Delete FMS_RESTART if it already exists'
PRINT *,'Processing terminated.'
GO TO 900
END IF
C Read in the restart data:
C LUA(25):
DO I=1,25
READ(LUREST,1000) LUA(I)
END DO
C LUF(25):
DO I=1,25
READ(LUREST,1000) LUF(I)
END DO
C LUB(25):
DO I=1,25
READ(LUREST,1000) LUB(I)
END DO
C LUX(25):
DO I=1,25
READ(LUREST,1000) LUX(I)
END DO
CLOSE(UNIT=LUREST,STATUS=K_OR_D,IOSTAT=ISTAT)
PRINT *,'Restart file ',LUREST,' was closed with status ',
1 K_OR_D
IF(ISTAT .NE. 0) THEN
PRINT *,'Error closing restart file FMS_RESTART'
PRINT *,'Status of closing file=',ISTAT
PRINT *,'Processing terminated.'
GO TO 900
END IF
C Print the information read from the restart file:
C WRITE(LUPR,2005)
C DO I=1,10
C WRITE(LUPR,2006) I, LUA(I), LUF(I), LUB(I), LUX(I)
C END DO
C DO I=11,15
C WRITE(LUPR,2006) I, LUA(I), LUF(I)
C END DO
N = LUA( 8)
IDTYPE = LUA(11)
ISTYPE = LUA(12)
NRHS = LUB( 6)
NUMEQ = N
IF( (IDTYPE.EQ.1) .AND. (ISTYPE.EQ.1) ) MOD = 1
IF( (IDTYPE.EQ.1) .AND. (ISTYPE.EQ.2) ) MOD = 2
IF( (IDTYPE.EQ.2) .AND. (ISTYPE.EQ.3) ) MOD = 3
IF( (IDTYPE.EQ.2) .AND. (ISTYPE.EQ.1) ) MOD = 4
IF( (IDTYPE.EQ.2) .AND. (ISTYPE.EQ.2) ) MOD = 5
END IF
C
C (2) Open FMS files:
DO I=11,25
LUB(I) = 0
LUX(I) = 0
END DO
CALL FMSOV (N, IDTYPE, NRHS, 'LUB', LUB)
CALL FMSOV (N, IDTYPE, NRHS, 'LUX', LUX)
C
C Allocate a vector of storage:
IF(IDTYPE .EQ. 1) THEN
LENX = LUX(4)
CALL FMSRMG (RMD, LX, LENX)
ELSE
LENX = LUX(4)/2
CALL FMSCMG (CMD, LX, LENX)
END IF
C NOTE: When these files are being reopened (due to restart),
C the file attribute lists must be defined. In addition,
C the FMS Parameter IEXIST must be 1.
IF(LOWASM .LE. N) THEN
IF(MOD.EQ.1) CALL RSDI (LOWEQ, N, 'LUA', LUA)
IF(MOD.EQ.2) CALL RNDI (LOWEQ, N, 'LUA', LUA)
IF(MOD.EQ.3) CALL CHDI (LOWEQ, N, 'LUA', LUA)
IF(MOD.EQ.4) CALL CSDI (LOWEQ, N, 'LUA', LUA)
IF(MOD.EQ.5) CALL CNDI (LOWEQ, N, 'LUA', LUA)
CALL FMSOM (LUA, 'LUF', LUF)
IF(RESTART) GO TO 400
ELSE
C Start with the solving:
IF(MOD.EQ.1) CALL RSDI (LOWEQ, N, 'LUF', LUF)
IF(MOD.EQ.2) CALL RNDI (LOWEQ, N, 'LUF', LUF)
IF(MOD.EQ.3) CALL CHDI (LOWEQ, N, 'LUF', LUF)
IF(MOD.EQ.4) CALL CSDI (LOWEQ, N, 'LUF', LUF)
IF(MOD.EQ.5) CALL CNDI (LOWEQ, N, 'LUF', LUF)
GO TO 450
END IF
C
C (3) Write data to FMS files:
C
C Matrix File:
CALL FMSIST ('MDATAU', 1)
IF(MOD.EQ.1) CALL RSDAF
1 (LUA0, RZERO, 0, LUS0, 0, LUA, LUA0, LUX0, LUX0, 0)
IF(MOD.EQ.2) CALL RNDAF
1 (LUA0, RZERO, 0, LUS0, 0, LUA, LUA0, LUX0, LUX0, 0)
IF(MOD.EQ.3) CALL CHDAF
1 (LUA0, RZERO, 0, LUS0, 0, LUA, LUA0, LUX0, LUX0, 0)
IF(MOD.EQ.4) CALL CSDAF
1 (LUA0, CZERO, 0, LUS0, 0, LUA, LUA0, LUX0, LUX0, 0)
IF(MOD.EQ.5) CALL CNDAF
1 (LUA0, CZERO, 0, LUS0, 0, LUA, LUA0, LUX0, LUX0, 0)
C
C R.H.S. Vector File:
C
C Populate test vector:
IF(IDTYPE .EQ. 1) THEN
DO 10 I = 1,LENX
RMD(LX-1+I) = 0.0D0
10 CONTINUE
RMD(LX) = 1.0D0
ELSE
DO 11 I = 1,LENX
CMD(LX-1+I) = (0.0D0,0.0D0)
11 CONTINUE
IF(MOD.EQ.3) THEN
CMD(LX) = (1.0D0,0.0D0)
ELSE
CMD(LX) = (1.0D0,1.0D0)
END IF
END IF
C
LDISK = 1
DO 30 NV = 1,NRHS
IF(IDTYPE .EQ. 1) THEN
CALL FMSWRT (LUB(1), LDISK, RMD(LX), LUB(4))
ELSE
CALL FMSWRT (LUB(1), LDISK, CMD(LX), LUB(4))
END IF
LDISK = LDISK + LUB(4)
30 CONTINUE
C
C Open the restart file:
OPEN (
1 UNIT = LUREST,
2 FILE = 'FMS_RESTART',
3 STATUS = 'NEW',
4 FORM = 'FORMATTED',
5 IOSTAT = ISTAT)
PRINT *,'Restart file ',LUREST,' was created and opened'
IF(ISTAT .NE. 0) THEN
PRINT *,'Error opening restart file FMS_RESTART'
PRINT *,'Status of OPEN=',ISTAT
PRINT *,'Processing terminated.'
GO TO 900
END IF
C
C Save file attribute arrays for restart:
C LUA(25):
DO I=1,25
WRITE(LUREST,2001) I, LUA(I)
END DO
C LUF(25):
DO I=1,25
WRITE(LUREST,2002) I, LUF(I)
END DO
C LUB(25):
DO I=1,25
WRITE(LUREST,2003) I, LUB(I)
END DO
C LUX(25):
DO I=1,25
WRITE(LUREST,2004) I, LUX(I)
END DO
C Print the information written to the restart file:
C WRITE(LUPR,2005)
C DO I=1,10
C WRITE(LUPR,2006) I, LUA(I), LUF(I), LUB(I), LUX(I)
C END DO
C DO I=11,25
C WRITE(LUPR,2006) I, LUA(I), LUF(I)
C END DO
CLOSE(UNIT=LUREST,STATUS='KEEP',IOSTAT=ISTAT)
PRINT *,'Restart file ',LUREST,' was closed with status KEEP'
IF(ISTAT .NE. 0) THEN
PRINT *,'Error closing restart file FMS_RESTART'
PRINT *,'Status of closing file=',ISTAT
PRINT *,'Processing terminated.'
GO TO 900
END IF
C
C (4) Perform matrix algebra:
400 CONTINUE
IF(MOD.EQ.1) CALL RSDF (LUA, LUF, LUB, LUX, NUMRED)
IF(MOD.EQ.2) CALL RNDF (LUA, LUF, LUB, LUX, NUMRED)
IF(MOD.EQ.3) CALL CHDF (LUA, LUF, LUB, LUX, NUMRED)
IF(MOD.EQ.4) CALL CSDF (LUA, LUF, LUB, LUX, NUMRED)
IF(MOD.EQ.5) CALL CNDF (LUA, LUF, LUB, LUX, NUMRED)
C
C Update the restart file:
IF(.NOT.DELETE) THEN
OPEN (
1 UNIT = LUREST,
2 FILE = 'FMS_RESTART',
3 STATUS = 'OLD',
4 FORM = 'FORMATTED',
5 IOSTAT = ISTAT)
IF(ISTAT .NE. 0) THEN
PRINT *,'Error opening restart file FMS_RESTART'
PRINT *,'Status of OPEN=',ISTAT
PRINT *,'Processing terminated.'
GO TO 900
END IF
C LUA(25):
DO I=1,25
WRITE(LUREST,2001) I, LUA(I)
END DO
C LUF(25):
DO I=1,25
WRITE(LUREST,2002) I, LUF(I)
END DO
C LUB(25):
DO I=1,25
WRITE(LUREST,2003) I, LUB(I)
END DO
C LUX(25):
DO I=1,25
WRITE(LUREST,2004) I, LUX(I)
END DO
CLOSE(UNIT=LUREST,STATUS='KEEP',IOSTAT=ISTAT)
PRINT *,'Restart file ',LUREST,' was updated'
IF(ISTAT .NE. 0) THEN
PRINT *,'Error closing restart file FMS_RESTART'
PRINT *,'Status of closing file=',ISTAT
PRINT *,'Processing terminated.'
GO TO 900
END IF
END IF
450 CONTINUE
IF(MOD.EQ.1) CALL RSDS (LUF, LUB, LUX, NRHS, ISKIP)
IF(MOD.EQ.2) CALL RNDS (LUF, LUB, LUX, NRHS, ISKIP)
IF(MOD.EQ.3) CALL CHDS (LUF, LUB, LUX, NRHS, ISKIP)
IF(MOD.EQ.4) CALL CSDS (LUF, LUB, LUX, NRHS, ISKIP)
IF(MOD.EQ.5) CALL CNDS (LUF, LUB, 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
IF(IDTYPE .EQ. 1) THEN
CALL FMSRED (LUX(1), LDISK, RMD(LX), N)
DO 50 I = 1,N
EI = ABS(RMD(LX-1+I) - 1.0D0)
IF(EI .GT. ERROR) ERROR = EI
50 CONTINUE
ELSE
CALL FMSRED (LUX(1), LDISK, CMD(LX), 2*N)
DO 51 I = 1,N
EI = ABS(CMD(LX-1+I) - 1.0D0)
IF(EI .GT. ERROR) ERROR = EI
51 CONTINUE
END IF
LDISK = LDISK + LUX(4)
60 CONTINUE
WRITE(6,*) 'MAXIMUM ERROR =', ERROR
C
C (6) Close FMS files:
IF(LOWASM .LE. N) CALL FMSCM (LUA)
CALL FMSCM (LUF)
CALL FMSCV (LUB)
CALL FMSCV (LUX)
IF(IDTYPE .EQ. 1) THEN
CALL FMSRMR (RMD, LX, LENX)
ELSE
CALL FMSCMR (CMD, LX, LENX)
END IF
IF(ASK('Do you want another solution?')) GO TO 100
CALL FMSPOP (MYNAME)
CALL FMSEND
900 CONTINUE
C ----------------------------------------------------------------
C FORMAT STATEMENTS
C ----------------------------------------------------------------
1000 FORMAT (9X,I15)
2000 FORMAT (/
1 ' Problem is being restarted:'/
2 ' FMS Module number..........=',I10/
3 ' Number of equations........=',I10/
4 ' Number of R.H.S. vectors...=',I10/
5 ' Restart equation...........=',I10/)
2001 FORMAT (' LUA(',I2,')=',I15)
2002 FORMAT (' LUF(',I2,')=',I15)
2003 FORMAT (' LUB(',I2,')=',I15)
2004 FORMAT (' LUX(',I2,')=',I15)
C2005 FORMAT (/
C 1 'Contents of restart file:'/
C 2 ' LUA LUF LUB LUX'/
C 3 ' ======= ======= ======= =======')
C2006 FORMAT ('(',I2,')',4I10)
END
C=======================================================================
LOGICAL FUNCTION ASK(QUESTION)
C=======================================================================
CHARACTER* (*) QUESTION
CHARACTER*1 IYN
WRITE(6,2000) QUESTION
READ (5,1000) IYN
WRITE(6,2001) IYN
IF( (IYN .EQ. 'Y') .OR. (IYN .EQ. 'y') ) THEN
ASK = .TRUE.
ELSE
ASK = .FALSE.
END IF
RETURN
1000 FORMAT (A)
2000 FORMAT (1X,A,' (y,n)>')
2001 FORMAT (4X,'You entered ',A)
END
C=======================================================================
INTEGER FUNCTION ASK_I(STRING)
C=======================================================================
CHARACTER* (*) STRING
WRITE(6,2000) STRING
READ (5,*) ASK_I
WRITE(6,2001) ASK_I
RETURN
2000 FORMAT (1X,A,'>')
2001 FORMAT (4X,'You entered',I10)
END
C=======================================================================
SUBROUTINE RSUBLK (A, D, LOWEQ, LOCEQ, IROW1, IROW2,
1 JCOL1, JCOL2, IJSTEP)
C=======================================================================
INTEGER IROW1, IROW2, JCOL1, JCOL2, IJSTEP
INTEGER LOWEQ(*), LOCEQ(*)
REAL*8 A(0:*), D(*), DIA, OFFDIA
PARAMETER (DIA = 1.0D0)
PARAMETER (OFFDIA=-1.0D0)
COMMON /MYDATA/MOD, N, 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) = DIA
10 CONTINUE
IF(IROW1 .EQ. 1) D(1) = DFLOAT(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) = OFFDIA
20 CONTINUE
RETURN
END
C=======================================================================
SUBROUTINE RNUBLK (A, D, LOWEQ, LOCEQ, IROW1, IROW2,
1 JCOL1, JCOL2, IJSTEP)
C=======================================================================
INTEGER IROW1, IROW2, JCOL1, JCOL2, IJSTEP
INTEGER LOWEQ(*), LOCEQ(*)
REAL*8 A(0:*), D(*), OFFDIA
PARAMETER (OFFDIA=-1.0D0)
C
C Populate profile of [AL] with test data:
C The term A(I,J) is addressed as A(LOCEQ(J)+IJSTEP*I)
DO 20 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) = OFFDIA
20 CONTINUE
RETURN
END
C=======================================================================
SUBROUTINE CHUBLK (A, D, LOWEQ, LOCEQ, IROW1, IROW2,
1 JCOL1, JCOL2, IJSTEP)
C=======================================================================
INTEGER IROW1, IROW2, JCOL1, JCOL2, IJSTEP
INTEGER LOWEQ(*), LOCEQ(*)
REAL*8 D(*), DIA
COMPLEX*16 A(0:*), OFFDIA
PARAMETER (DIA = 1.0D0 )
PARAMETER (OFFDIA=(-1.0D0,0.0D0))
COMMON /MYDATA/MOD, N, 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) = DIA
10 CONTINUE
IF(IROW1 .EQ. 1) D(1) = DFLOAT(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) = OFFDIA
20 CONTINUE
RETURN
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(*), DIA, OFFDIA
PARAMETER (DIA =( 1.0D0, 1.0D0))
PARAMETER (OFFDIA=(-1.0D0,-1.0D0))
COMMON /MYDATA/MOD, N, 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) = DIA
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) = OFFDIA
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(*), OFFDIA
PARAMETER (OFFDIA=(-1.0D0,-1.0D0))
C
C Populate profile of [AL] with test data:
C The term A(I,J) is addressed as A(LOCEQ(J)+IJSTEP*I)
DO 20 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) = OFFDIA
20 CONTINUE
RETURN
END
C=======================================================================
SUBROUTINE RNUSLB (A, LOCI, LOCJ, LUFLAG, JEQN1, JEQN2, NUMEQ)
C=======================================================================
INTEGER JEQN1, JEQN2, NUMEQ
INTEGER LOCI(NUMEQ), LOCJ(JEQN1:JEQN2,2), LUFLAG(NUMEQ)
REAL*8 A(0:*), DIA, OFFDIA
PARAMETER (DIA = 1.0D0)
PARAMETER (OFFDIA=-1.0D0)
C
C Fill first row and diagonal:
LU1 = LUFLAG(1)
LI1 = LOCI(1)
DO 10 J = JEQN1,JEQN2
LIJ = LI1 + LOCJ(J,LU1)
A(LIJ) = OFFDIA
LD = LOCI(J) + LOCJ(J,LUFLAG(J))
A(LD ) = DIA
10 CONTINUE
C
C Fill column 1 if this is the first slab:
IF(JEQN1 .EQ. 1) THEN
DO 20 I = 1,NUMEQ
LIJ = LOCI(I) + LOCJ(1,LUFLAG(I))
A(LIJ) = OFFDIA
20 CONTINUE
LD = LOCI(1) + LOCJ(1,LUFLAG(1))
A(LD) = DFLOAT(NUMEQ)
END IF
RETURN
END
C=======================================================================
SUBROUTINE CNUSLB (A, LOCI, LOCJ, LUFLAG, JEQN1, JEQN2, NUMEQ)
C=======================================================================
INTEGER JEQN1, JEQN2, NUMEQ
INTEGER LOCI(NUMEQ), LOCJ(JEQN1:JEQN2,2), LUFLAG(NUMEQ)
COMPLEX*16 A(0:*), DIA, OFFDIA
PARAMETER (DIA =( 1.0D0, 1.0D0))
PARAMETER (OFFDIA=(-1.0D0,-1.0D0))
C
C Fill first row and diagonal:
LU1 = LUFLAG(1)
LI1 = LOCI(1)
DO 10 J = JEQN1,JEQN2
LIJ = LI1 + LOCJ(J,LU1)
A(LIJ) = OFFDIA
LD = LOCI(J) + LOCJ(J,LUFLAG(J))
A(LD ) = DIA
10 CONTINUE
C
C Fill column 1 if this is the first slab:
IF(JEQN1 .EQ. 1) THEN
DO 20 I = 1,NUMEQ
LIJ = LOCI(I) + LOCJ(1,LUFLAG(I))
A(LIJ) = OFFDIA
20 CONTINUE
LD = LOCI(1) + LOCJ(1,LUFLAG(1))
A(LD) = DCMPLX(NUMEQ,NUMEQ)
END IF
RETURN
END