#! /usr/bin/env python3
#
def i4_log_10 ( i ):
#*****************************************************************************80
#
## i4_log_10() returns the integer part of the logarithm base 10 of ABS(X).
#
# Example:
#
# I VALUE
# ----- --------
# 0 0
# 1 0
# 2 0
# 9 0
# 10 1
# 11 1
# 99 1
# 100 2
# 101 2
# 999 2
# 1000 3
# 1001 3
# 9999 3
# 10000 4
#
# Discussion:
#
# I4_LOG_10 ( I ) + 1 is the number of decimal digits in I.
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 08 May 2013
#
# Author:
#
# John Burkardt
#
# Input:
#
# integer I, the number whose logarithm base 10 is desired.
#
# Output:
#
# integer VALUE, the integer part of the logarithm base 10 of
# the absolute value of X.
#
if ( i == 0 ):
value = 0
else:
value = 0
ten_pow = 10
i_abs = abs ( i )
while ( ten_pow <= i_abs ):
value = value + 1
ten_pow = ten_pow * 10
return value
def i4vec_search_binary_a ( n, a, b ):
#*****************************************************************************80
#
## i4vec_search_binary_a() searches an ascending sorted I4VEC.
#
# Discussion:
#
# Binary search is used.
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 30 September 2015
#
# Author:
#
# John Burkardt
#
# Reference:
#
# Donald Kreher, Douglas Simpson,
# Combinatorial Algorithms,
# CRC Press, 1998, page 26.
#
# Input:
#
# integer N, the number of elements in the vector.
#
# integer A(N), the array to be searched. A must
# be sorted in ascending order.
#
# integer B, the value to be searched for.
#
# Output:
#
# integer INDX, the result of the search.
# -1, B does not occur in A.
# I, A(I) = B.
#
indx = - 1
low = 0
high = n - 1
while ( low <= high ):
mid = ( ( low + high ) // 2 )
if ( a[mid] == b ):
indx = mid
break
elif ( a[mid] < b ):
low = mid + 1
elif ( b < a[mid] ):
high = mid - 1
return indx
def r8ccs_test ( ):
#*****************************************************************************80
#
## r8ccs_test() tests r8ccs().
#
# Discussion:
#
# CCS is the compressed column storage format.
#
# Associated with this format, we have an M by N matrix
# with NZ_NUM nonzero entries. We construct the column pointer
# vector COL of length N+1, such that entries of column J will be
# stored in positions COL(J) through COL(J+1)-1. This indexing
# refers to both the ROW and A vectors, which store the row indices
# and the values of the nonzero entries. The entries of the
# ROW vector corresponding to each column are assumed to be
# ascending sorted.
#
# The CCS format is equivalent to the MATLAB "sparse" format,
# and the Harwell Boeing "real unsymmetric assembled" (RUA) format.
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 21 June 2022
#
# Author:
#
# John Burkardt
#
import platform
print ( '' )
print ( 'r8ccs_test():' )
print ( ' Python version: ' + platform.python_version ( ) )
print ( ' Test the r8ccs() library for real matrices in CCS format.' )
r8ccs_dif2_test ( )
r8ccs_get_test ( )
r8ccs_ijk_test ( )
r8ccs_inc_test ( )
r8ccs_indicator_test ( )
r8ccs_kij_test ( )
r8ccs_mtv_test ( )
r8ccs_mv_test ( )
r8ccs_print_test ( )
r8ccs_print_some_test ( )
r8ccs_random_test ( )
r8ccs_read_test ( )
r8ccs_set_test ( )
r8ccs_to_r8ge_test ( )
r8ccs_write_test ( )
r8ccs_zeros_test ( )
return
def r8ccs_dif2 ( m, n, nz_num ):
#*****************************************************************************80
#
## r8ccs_dif2() sets the second difference as a CCS matrix.
#
# Discussion:
#
# CCS is the compressed column storage format.
#
# Associated with this format, we have an M by N matrix
# with NZ_NUM nonzero entries. We construct the column pointer
# vector COL of length N+1, such that entries of column J will be
# stored in positions COL(J) through COL(J+1)-1. This indexing
# refers to both the ROW and A vectors, which store the row indices
# and the values of the nonzero entries. The entries of the
# ROW vector corresponding to each column are assumed to be
# ascending sorted.
#
# The CCS format is equivalent to the MATLAB "sparse" format,
# and the Harwell Boeing "real unsymmetric assembled" (RUA) format.
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 02 October 2015
#
# Author:
#
# John Burkardt
#
# Input:
#
# integer M, the number of rows of the matrix.
#
# integer N, the number of columns of the matrix.
#
# integer NZ_NUM, the number of nonzero entries.
#
# Output:
#
# integer COLPTR(N+1), indicate where each column's data begins.
#
# integer ROW(NZ_NUM), the row indices.
#
# real A(NZ_NUM), the nonzero entries.
#
import numpy as np
#
# Column pointers
#
colptr = np.zeros ( n + 1, dtype = np.int32 )
colptr[0] = 0
colptr[1] = 2
for j in range ( 2, n ):
colptr[j] = colptr[j-1] + 3
colptr[n] = colptr[n-1] + 2
#
# Row indices
#
rowind = np.zeros ( nz_num, dtype = np.int32 )
k = 0
rowind[k] = 0
k = k + 1
rowind[k] = 1
k = k + 1
for j in range ( 1, n - 1 ):
for i in range ( j - 1, j + 2 ):
rowind[k] = i
k = k + 1
rowind[k] = m - 2
k = k + 1
rowind[k] = m - 1
k = k + 1
#
# Values
#
a = np.zeros ( nz_num, dtype = np.float64 )
k = 0
j = 0
i = 0
a[k] = 2.0
k = k + 1
i = 1
a[k] = -1.0
k = k + 1
for j in range ( 1, n - 1 ):
i = j - 1
a[k] = -1.0
k = k + 1
i = j
a[k] = 2.0
k = k + 1
i = j + 1
a[k] = -1.0
k = k + 1
j = n - 1
i = m - 2
a[k] = -1.0
k = k + 1
i = m - 1
a[k] = 2.0
k = k + 1
return colptr, rowind, a
def r8ccs_dif2_test ( ):
#*****************************************************************************80
#
## r8ccs_dif2_test() tests r8ccs_dif2().
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 02 October 2015
#
# Author:
#
# John Burkardt
#
import numpy as np
print ( '' )
print ( 'r8ccs_dif2_test():' )
print ( ' r8ccs_dif2() returns the second difference as a CCS matrix.' )
m = 5
n = 5
nz_num = 13
#
# Set the matrix.
#
colptr, rowind, a = r8ccs_dif2 ( m, n, nz_num )
#
# Print the matrix.
#
title = 'The second difference matrix in CCS format:'
r8ccs_print_header ( m, n, nz_num, colptr, rowind, title )
r8ccs_print ( m, n, nz_num, colptr, rowind, a, title )
return
def r8ccs_get ( m, n, nz_num, colptr, rowind, a, i, j ):
#*****************************************************************************80
#
## r8ccs_get() gets a value of a CCS matrix.
#
# Discussion:
#
# It is legal to request entries of the matrix for which no storage
# was set aside. In that case, a zero value will be returned.
#
# CCS is the compressed column storage format.
#
# Associated with this format, we have an M by N matrix
# with NZ_NUM nonzero entries. We construct the column pointer
# vector COL of length N+1, such that entries of column J will be
# stored in positions COL(J) through COL(J+1)-1. This indexing
# refers to both the ROW and A vectors, which store the row indices
# and the values of the nonzero entries. The entries of the
# ROW vector corresponding to each column are assumed to be
# ascending sorted.
#
# The CCS format is equivalent to the MATLAB "sparse" format,
# and the Harwell Boeing "real unsymmetric assembled" (RUA) format.
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 29 September 2015
#
# Author:
#
# John Burkardt
#
# Reference:
#
# Iain Duff, Roger Grimes, John Lewis,
# User's Guide for the Harwell-Boeing Sparse Matrix Collection,
# October 1992
#
# Input:
#
# integer M, the number of rows of the matrix.
#
# integer N, the number of columns of the matrix.
#
# integer NZ_NUM, the number of nonzero entries.
#
# integer COLPTR(N+1), indicate where each column's data begins.
#
# integer ROWIND(NZ_NUM), the row indices.
#
# real A(NZ_NUM), the nonzero entries.
#
# integer I, J, the indices of the value to retrieve.
#
# Output:
#
# real AIJ, the value of A(I,J).
#
#
# Seek sparse index K corresponding to full index (I,J).
#
k = r8ccs_ijk ( m, n, nz_num, colptr, rowind, i, j )
#
# If no K was found, then be merciful, and simply return 0.
#
if ( k == -1 ):
aij = 0.0
else:
aij = a[k]
return aij
def r8ccs_get_test ( ):
#*****************************************************************************80
#
## r8ccs_get_test() tests r8ccs_get().
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 03 October 2015
#
# Author:
#
# John Burkardt
#
import numpy as np
m = 5
n = 5
nz_num = 12
colptr = np.array ( [ 0, 3, 5, 7, 9, 12 ], dtype = np.int32 )
rowind = np.array ( [ 0, 1, 3, 0, 1, 2, 4, 3, 4, 0, 1, 4 ], dtype = np.int32 )
print ( '' )
print ( 'r8ccs_get_test():' )
print ( ' r8ccs_get() gets an entry of a matrix in the CCS format.' )
#
# Randomize the matrix.
#
a = r8ccs_random ( m, n, nz_num, colptr, rowind )
#
# Print the matrix.
#
title = 'A random CCS matrix:'
r8ccs_print_header ( m, n, nz_num, colptr, rowind, title )
r8ccs_print ( m, n, nz_num, colptr, rowind, a, title )
print ( '' )
print ( ' r8ccs_get() retrieves 10 entries.' )
print ( '' )
print ( ' I J K VALUE' )
print ( '' )
for test in range ( 0, 10 ):
k = np.random.randint ( 0, nz_num )
i, j = r8ccs_kij ( m, n, nz_num, colptr, rowind, k )
value = r8ccs_get ( m, n, nz_num, colptr, rowind, a, i, j )
print ( ' %8d %8d %8d %14.6g' % ( i, j, k, value ) )
return
def r8ccs_ijk ( m, n, nz_num, colptr, rowind, i, j ):
#*****************************************************************************80
#
## r8ccs_ijk() seeks K, the sparse index of (I,J), the full index of a CCS matrix.
#
# Discussion:
#
# CCS is the compressed column storage format.
#
# Associated with this format, we have an M by N matrix
# with NZ_NUM nonzero entries. We construct the column pointer
# vector COL of length N+1, such that entries of column J will be
# stored in positions COL(J) through COL(J+1)-1. This indexing
# refers to both the ROW and A vectors, which store the row indices
# and the values of the nonzero entries. The entries of the
# ROW vector corresponding to each column are assumed to be
# ascending sorted.
#
# The CCS format is equivalent to the MATLAB "sparse" format,
# and the Harwell Boeing "real unsymmetric assembled" (RUA) format.
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 29 September 2015
#
# Author:
#
# John Burkardt
#
# Reference:
#
# Iain Duff, Roger Grimes, John Lewis,
# User's Guide for the Harwell-Boeing Sparse Matrix Collection,
# October 1992
#
# Input:
#
# integer M, the number of rows of the matrix.
#
# integer N, the number of columns of the matrix.
#
# integer NZ_NUM, the number of nonzero entries.
#
# integer COLPTR(N+1), indicate where each column's data begins.
#
# integer ROWIND(NZ_NUM), the row indices.
#
# integer I, J, the indices of the value to retrieve.
#
# Output:
#
# integer K, the index of the sparse matrix in which entry
# (I,J) is stored, or -1 if no such entry exists.
#
import numpy as np
#
# Determine the part of ROW containing row indices of entries in column J.
#
k1 = colptr[j]
k2 = colptr[j+1] - 1
#
# Seek the location K for which ROWIND(K) = I.
#
rowsub = np.zeros ( k2 + 1 - k1 )
for i2 in range ( 0, k2 + 1 - k1 ):
rowsub[i2] = rowind[k1+i2]
k = i4vec_search_binary_a ( k2+1-k1, rowsub, i )
if ( k != -1 ):
k = k + k1
return k
def r8ccs_ijk_test ( ):
#*****************************************************************************80
#
## r8ccs_ijk_test() tests r8ccs_ijk().
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 03 October 2015
#
# Author:
#
# John Burkardt
#
import numpy as np
m = 5
n = 5
nz_num = 12
colptr = np.array ( [ 0, 3, 5, 7, 9, 12 ], dtype = np.int32 )
rowind = np.array ( [ 0, 1, 3, 0, 1, 2, 4, 3, 4, 0, 1, 4 ], dtype = np.int32 )
test_num = 20
print ( '' )
print ( 'r8ccs_ijk_test():' )
print ( ' r8ccs_ijk() gets K from (I,J) for a CCS matrix.' )
#
# Randomize the matrix.
#
a = r8ccs_random ( m, n, nz_num, colptr, rowind )
#
# Print the matrix.
#
title = 'A random CCS matrix:'
r8ccs_print_header ( m, n, nz_num, colptr, rowind, title )
r8ccs_print ( m, n, nz_num, colptr, rowind, a, title )
print ( '' )
print ( ' r8ccs_ijk() locates some (I,J) entries.' )
print ( '' )
print ( ' I J K' )
print ( '' )
for test in range ( 0, test_num ):
i = np.random.randint ( 0, m )
j = np.random.randint ( 0, n )
k = r8ccs_ijk ( m, n, nz_num, colptr, rowind, i, j )
print ( ' %8d %8d %8d' % ( i, j, k ) )
return
def r8ccs_inc ( m, n, nz_num, colptr, rowind, a, i, j, aij ):
#*****************************************************************************80
#
## r8ccs_inc() increments a value of a CCS matrix.
#
# Discussion:
#
# CCS is the compressed column storage format.
#
# Associated with this format, we have an M by N matrix
# with NZ_NUM nonzero entries. We construct the column pointer
# vector COL of length N+1, such that entries of column J will be
# stored in positions COL(J) through COL(J+1)-1. This indexing
# refers to both the ROW and A vectors, which store the row indices
# and the values of the nonzero entries. The entries of the
# ROW vector corresponding to each column are assumed to be
# ascending sorted.
#
# The CCS format is equivalent to the MATLAB "sparse" format,
# and the Harwell Boeing "real unsymmetric assembled" (RUA) format.
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 29 September 2015
#
# Author:
#
# John Burkardt
#
# Reference:
#
# Iain Duff, Roger Grimes, John Lewis,
# User's Guide for the Harwell-Boeing Sparse Matrix Collection,
# October 1992
#
# Input:
#
# integer M, the number of rows of the matrix.
#
# integer N, the number of columns of the matrix.
#
# integer NZ_NUM, the number of nonzero entries.
#
# integer COLPTR(N+1), indicate where each column's data begins.
#
# integer ROWIND(NZ_NUM), the row indices.
#
# real A(NZ_NUM), the nonzero entries.
#
# integer I, J, the indices of the value to retrieve.
#
# real AIJ, the value to be added to A(I,J).
#
# Output:
#
# real A(NZ_NUM), entry (I,J) has been incremented.
#
#
# Seek sparse index K corresponding to full index (I,J).
#
k = r8ccs_ijk ( m, n, nz_num, colptr, rowind, i, j )
#
# If no K was found, we fail.
#
if ( k == -1 ):
print ( '' )
print ( 'r8ccs_inc(): Fatal error!' )
print ( ' r8ccs_ijk() could not find the entry.' )
print ( ' Row I = ', i )
print ( ' Col J = ', j )
raise Exception ( 'r8ccs_inc(): Fatal error!' )
a[k] = a[k] + aij
return a
def r8ccs_inc_test ( ):
#*****************************************************************************80
#
## r8ccs_inc_test() tests r8ccs_inc().
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 03 October 2015
#
# Author:
#
# John Burkardt
#
import numpy as np
m = 5
n = 5
nz_num = 12
colptr = np.array ( [ 0, 3, 5, 7, 9, 12 ], dtype = np.int32 )
rowind = np.array ( [ 0, 1, 3, 0, 1, 2, 4, 3, 4, 0, 1, 4 ], dtype = np.int32 )
print ( '' )
print ( 'r8ccs_inc_test():' )
print ( ' r8ccs_inc() increments entries in a CCS matrix.' )
#
# Randomize the matrix.
#
a = r8ccs_random ( m, n, nz_num, colptr, rowind )
#
# Print the matrix.
#
title = 'A random CCS matrix:'
r8ccs_print_header ( m, n, nz_num, colptr, rowind, title )
r8ccs_print ( m, n, nz_num, colptr, rowind, a, title )
print ( '' )
print ( ' r8ccs_inc() increments 10 entries at random.' )
print ( '' )
print ( ' I J K NEW_VALUE' )
print ( '' )
for test in range ( 0, 10 ):
k = np.random.randint ( 0, nz_num )
i, j = r8ccs_kij ( m, n, nz_num, colptr, rowind, k )
value = 20.0 + test
a = r8ccs_inc ( m, n, nz_num, colptr, rowind, a, i, j, value )
value = r8ccs_get ( m, n, nz_num, colptr, rowind, a, i, j )
print ( ' %8d %8d %8d %14f' % ( i, j, k, value ) )
r8ccs_print ( m, n, nz_num, colptr, rowind, a, ' The final CCS matrix:' )
return
def r8ccs_indicator ( m, n, nz_num, colptr, rowind ):
#*****************************************************************************80
#
## r8ccs_indicator() sets up a CCS indicator matrix.
#
# Discussion:
#
# CCS is the compressed column storage format.
#
# Associated with this format, we have an M by N matrix
# with NZ_NUM nonzero entries. We construct the column pointer
# vector COL of length N+1, such that entries of column J will be
# stored in positions COL(J) through COL(J+1)-1. This indexing
# refers to both the ROW and A vectors, which store the row indices
# and the values of the nonzero entries. The entries of the
# ROW vector corresponding to each column are assumed to be
# ascending sorted.
#
# The CCS format is equivalent to the MATLAB "sparse" format,
# and the Harwell Boeing "real unsymmetric assembled" (RUA) format.
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 01 October 2015
#
# Author:
#
# John Burkardt
#
# Reference:
#
# Iain Duff, Roger Grimes, John Lewis,
# User's Guide for the Harwell-Boeing Sparse Matrix Collection,
# October 1992
#
# Input:
#
# integer M, the number of rows of the matrix.
#
# integer N, the number of columns of the matrix.
#
# integer NZ_NUM, the number of nonzero elements in A.
#
# integer COLPTR(N+1), points to the first element of each column.
#
# integer ROWIND(NZ_NUM), contains the row indices of the elements.
#
# Output:
#
# real A(NZ_NUM), the matrix.
#
import numpy as np
fac = 10 ** ( i4_log_10 ( n ) + 1 )
a = np.zeros ( nz_num )
for j in range ( 0, n ):
for k in range ( colptr[j], colptr[j+1] ):
i = rowind[k]
a[k] = fac * ( i + 1 ) + ( j + 1 )
return a
def r8ccs_indicator_test ( ):
#*****************************************************************************80
#
## r8ccs_indicator_test() tests r8ccs_indicator().
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 01 October 2015
#
# Author:
#
# John Burkardt
#
import numpy as np
m = 5
n = 5
nz_num = 12
colptr = np.array ( [ 0, 3, 5, 7, 9, 12 ] )
rowind = np.array ( [ 0, 1, 3, 0, 1, 2, 4, 3, 4, 0, 1, 4 ] )
print ( '' )
print ( 'r8ccs_indicator_test():' )
print ( ' r8ccs_indicator() sets up a CCS indicator matrix;' )
#
# Set the matrix.
#
a = r8ccs_indicator ( m, n, nz_num, colptr, rowind )
#
# Print the matrix.
#
title = 'An indicator matrix in CCS format:'
r8ccs_print_header ( m, n, nz_num, colptr, rowind, title )
r8ccs_print ( m, n, nz_num, colptr, rowind, a, title )
return
def r8ccs_kij ( m, n, nz_num, colptr, rowind, k ):
#*****************************************************************************80
#
## r8ccs_kij() seeks (I,J), the full index of K, the sparse index of a CCS matrix.
#
# Discussion:
#
# CCS is the compressed column storage format.
#
# Associated with this format, we have an M by N matrix
# with NZ_NUM nonzero entries. We construct the column pointer
# vector COL of length N+1, such that entries of column J will be
# stored in positions COL(J) through COL(J+1)-1. This indexing
# refers to both the ROW and A vectors, which store the row indices
# and the values of the nonzero entries. The entries of the
# ROW vector corresponding to each column are assumed to be
# ascending sorted.
#
# The CCS format is equivalent to the MATLAB "sparse" format,
# and the Harwell Boeing "real unsymmetric assembled" (RUA) format.
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 29 September 2015
#
# Author:
#
# John Burkardt
#
# Reference:
#
# Iain Duff, Roger Grimes, John Lewis,
# User's Guide for the Harwell-Boeing Sparse Matrix Collection,
# October 1992
#
# Input:
#
# integer M, the number of rows of the matrix.
#
# integer N, the number of columns of the matrix.
#
# integer NZ_NUM, the number of nonzero entries.
#
# integer COLPTR(N+1), indicate where each column's data begins.
#
# integer ROWIND(NZ_NUM), the row indices.
#
# integer K, the sparse index of an entry of the matrix.
# 1 <= K <= NZ_NUM.
#
# Output:
#
# integer I, J, the full indices corresponding to the sparse
# index K.
#
i = -1
j = -1
if ( k < 0 or nz_num <= k ):
return i, j
#
# The row index is easy.
#
i = rowind[k]
#
# Determine the column by bracketing in COLPTR.
#
for jj in range ( 0, n ):
k1 = colptr[jj]
k2 = colptr[jj+1] - 1
if ( k1 <= k and k <= k2 ):
j = jj
break
if ( j == -1 ):
return i, j
return i, j
def r8ccs_kij_test ( ):
#*****************************************************************************80
#
## r8ccs_kij_test() tests r8ccs_kij().
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 03 October 2015
#
# Author:
#
# John Burkardt
#
import numpy as np
m = 5
n = 5
nz_num = 12
colptr = np.array ( [ 0, 3, 5, 7, 9, 12 ] )
rowind = np.array ( [ 0, 1, 3, 0, 1, 2, 4, 3, 4, 0, 1, 4 ] )
print ( '' )
print ( 'r8ccs_kij_test():' )
print ( ' r8ccs_kij() gets (I,J) from K in a CCS matrix.' )
#
# Randomize the matrix.
#
a = r8ccs_random ( m, n, nz_num, colptr, rowind )
#
# Print the matrix.
#
title = 'A random CCS matrix:'
r8ccs_print_header ( m, n, nz_num, colptr, rowind, title )
r8ccs_print ( m, n, nz_num, colptr, rowind, a, title )
print ( '' )
print ( ' r8ccs_kij() locates some K entries.' )
print ( '' )
print ( ' K I J' )
print ( '' )
for test in range ( 0, 20 ):
k = np.random.randint ( 0, nz_num )
i, j = r8ccs_kij ( m, n, nz_num, colptr, rowind, k )
print ( ' %8d %8d %8d' % ( k, i, j ) )
return
def r8ccs_mtv ( m, n, nz_num, colptr, rowind, a, x ):
#*****************************************************************************80
#
## r8ccs_mtv() multiplies a vector times a CCS matrix.
#
# Discussion:
#
# CCS is the compressed column storage format.
#
# Associated with this format, we have an M by N matrix
# with NZ_NUM nonzero entries. We construct the column pointer
# vector COL of length N+1, such that entries of column J will be
# stored in positions COL(J) through COL(J+1)-1. This indexing
# refers to both the ROW and A vectors, which store the row indices
# and the values of the nonzero entries. The entries of the
# ROW vector corresponding to each column are assumed to be
# ascending sorted.
#
# The CCS format is equivalent to the MATLAB "sparse" format,
# and the Harwell Boeing "real unsymmetric assembled" (RUA) format.
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 29 September 2015
#
# Author:
#
# John Burkardt
#
# Reference:
#
# Iain Duff, Roger Grimes, John Lewis,
# User's Guide for the Harwell-Boeing Sparse Matrix Collection,
# October 1992
#
# Input:
#
# integer M, the number of rows of the matrix.
#
# integer N, the number of columns of the matrix.
#
# integer NZ_NUM, the number of nonzero elements in A.
#
# integer COLPTR(N+1), points to the first element of each column.
#
# integer ROWIND(NZ_NUM), contains the row indices of the elements.
#
# real A(NZ_NUM), the matrix.
#
# real X(M), the vector to be multiplied.
#
# Output:
#
# real B(N), the product A'*X.
#
import numpy as np
b = np.zeros ( n )
for j in range ( 0, n ):
for k in range ( colptr[j], colptr[j+1] ):
i = rowind[k]
b[j] = b[j] + a[k] * x[i]
return b
def r8ccs_mtv_test ( ):
#*****************************************************************************80
#
## r8ccs_mtv_test() tests r8ccs_mtv().
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 04 October 2015
#
# Author:
#
# John Burkardt
#
import numpy as np
m = 5
n = 5
nz_num = 12
colptr = np.array ( [ 0, 3, 5, 7, 9, 12 ] )
rowind = np.array ( [ 0, 1, 3, 0, 1, 2, 4, 3, 4, 0, 1, 4 ] )
print ( '' )
print ( 'r8ccs_mtv_test():' )
print ( ' r8ccs_mtv() computes b=A\'*x, where A is a CCS matrix.' )
#
# Randomize the matrix.
#
a = r8ccs_random ( m, n, nz_num, colptr, rowind )
#
# Print the matrix.
#
title = 'A random CCS matrix:'
r8ccs_print_header ( m, n, nz_num, colptr, rowind, title )
r8ccs_print ( m, n, nz_num, colptr, rowind, a, title )
x = np.zeros ( n )
x[0] = 1.0
x[n-1] = -1.0
print ( ' x:' )
print ( x )
b = r8ccs_mtv ( m, n, nz_num, colptr, rowind, a, x )
print ( ' b=A\'*x:' )
print ( b )
return
def r8ccs_mv ( m, n, nz_num, colptr, rowind, a, x ):
#*****************************************************************************80
#
## r8ccs_mv() multiplies a CCS matrix times a vector.
#
# Discussion:
#
# CCS is the compressed column storage format.
#
# Associated with this format, we have an M by N matrix
# with NZ_NUM nonzero entries. We construct the column pointer
# vector COL of length N+1, such that entries of column J will be
# stored in positions COL(J) through COL(J+1)-1. This indexing
# refers to both the ROW and A vectors, which store the row indices
# and the values of the nonzero entries. The entries of the
# ROW vector corresponding to each column are assumed to be
# ascending sorted.
#
# The CCS format is equivalent to the MATLAB "sparse" format,
# and the Harwell Boeing "real unsymmetric assembled" (RUA) format.
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 29 September 2015
#
# Author:
#
# John Burkardt
#
# Reference:
#
# Iain Duff, Roger Grimes, John Lewis,
# User's Guide for the Harwell-Boeing Sparse Matrix Collection,
# October 1992
#
# Input:
#
# integer M, the number of rows of the matrix.
#
# integer N, the number of columns of the matrix.
#
# integer NZ_NUM, the number of nonzero elements in A.
#
# integer COLPTR(N+1), points to the first element of each column.
#
# integer ROWIND(NZ_NUM), contains the row indices of the elements.
#
# real A(NZ_NUM), the matrix.
#
# real X(N), the vector to be multiplied.
#
# Output:
#
# real B(M), the product A*X.
#
import numpy as np
b = np.zeros ( m )
for j in range ( 0, n ):
for k in range ( colptr[j], colptr[j+1] ):
i = rowind[k]
b[i] = b[i] + a[k] * x[j]
return b
def r8ccs_mv_test ( ):
#*****************************************************************************80
#
## r8ccs_mv_test() tests r8ccs_mv().
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 04 October 2015
#
# Author:
#
# John Burkardt
#
import numpy as np
m = 5
n = 5
nz_num = 12
colptr = np.array ( [ 0, 3, 5, 7, 9, 12 ] )
rowind = np.array ( [ 0, 1, 3, 0, 1, 2, 4, 3, 4, 0, 1, 4 ] )
print ( '' )
print ( 'r8ccs_mv_test():' )
print ( ' r8ccs_mv() computes b=A*x, where A is a CCS matrix.' )
#
# Randomize the matrix.
#
a = r8ccs_random ( m, n, nz_num, colptr, rowind )
#
# Print the matrix.
#
title = 'A random CCS matrix:'
r8ccs_print_header ( m, n, nz_num, colptr, rowind, title )
r8ccs_print ( m, n, nz_num, colptr, rowind, a, title )
x = np.zeros ( n )
x[0] = 1.0
x[n-1] = -1.0
print ( ' x:' )
print ( x )
b = r8ccs_mv ( m, n, nz_num, colptr, rowind, a, x )
print ( ' b=A*x:' )
print ( b )
return
def r8ccs_print ( m, n, nz_num, colptr, rowind, a, title ):
#*****************************************************************************80
#
## r8ccs_print() prints a CCS matrix.
#
# Discussion:
#
# CCS is the compressed column storage format.
#
# Associated with this format, we have an M by N matrix
# with NZ_NUM nonzero entries. We construct the column pointer
# vector COL of length N+1, such that entries of column J will be
# stored in positions COL(J) through COL(J+1)-1. This indexing
# refers to both the ROW and A vectors, which store the row indices
# and the values of the nonzero entries. The entries of the
# ROW vector corresponding to each column are assumed to be
# ascending sorted.
#
# The CCS format is equivalent to the MATLAB "sparse" format,
# and the Harwell Boeing "real unsymmetric assembled" (RUA) format.
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 29 September 2015
#
# Author:
#
# John Burkardt
#
# Reference:
#
# Iain Duff, Roger Grimes, John Lewis,
# User's Guide for the Harwell-Boeing Sparse Matrix Collection,
# October 1992
#
# Input:
#
# integer M, the number of rows of the matrix.
#
# integer N, the number of columns of the matrix.
#
# integer NZ_NUM, the number of nonzero elements in A.
#
# integer COLPTR(N+1), points to the first element of each column.
#
# integer ROWIND(NZ_NUM), contains the row indices of the elements.
#
# real A(NZ_NUM), the matrix.
#
# string TITLE, a title.
#
r8ccs_print_some ( m, n, nz_num, colptr, rowind, a, 0, 0, m - 1, n - 1, title )
return
def r8ccs_print_test ( ):
#*****************************************************************************80
#
## r8ccs_print_test() tests r8ccs_print().
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 04 October 2015
#
# Author:
#
# John Burkardt
#
import numpy as np
m = 5
n = 5
nz_num = 12
colptr = np.array ( [ 0, 3, 5, 7, 9, 12 ] )
rowind = np.array ( [ 0, 1, 3, 0, 1, 2, 4, 3, 4, 0, 1, 4 ] )
print ( '' )
print ( 'r8ccs_print_test():' )
print ( ' r8ccs_print() prints a CCS matrix.' )
#
# Randomize the matrix.
#
a = r8ccs_random ( m, n, nz_num, colptr, rowind )
#
# Print the matrix.
#
title = 'A random CCS matrix:'
r8ccs_print_header ( m, n, nz_num, colptr, rowind, title )
r8ccs_print ( m, n, nz_num, colptr, rowind, a, title )
return
def r8ccs_print_header ( m, n, nz_num, colptr, rowind, title ):
#*****************************************************************************80
#
## r8ccs_print_header() prints the header of a CCS matrix.
#
# Discussion:
#
# CCS is the compressed column storage format.
#
# Associated with this format, we have an M by N matrix
# with NZ_NUM nonzero entries. We construct the column pointer
# vector COL of length N+1, such that entries of column J will be
# stored in positions COL(J) through COL(J+1)-1. This indexing
# refers to both the ROW and A vectors, which store the row indices
# and the values of the nonzero entries. The entries of the
# ROW vector corresponding to each column are assumed to be
# ascending sorted.
#
# The CCS format is equivalent to the MATLAB "sparse" format,
# and the Harwell Boeing "real unsymmetric assembled" (RUA) format.
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 22 June 2022
#
# Author:
#
# John Burkardt
#
# Reference:
#
# Iain Duff, Roger Grimes, John Lewis,
# User's Guide for the Harwell-Boeing Sparse Matrix Collection,
# October 1992
#
# Input:
#
# integer M, the number of rows of the matrix.
#
# integer N, the number of columns of the matrix.
#
# integer NZ_NUM, the number of nonzero elements in A.
#
# integer COLPTR(N+1), points to the first element of each column.
#
# integer ROWIND(NZ_NUM), contains the row indices of the elements.
#
# string TITLE, a title.
#
print ( '' )
print ( title )
print ( ' Matrix rows M = ', m )
print ( ' Matrix columns N = ', n )
print ( ' Matrix nonzeros NZ_NUM = ', nz_num )
print ( ' The COLPTR vector:' )
print ( colptr )
print ( ' The ROWIND vector:' )
print ( rowind )
return
def r8ccs_print_some ( m, n, nz_num, colptr, rowind, a, ilo, jlo, ihi, jhi, title ):
#*****************************************************************************80
#
## r8ccs_print_some() prints some of a CCS matrix.
#
# Discussion:
#
# CCS is the compressed column storage format.
#
# Associated with this format, we have an M by N matrix
# with NZ_NUM nonzero entries. We construct the column pointer
# vector COL of length N+1, such that entries of column J will be
# stored in positions COL(J) through COL(J+1)-1. This indexing
# refers to both the ROW and A vectors, which store the row indices
# and the values of the nonzero entries. The entries of the
# ROW vector corresponding to each column are assumed to be
# ascending sorted.
#
# The CCS format is equivalent to the MATLAB "sparse" format,
# and the Harwell Boeing "real unsymmetric assembled" (RUA) format.
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 29 September 2015
#
# Author:
#
# John Burkardt
#
# Reference:
#
# Iain Duff, Roger Grimes, John Lewis,
# User's Guide for the Harwell-Boeing Sparse Matrix Collection,
# October 1992
#
# Input:
#
# integer M, the number of rows of the matrix.
#
# integer N, the number of columns of the matrix.
#
# integer NZ_NUM, the number of nonzero elements in A.
#
# integer COLPTR(N+1), points to the first element of each column.
#
# integer ROWIND(NZ_NUM), contains the row indices of the elements.
#
# real A(NZ_NUM), the matrix.
#
# integer ILO, JLO, IHI, JHI, the first row and
# column, and the last row and column to be printed.
#
# string TITLE, a title.
#
incx = 5
print ( '' )
print ( title )
#
# Print the columns of the matrix, in strips of 5.
#
for j2lo in range ( jlo, jhi + 1, incx ):
j2hi = j2lo + incx - 1
j2hi = min ( j2hi, n - 1 )
j2hi = min ( j2hi, jhi )
inc = j2hi + 1 - j2lo
print ( '' )
print ( ' Col: ', end = '' )
for j in range ( j2lo, j2hi + 1 ):
print ( '%7d ' % ( j ), end = '' )
print ( '' )
print ( ' Row' )
print ( ' ---' )
#
# Determine the range of the rows in this strip.
#
i2lo = max ( ilo, 0 )
i2hi = min ( ihi, m - 1 )
for i in range ( i2lo, i2hi + 1 ):
print ( '%4d' % ( i ), end = '' )
#
# Print out (up to) 5 entries in row I, that lie in the current strip.
#
for j in range ( j2lo, j2hi + 1 ):
aij = 0.0
for k in range ( colptr[j], colptr[j+1] ):
if ( rowind[k] == i ):
aij = a[k]
break
print ( ' %12g' % ( aij ), end = '' )
print ( '' )
return
def r8ccs_print_some_test ( ):
#*****************************************************************************80
#
## r8ccs_print_some_test() tests r8ccs_print_some().
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 04 October 2015
#
# Author:
#
# John Burkardt
#
import numpy as np
m = 10
n = 10
nz_num = 28
colptr = np.array ( [ 0, 2, 5, 8, 11, 14, 17, 20, 23, 26, 28 ] )
rowind = np.array ( [ \
1, 2, \
1, 2, 3, \
2, 3, 4, \
3, 4, 5, \
4, 5, 6, \
5, 6, 7, \
6, 7, 8, \
7, 8, 9, \
8, 9, 10, \
9, 10 ] )
print ( '' )
print ( 'r8ccs_print_some_test():' )
print ( ' r8ccs_print_some() prints some of a CCS matrix.' )
#
# Set the matrix.
#
a = r8ccs_indicator ( m, n, nz_num, colptr, rowind )
#
# Print the matrix.
#
title = 'The indicator matrix in CCS format:'
r8ccs_print_header ( m, n, nz_num, colptr, rowind, title )
title = 'Rows 1-5, Cols 4-7:'
r8ccs_print_some ( m, n, nz_num, colptr, rowind, a, 1, 4, 5, 7, title )
return
def r8ccs_random ( m, n, nz_num, colptr, rowind ):
#*****************************************************************************80
#
## r8ccs_random() randomizes a CCS matrix.
#
# Discussion:
#
# CCS is the compressed column storage format.
#
# Associated with this format, we have an M by N matrix
# with NZ_NUM nonzero entries. We construct the column pointer
# vector COL of length N+1, such that entries of column J will be
# stored in positions COL(J) through COL(J+1)-1. This indexing
# refers to both the ROW and A vectors, which store the row indices
# and the values of the nonzero entries. The entries of the
# ROW vector corresponding to each column are assumed to be
# ascending sorted.
#
# The CCS format is equivalent to the MATLAB "sparse" format,
# and the Harwell Boeing "real unsymmetric assembled" (RUA) format.
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 29 September 2015
#
# Author:
#
# John Burkardt
#
# Reference:
#
# Iain Duff, Roger Grimes, John Lewis,
# User's Guide for the Harwell-Boeing Sparse Matrix Collection,
# October 1992
#
# Input:
#
# integer M, the number of rows of the matrix.
#
# integer N, the number of columns of the matrix.
#
# integer NZ_NUM, the number of nonzero elements in A.
#
# integer COLPTR(N+1), points to the first element of each column.
#
# integer ROWIND(NZ_NUM), contains the row indices of the elements.
#
# Output:
#
# real A(NZ_NUM), the matrix.
#
import numpy as np
a = np.random.normal ( loc = 0.0, scale = 1.0, size = nz_num )
return a
def r8ccs_random_test ( ):
#*****************************************************************************80
#
## r8ccs_random_test() tests r8ccs_random().
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 03 October 2015
#
# Author:
#
# John Burkardt
#
import numpy as np
m = 5
n = 5
nz_num = 12
colptr = np.array ( [ 0, 3, 5, 7, 9, 12 ] )
rowind = np.array ( [ 0, 1, 3, 0, 1, 2, 4, 3, 4, 0, 1, 4 ] )
print ( '' )
print ( 'r8ccs_random_test():' )
print ( ' r8ccs_random() randomizes a CCS matrix;' )
#
# Set the matrix.
#
a = r8ccs_random ( m, n, nz_num, colptr, rowind )
#
# Print the matrix.
#
title = 'A random matrix in CCS format:'
r8ccs_print_header ( m, n, nz_num, colptr, rowind, title )
r8ccs_print ( m, n, nz_num, colptr, rowind, a, title )
return
def r8ccs_read ( col_file, row_file, a_file, m, n, nz_num ):
#*****************************************************************************80
#
## r8ccs_read() reads a CCS matrix from three files.
#
# Discussion:
#
# This routine needs the values of M, N, and NZ_NUM, which can be
# determined by a call to r8ccs_read_size().
#
# CCS is the compressed column storage format.
#
# Associated with this format, we have an M by N matrix
# with NZ_NUM nonzero entries. We construct the column pointer
# vector COL of length N+1, such that entries of column J will be
# stored in positions COL(J) through COL(J+1)-1. This indexing
# refers to both the ROW and A vectors, which store the row indices
# and the values of the nonzero entries. The entries of the
# ROW vector corresponding to each column are assumed to be
# ascending sorted.
#
# The CCS format is equivalent to the MATLAB "sparse" format,
# and the Harwell Boeing "real unsymmetric assembled" (RUA) format.
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 29 September 2015
#
# Author:
#
# John Burkardt
#
# Reference:
#
# Iain Duff, Roger Grimes, John Lewis,
# User's Guide for the Harwell-Boeing Sparse Matrix Collection,
# October 1992
#
# Input:
#
# string COL_FILE, ROW_FILE, A_FILE, the names of the
# files containing the column pointers, row indices, and matrix entries.
#
# integer M, N, the number of rows and columns in the matrix.
#
# integer NZ_NUM, the number of nonzero elements in the matrix.
#
# Output:
#
# integer COLPTR(N+1), the column pointers.
#
# integer ROWIND(NZ_NUM), the row indices.
#
# real A(NZ_NUM), the nonzero elements
# of the matrix.
#
import numpy as np
colptr = np.zeros ( nz_num, dtype = np.int32 )
rowind = np.zeros ( nz_num, dtype = np.int32 )
a = np.zeros ( nz_num, dtype = np.float64 )
#
# Read the column information.
#
input_unit = open ( col_file, 'r' )
k = 0
for line in input_unit:
words = line.split ( )
colptr[k] = int ( words[0] )
k = k + 1
input_unit.close ( )
#
# Read the row information.
#
input_unit = open ( row_file, 'r' )
k = 0
for line in input_unit:
words = line.split ( )
rowind[k] = int ( words[0] )
k = k + 1
input_unit.close ( )
#
# Read the value information.
#
input_unit = open ( a_file, 'r' )
k = 0
for line in input_unit:
words = line.split ( )
a[k] = float ( words[0] )
k = k + 1
input_unit.close ( )
return colptr, rowind, a
def r8ccs_read_size ( col_file, row_file ):
#*****************************************************************************80
#
## r8ccs_read_size() reads the sizes of a CCS matrix from a file.
#
# Discussion:
#
# The value of M is "guessed" to be the largest value that occurs in
# the ROW file. However, if a row index of 0 is encountered, then
# the value of M is incremented by 1.
#
# The value of N is the number of records in the COL file minus 1.
#
# The value of NZ_NUM is simply the number of records in the ROW file.
#
# The value of BASE is 0 or 1, depending on whether the program
# "guesses" that the row and column indices are 0-based or 1-based.
# Although the first entry of the COL array might be used as evidence,
# this program makes its determination based on whether it encounters
# a 0 index in the ROW file.
#
# CCS is the compressed column storage format.
#
# Associated with this format, we have an M by N matrix
# with NZ_NUM nonzero entries. We construct the column pointer
# vector COL of length N+1, such that entries of column J will be
# stored in positions COL(J) through COL(J+1)-1. This indexing
# refers to both the ROW and A vectors, which store the row indices
# and the values of the nonzero entries. The entries of the
# ROW vector corresponding to each column are assumed to be
# ascending sorted.
#
# The CCS format is equivalent to the MATLAB "sparse" format,
# and the Harwell Boeing "real unsymmetric assembled" (RUA) format.
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 29 September 2015
#
# Author:
#
# John Burkardt
#
# Reference:
#
# Iain Duff, Roger Grimes, John Lewis,
# User's Guide for the Harwell-Boeing Sparse Matrix Collection,
# October 1992
#
# Input:
#
# string COL_FILE, ROW_FILE, the names of the
# column and row files that describe the structure of the matrix.
#
# Output:
#
# integer M, N, the inferred number of rows and columns
# in the sparse matrix.
#
# integer NZ_NUM, the number of nonzero entries in the
# sparse matrix.
#
# integer BASE, is 0 if the row indexing is believed
# to be 0-based, and 1 if the row-index is believed to be
# 1-based. In uncertain cases, BASE = 1 is the default.
#
#
# Default values.
#
m = -1
n = -1
nz_num = 0
base = 1
#
# Check the COL file first.
#
input_unit = open ( col_file, 'r' )
n = 0
for line in input_unit:
words = line.split ( )
k = int ( words[0] )
n = n + 1
n = n - 1
input_unit.close ( )
#
# Check the ROW file.
#
input_unit = open ( row_file, 'r' )
base = 1
m = 0
for line in input_unit:
words = line.split ( )
i = int ( words[0] )
nz_num = nz_num + 1
m = max ( m, i )
if ( i == 0 ):
base = 0
input_unit.close ( )
return m, n, nz_num, base
def r8ccs_read_test ( ):
#*****************************************************************************80
#
## r8ccs_read_test() tests r8ccs_read().
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 04 October 2015
#
# Author:
#
# John Burkardt
#
a_file = 'r8ccs_a.txt';
col_file = 'r8ccs_col.txt';
row_file = 'r8ccs_row.txt';
print ( '' )
print ( 'r8ccs_read_test():' )
print ( ' r8ccs_read() reads a CCS matrix from a file.' )
#
# Read the matrix sizes, and the matrix values.
#
m, n, nz_num, base = r8ccs_read_size ( col_file, row_file )
colptr, rowind, a = r8ccs_read ( col_file, row_file, a_file, m, n, nz_num )
if ( base == 1 ):
print ( '' )
print ( ' ROWIND and COLPTR indexing is 1-based.' )
colptr = colptr - 1
rowind = rowind - 1
print ( ' ROWIND and COLPTR indexing rebased at 0.' )
#
# Print the matrix.
#
title = 'The CCS matrix read from the file:'
r8ccs_print_header ( m, n, nz_num, colptr, rowind, title )
r8ccs_print ( m, n, nz_num, colptr, rowind, a, title )
return
def r8ccs_set ( m, n, nz_num, colptr, rowind, a, i, j, aij ):
#*****************************************************************************80
#
## r8ccs_set() sets a value of a CCS matrix.
#
# Discussion:
#
# CCS is the compressed column storage format.
#
# Associated with this format, we have an M by N matrix
# with NZ_NUM nonzero entries. We construct the column pointer
# vector COL of length N+1, such that entries of column J will be
# stored in positions COL(J) through COL(J+1)-1. This indexing
# refers to both the ROW and A vectors, which store the row indices
# and the values of the nonzero entries. The entries of the
# ROW vector corresponding to each column are assumed to be
# ascending sorted.
#
# The CCS format is equivalent to the MATLAB "sparse" format,
# and the Harwell Boeing "real unsymmetric assembled" (RUA) format.
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 29 September 2015
#
# Author:
#
# John Burkardt
#
# Reference:
#
# Iain Duff, Roger Grimes, John Lewis,
# User's Guide for the Harwell-Boeing Sparse Matrix Collection,
# October 1992
#
# Input:
#
# integer M, the number of rows of the matrix.
#
# integer N, the number of columns of the matrix.
#
# integer NZ_NUM, the number of nonzero entries.
#
# integer COLPTR(N+1), indicate where each column's data begins.
#
# integer ROWIND(NZ_NUM), the row indices.
#
# real A(NZ_NUM), the nonzero entries.
#
# integer I, J, the indices of the value to retrieve.
#
# real AIJ, the new value of A(I,J).
#
# Output:
#
# real A(NZ_NUM), the entry of A corresponding to (I,J)
# has been reset.
#
#
# Seek sparse index K corresponding to full index (I,J).
#
k = r8ccs_ijk ( m, n, nz_num, colptr, rowind, i, j )
#
# If no K was found, we fail.
#
if ( k == -1 ):
print ( '' )
print ( 'r8ccs_set(): Fatal error!' )
print ( ' r8ccs_ijk() could not find the entry.' )
print ( ' Row I = ', i )
print ( ' Col J = ', j )
raise Exception ( 'r8ccs_set(): Fatal error!' )
a[k] = aij
return a
def r8ccs_set_test ( ):
#*****************************************************************************80
#
## r8ccs_set_test() tests r8ccs_set().
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 04 October 2015
#
# Author:
#
# John Burkardt
#
import numpy as np
m = 5
n = 5
nz_num = 12
colptr = np.array ( [ 0, 3, 5, 7, 9, 12 ] )
rowind = np.array ( [ 0, 1, 3, 0, 1, 2, 4, 3, 4, 0, 1, 4 ] )
print ( '' )
print ( 'r8ccs_set_test():' )
print ( ' r8ccs_set() sets entries in a CCS matrix' )
#
# Randomize the matrix.
#
a = r8ccs_random ( m, n, nz_num, colptr, rowind )
#
# Print the matrix.
#
title = 'A random CCS matrix:'
r8ccs_print_header ( m, n, nz_num, colptr, rowind, title )
r8ccs_print ( m, n, nz_num, colptr, rowind, a, title )
print ( '' )
print ( ' r8ccs_set() sets 10 entries at random.' )
print ( '' )
print ( ' I J K NEW_VALUE' )
print ( '' )
for test in range ( 0, 10 ):
k = np.random.randint ( 0, nz_num )
i, j = r8ccs_kij ( m, n, nz_num, colptr, rowind, k )
value = 100.0 + test
a = r8ccs_set ( m, n, nz_num, colptr, rowind, a, i, j, value )
print ( ' %8d %8d %8d %14f' % ( i, j, k, value ) )
r8ccs_print ( m, n, nz_num, colptr, rowind, a, ' The final CCS matrix:' )
return
def r8ccs_to_r8ge ( m, n, nz_num, colptr, rowind, a ):
#*****************************************************************************80
#
## r8ccs_to_r8ge() converts a CCS matrix to GE format.
#
# Discussion:
#
# CCS is the compressed column storage format.
#
# Associated with this format, we have an M by N matrix
# with NZ_NUM nonzero entries. We construct the column pointer
# vector COL of length N+1, such that entries of column J will be
# stored in positions COL(J) through COL(J+1)-1. This indexing
# refers to both the ROW and A vectors, which store the row indices
# and the values of the nonzero entries. The entries of the
# ROW vector corresponding to each column are assumed to be
# ascending sorted.
#
# The CCS format is equivalent to the MATLAB "sparse" format,
# and the Harwell Boeing "real unsymmetric assembled" (RUA) format.
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 29 September 2015
#
# Author:
#
# John Burkardt
#
# Reference:
#
# Iain Duff, Roger Grimes, John Lewis,
# User's Guide for the Harwell-Boeing Sparse Matrix Collection,
# October 1992
#
# Input:
#
# integer M, the number of rows of the matrix.
#
# integer N, the number of columns of the matrix.
#
# integer NZ_NUM, the number of nonzero elements in A.
#
# integer COLPTR(N+1), points to the first element of each column.
#
# integer ROWIND(NZ_NUM), contains the row indices of the elements.
#
# real A(NZ_NUM), the CCS matrix.
#
# Output:
#
# real B(M,N), the R8GE matrix.
#
import numpy as np
b = np.zeros ( [ m, n ] )
for j in range ( 0, n ):
for k in range ( colptr[j], colptr[j+1] ):
i = rowind[k]
b[i,j] = a[k]
return b
def r8ccs_to_r8ge_test ( ):
#*****************************************************************************80
#
## r8ccs_to_r8ge_test() tests r8ccs_to_r8ge().
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 04 October 2015
#
# Author:
#
# John Burkardt
#
import numpy as np
m = 5
n = 5
nz_num = 12
colptr = np.array ( [ 0, 3, 5, 7, 9, 12 ] )
rowind = np.array ( [ 0, 1, 3, 0, 1, 2, 4, 3, 4, 0, 1, 4 ] )
print ( '' )
print ( 'r8ccs_to_r8ge_test():' )
print ( ' r8ccs_to_r8ge() converts a matrix from CCS to GE format;' )
#
# Set the matrix.
#
a_r8cc = r8ccs_indicator ( m, n, nz_num, colptr, rowind )
#
# Print the matrix.
#
title = 'An indicator matrix in CCS format:'
r8ccs_print_header ( m, n, nz_num, colptr, rowind, title )
r8ccs_print ( m, n, nz_num, colptr, rowind, a_r8cc, title )
#
# Make a GE copy of the CCS matrix.
#
a_r8ge = r8ccs_to_r8ge ( m, n, nz_num, colptr, rowind, a_r8cc )
print ( ' The GE matrix:' )
print ( a_r8ge )
return
def r8ccs_write ( col_file, row_file, a_file, m, n, nz_num, colptr, rowind, a ):
#*****************************************************************************80
#
## r8ccs_write() writes a CCS matrix to three files.
#
# Discussion:
#
# CCS is the compressed column storage format.
#
# Associated with this format, we have an M by N matrix
# with NZ_NUM nonzero entries. We construct the column pointer
# vector COL of length N+1, such that entries of column J will be
# stored in positions COL(J) through COL(J+1)-1. This indexing
# refers to both the ROW and A vectors, which store the row indices
# and the values of the nonzero entries. The entries of the
# ROW vector corresponding to each column are assumed to be
# ascending sorted.
#
# The CCS format is equivalent to the MATLAB "sparse" format,
# and the Harwell Boeing "real unsymmetric assembled" (RUA) format.
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 29 September 2015
#
# Author:
#
# John Burkardt
#
# Reference:
#
# Iain Duff, Roger Grimes, John Lewis,
# User's Guide for the Harwell-Boeing Sparse Matrix Collection,
# October 1992
#
# Input:
#
# string COL_FILE, ROW_FILE, A_FILE, the names of the
# files containing the column pointers, row entries, and matrix entries.
#
# integer M, N, the number of rows and columns in the matrix.
#
# integer NZ_NUM, the number of nonzero elements in the matrix.
#
# integer COLPTR(N+1), the column pointers.
#
# integer ROWIND(NZ_NUM), the row indices.
#
# real A(NZ_NUM), the nonzero elements of the matrix.
#
#
# Write the column information.
#
output_unit = open ( col_file, 'w' )
for k in range ( 0, n + 1 ):
s = '%d\n' % ( colptr[k] )
output_unit.write ( s )
output_unit.close ( )
#
# Write the row information.
#
output_unit = open ( row_file, 'w' )
for k in range ( 0, nz_num ):
s = '%d\n' % ( rowind[k] )
output_unit.write ( s )
output_unit.close ( )
#
# Write the value information.
#
output_unit = open ( a_file, 'w' )
for k in range ( 0, nz_num ):
s = '%d\n' % ( a[k] )
output_unit.write ( s )
output_unit.close ( )
return
def r8ccs_write_test ( ):
#*****************************************************************************80
#
## r8ccs_write_test() tests r8ccs_write().
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 04 October 2015
#
# Author:
#
# John Burkardt
#
import numpy as np
m = 5
n = 5
nz_num = 12
a_file = 'r8ccs_a.txt'
colptr = np.array ( [ 0, 3, 5, 7, 9, 12 ] )
col_file = 'r8ccs_col.txt'
rowind = np.array ( [ 0, 1, 3, 0, 1, 2, 4, 3, 4, 0, 1, 4 ] )
row_file = 'r8ccs_row.txt'
print ( '' )
print ( 'r8ccs_write_test():' )
print ( ' r8ccs_write() writes a CCS matrix to 3 files' )
#
# Set the matrix.
#
a = r8ccs_indicator ( m, n, nz_num, colptr, rowind )
#
# Print the matrix.
#
title = 'An indicator matrix in CCS format:'
r8ccs_print_header ( m, n, nz_num, colptr, rowind, title )
r8ccs_print ( m, n, nz_num, colptr, rowind, a, title )
#
# Write the matrix to files.
#
r8ccs_write ( col_file, row_file, a_file, m, n, nz_num, colptr, rowind, a )
return
def r8ccs_zeros ( m, n, nz_num, colptr, rowind ):
#*****************************************************************************80
#
## r8ccs_zeros() zeros a CCS matrix.
#
# Discussion:
#
# CCS is the compressed column storage format.
#
# Associated with this format, we have an M by N matrix
# with NZ_NUM nonzero entries. We construct the column pointer
# vector COL of length N+1, such that entries of column J will be
# stored in positions COL(J) through COL(J+1)-1. This indexing
# refers to both the ROW and A vectors, which store the row indices
# and the values of the nonzero entries. The entries of the
# ROW vector corresponding to each column are assumed to be
# ascending sorted.
#
# The CCS format is equivalent to the MATLAB "sparse" format,
# and the Harwell Boeing "real unsymmetric assembled" (RUA) format.
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 29 September 2015
#
# Author:
#
# John Burkardt
#
# Reference:
#
# Iain Duff, Roger Grimes, John Lewis,
# User's Guide for the Harwell-Boeing Sparse Matrix Collection,
# October 1992
#
# Input:
#
# integer M, the number of rows of the matrix.
#
# integer N, the number of columns of the matrix.
#
# integer NZ_NUM, the number of nonzero elements in A.
#
# integer COLPTR(N+1), points to the first element of each column.
#
# integer ROWIND(NZ_NUM), contains the row indices of the elements.
#
# Output:
#
# real A(NZ_NUM), the matrix.
#
import numpy as np
a = np.zeros ( nz_num )
return a
def r8ccs_zeros_test ( ):
#*****************************************************************************80
#
## r8ccs_zeros_test() tests r8ccs_zeros().
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 03 October 2015
#
# Author:
#
# John Burkardt
#
import numpy as np
m = 5
n = 5
nz_num = 12
colptr = np.array ( [ 0, 3, 5, 7, 9, 12 ] )
rowind = np.array ( [ 0, 1, 3, 0, 1, 2, 4, 3, 4, 0, 1, 4 ] )
print ( '' )
print ( 'r8ccs_zeros_test():' )
print ( ' r8ccs_zeros() zeros a CCS matrix;' )
#
# Set the matrix.
#
a = r8ccs_zeros ( m, n, nz_num, colptr, rowind )
#
# Print the matrix.
#
title = 'Zero matrix in CCS format:'
r8ccs_print_header ( m, n, nz_num, colptr, rowind, title )
r8ccs_print ( m, n, nz_num, colptr, rowind, a, title )
return
def timestamp ( ):
#*****************************************************************************80
#
## timestamp() prints the date as a timestamp.
#
# Licensing:
#
# This code is distributed under the MIT license.
#
# Modified:
#
# 21 August 2019
#
# Author:
#
# John Burkardt
#
import time
t = time.time ( )
print ( time.ctime ( t ) )
return
if ( __name__ == '__main__' ):
timestamp ( )
r8ccs_test ( )
timestamp ( )