Next: G.2.12 Residual_ELL_Matrix Procedure Up: G.2 ELL_Matrix Class Code Previous: G.2.10 Output_ELL_Matrix Procedure
The main documentation of the Read_Harwell_Boeing_ELL_Matrix Procedure contains additional explanation of this code listing.
subroutine Read_Harwell_Boeing_ELL_Matrix (ELLM, RHS_MV, Solution_MV, &
Guess_MV, Row_Structure, &
Column_Structure, Unit, status)
! Input variables.
type(integer), intent(in), optional :: Unit ! Input unit.
! Output variables.
! ELL_Matrix to be initialized.
type(ELL_Matrix_type), intent(out) :: ELLM
type(Status_type), intent(out), optional :: status ! Exit status.
! Row and Column Base_Structures.
type(Base_Structure_type), intent(inout), target :: Row_Structure, &
Column_Structure
! RHS mathematic vector (RHSVAL).
type(Mathematic_Vector_type), optional :: RHS_MV
! Guess mathematic vector (SGUESS).
type(Mathematic_Vector_type), optional :: Guess_MV
! Solution mathematic vector (XEXACT).
type(Mathematic_Vector_type), optional :: Solution_MV
! Internal variables.
type(integer) :: A_Unit ! Actual output unit.
type(Status_type), dimension(33) :: allocate_status ! Allocation Status.
type(Status_type) :: consolidated_status ! Consolidated Status.
type(integer) :: Max_Nonzeros ! Maximum number of nonzero elements/row.
type(integer) :: Max_Nonzeros_PE ! Maximum number of nonzero elements/row,
! set to zero on non-IO PEs.
! Line 1:
type(character,name_length) :: Name ! Matrix name (TITLE).
type(character,8) :: Key ! Key for matrix (KEY).
! Line 2: Number of lines in the file, for...
type(integer) :: NTotal_Lines ! Total, not including
! header (TOTCRD).
type(integer) :: NEntries_of_Columns_Lines ! Entries of Columns
! (PTRCRD).
type(integer) :: NRows_of_Entries_Lines ! Rows of Entries (INDCRD).
type(integer) :: NValues_Lines ! Numerical values (VALCRD).
type(integer) :: NRHS_Lines ! Right-hand sides,
! including starting
! guesses and solution
! vectors (RHSCRD).
! Line 3:
type(character,3) :: Matrix_Type ! Matrix type (MXTYPE).
type(integer) :: NRows ! Number of rows (NROW).
type(integer) :: NColumns ! Number of Columns (NCOL).
type(integer) :: NEntries ! Number of nonzero entries in
! the matrix (NNZERO).
type(integer) :: NElemental_Matrices ! Number of elemental matrices
! (NELTVL).
! Line 4: Format strings for...
type(character,16) :: Entries_of_Columns_Format ! Entries of Columns
! (PTRFMT).
type(character,16) :: Rows_of_Entries_Format ! Rows of Entries (INDFMT).
type(character,20) :: Values_Format ! Values (VALFMT).
type(character,20) :: RHS_Format ! RHS values (RHSFMT).
! Line 5:
type(character,3) :: RHS_Type ! RHS type (RHSTYP).
type(integer) :: NRHS ! Number of RHS vectors, not
! including solution or guess
! vectors (NRHS).
type(integer) :: NRHS_Nonzero_Rows ! Number of RHS nonzero rows
! (NRHSIX).
! Data:
type(integer,1) :: Entries_of_Columns ! Entry numbers of each column
! (POINTR).
type(integer,1) :: Rows_of_Entries ! Row numbers of each Entry
! (ROWIND).
type(real,1) :: Values_CSC ! Numerical values for the
! Compressed Sparse Column
! (CSC) matrix (VALUES).
type(real,2) :: Values ! Numerical values for the
! ELL matrix (VALUES).
type(integer,2) :: Columns ! Column indices for the
! ELL matrix ().
type(real,2) :: Values_PE ! Numerical values for the
! ELL matrix on this PE
! (VALUES).
type(integer,2) :: Columns_PE ! Column indices for the
! ELL matrix on this PE
! ().
type(real,1) :: RHS_Values_PE ! RHS vector for each PE
! (RHSVAL).
type(real,1) :: Solution_PE ! Solution vector for each PE
! (XEXACT).
type(real,1) :: Guess_PE ! Initial guess vector for each
! PE (SGUESS).
! Temporary structures for input.
type(real,1) :: Temporary_BNV
type(Assembled_Vector_type) :: Temporary_AV
type(Distributed_Vector_type) :: Temporary_DV
! Other internal variables.
type(integer,1) :: Row_Length_Vector ! Row lengths on all the PEs.
type(integer,1) :: Column_Length_Vector ! Column lengths on all the PEs.
type(character,name_length) :: Locus_Name ! Structure locus name.
type(integer) :: column, column_PE, entry, row ! Loop parameters.
type(integer,1) :: ELL_Column ! Current ELL column number.
type(integer,1) :: Nonzeros_per_Row ! Number of nonzero entries
! per row.
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
! Verify requirements.
! ELLM, Row_Structure, Column_Structure have not been initialized.
VERIFY(.not.Initialized(ELLM),5)
VERIFY(.not.Initialized(Row_Structure),5)
VERIFY(.not.Initialized(Column_Structure),5)
! Set unit number.
if (PRESENT(Unit)) then
A_Unit = Unit
else
A_Unit = 5
end if
! Allocations and initializations.
call Initialize (allocate_status)
call Initialize (consolidated_status)
! Read in the header block.
if (this_is_IO_PE) then
Name = ' '
read (A_Unit,100) Name(1:72), Key
read (A_Unit,101) NTotal_Lines, NEntries_of_Columns_Lines, &
NRows_of_Entries_Lines, NValues_Lines, NRHS_Lines
read (A_Unit,102) Matrix_Type, NRows, NColumns, NEntries, &
NElemental_Matrices
read (A_Unit,103) Entries_of_Columns_Format, Rows_of_Entries_Format,&
Values_Format, RHS_Format
if (NRHS_Lines > 0) then
read (A_Unit,104) RHS_Type, NRHS, NRHS_Nonzero_Rows
else
RHS_Type = 'F '
NRHS = 0
NRHS_Nonzero_Rows = 0
end if
end if
! Format statements.
100 format (a72, a8)
101 format (5i14)
102 format (a3, 11x, 4i14)
103 format (2a16, 2a20)
104 format ( a3, 11x, 2i14 )
! Broadcast needed information (could put this into a
! smaller number of broadcasts (by type) later).
call Broadcast (Matrix_Type)
call Broadcast (NColumns)
call Broadcast (NElemental_Matrices)
call Broadcast (NEntries)
call Broadcast (NRHS)
call Broadcast (NRHS_Lines)
call Broadcast (NRows)
call Broadcast (NValues_Lines)
call Broadcast (RHS_Type)
! Input checks.
! Can only read real, unsymmetric, assembled matrices.
VERIFY(Matrix_Type=='RUA' .or. Matrix_Type=='rua',2)
! Can only read assembled matrices (warn on this only because
! often-used SPARSKIT routines do this incorrectly).
WARN_IF(NElemental_Matrices/=0,2)
! Can read at most one right-hand-side vector.
VERIFY(NRHS<=1,2)
! Cannot read sparse-storage RHS vectors.
VERIFY(.not.(NRHS_Lines>0 .and. RHS_Type(1:1)/='F'),2)
! Gots to be a matrix.
VERIFY(NValues_Lines>0,2)
! Set Row and Column Structures (divide evenly among the PEs).
call Initialize (Row_Length_Vector, NPEs)
call Generate_Even_Distribution (Row_Length_Vector, NRows)
Locus_Name = 'Rows'
call Initialize (Row_Structure, Row_Length_Vector, Locus_Name, &
allocate_status(1))
call Initialize (Column_Length_Vector, NPEs)
call Generate_Even_Distribution (Column_Length_Vector, NColumns)
Locus_Name = 'Columns'
call Initialize (Column_Structure, Column_Length_Vector, Locus_Name, &
allocate_status(2))
! Read in matrix structure.
call Initialize (Entries_of_Columns, NColumns+1, allocate_status(3))
call Initialize (Rows_of_Entries, NEntries, allocate_status(4))
if (this_is_IO_PE) then
read (A_Unit,Entries_of_Columns_Format) Entries_of_Columns
read (A_Unit,Rows_of_Entries_Format) Rows_of_Entries
end if
! Determine Max_Nonzeros.
call Initialize (Nonzeros_per_Row, NRows, allocate_status(5))
if (this_is_IO_PE) then
do entry = 1, NEntries
Nonzeros_per_Row(Rows_of_Entries(entry)) = &
Nonzeros_per_Row(Rows_of_Entries(entry)) + 1
end do
Max_Nonzeros = MAXVAL(Nonzeros_per_Row)
end if
call Broadcast (Max_Nonzeros)
call Finalize (Nonzeros_per_Row, allocate_status(6))
! Initialize temporary BNA ELL structure - on non-IO PEs, this is a
! dummy, unused structure.
if (this_is_IO_PE) then
Max_Nonzeros_PE = Max_Nonzeros
else
Max_Nonzeros_PE = 1
end if
call Initialize (Values, NRows, Max_Nonzeros_PE, allocate_status(7))
call Initialize (Columns, NRows, Max_Nonzeros_PE, allocate_status(8))
! Read matrix values.
call Initialize (Values_CSC, NEntries, allocate_status(9))
if (this_is_IO_PE) then
read (A_Unit,Values_Format) Values_CSC
end if
! Convert to ELL structure.
call Initialize (ELL_Column, NRows, allocate_status(10))
if (this_is_IO_PE) then
ELL_Column = 0
do column = 1, NColumns
do entry = Entries_of_Columns(column), Entries_of_Columns(column+1)-1
row = Rows_of_Entries(entry)
ELL_Column(row) = ELL_Column(row) + 1
Values(row,ELL_Column(row)) = Values_CSC(entry)
Columns(row,ELL_Column(row)) = column
end do
end do
else
ELL_Column = 1
end if
! Max_Nonzeros matches ELL_Column.
VERIFY(MAXVAL(ELL_Column)==Max_Nonzeros_PE,7)
call Finalize (ELL_Column, allocate_status(11))
call Finalize (Values_CSC, allocate_status(12))
call Finalize (Entries_of_Columns, allocate_status(13))
call Finalize (Rows_of_Entries, allocate_status(14))
! Distribute matrix to PEs.
call Initialize (Values_PE, Row_Length_Vector(this_PE), Max_Nonzeros, &
allocate_status(15))
call Initialize (Columns_PE, Row_Length_Vector(this_PE), Max_Nonzeros, &
allocate_status(16))
do column = 1, Max_Nonzeros
if (this_is_IO_PE) then
column_PE = column
else
column_PE = 1
end if
call Distribute (Values_PE(:,column), Values(:,column_PE), &
Row_Length_Vector)
call Distribute (Columns_PE(:,column), Columns(:,column_PE), &
Row_Length_Vector)
end do
call Finalize (Values, allocate_status(17))
call Finalize (Columns, allocate_status(18))
! Initialize matrix structure and set.
call Initialize (ELLM, Max_Nonzeros, Row_Structure, &
Column_Structure, Name, allocate_status(19))
call Set_Values (ELLM, Values_PE, Columns_PE)
call Finalize (Values_PE, allocate_status(20))
call Finalize (Columns_PE, allocate_status(21))
! Read Right-Hand Sides.
if (NRHS > 0) then
if (PRESENT(RHS_MV)) then
call Initialize (Temporary_BNV, NRows, allocate_status(22))
call Initialize (Temporary_AV, Row_Structure, 1, &
status=allocate_status(23))
call Initialize (Temporary_DV, Row_Structure, 1, &
status=allocate_status(24))
call Initialize (RHS_Values_PE, Row_Length_Vector(this_PE), &
status=allocate_status(25))
call Initialize (RHS_MV, Row_Structure, 'RHS: '//Name, &
allocate_status(26))
if (this_is_IO_PE) then
read (A_Unit,RHS_Format) Temporary_BNV
end if
Temporary_AV = Temporary_BNV
Temporary_DV = Temporary_AV
RHS_Values_PE = Temporary_DV
RHS_MV = RHS_Values_PE
end if
! Read starting guesses.
if (RHS_Type(2:2) == 'G' .and. PRESENT(Guess_MV)) then
call Initialize (Guess_PE, Column_Length_Vector(this_PE), &
status=allocate_status(27))
call Initialize (Guess_MV, Column_Structure, 'Guess: '//Name, &
allocate_status(28))
if (this_is_IO_PE) then
read (A_Unit,RHS_Format) Temporary_BNV
end if
Temporary_AV = Temporary_BNV
Temporary_DV = Temporary_AV
Guess_PE = Temporary_DV
Guess_MV = Guess_PE
end if
! Read solution vectors.
if (RHS_Type(3:3) == 'X' .and. PRESENT(Solution_MV)) then
call Initialize (Solution_PE, Column_Length_Vector(this_PE), &
status=allocate_status(29))
call Initialize (Solution_MV, Column_Structure, 'Solution: '//Name, &
allocate_status(30))
if (this_is_IO_PE) then
read (A_Unit,RHS_Format) Temporary_BNV
end if
Temporary_AV = Temporary_BNV
Temporary_DV = Temporary_AV
Solution_PE = Temporary_DV
Solution_MV = Solution_PE
end if
! Clean up.
call Finalize (Temporary_BNV, allocate_status(31))
call Finalize (Temporary_AV, allocate_status(32))
call Finalize (Temporary_DV, allocate_status(33))
end if
! Process status variables.
consolidated_status = allocate_status
if (PRESENT(status)) then
WARN_IF(Error(consolidated_status), 5)
status = consolidated_status
else
VERIFY(Normal(consolidated_status), 5)
end if
call Finalize (consolidated_status)
call Finalize (allocate_status)
! Verify guarantees.
VERIFY(Valid_State(ELLM),5) ! ELL_Matrix is now valid.
return
end subroutine Read_Harwell_Boeing_ELL_Matrix