P573 Scientific Computing Lecture 5 Matrix Operations 1

1
P573Scientific ComputingLecture 5 - Matrix
Operations 1

• Peter Gottschling
• pgottsch_at_cs.indiana.edu
• www.osl.iu.edu/pgottsch/courses/p573-06

2
Overview

• Definitions
• Addition
• Norms
• Matrix vector product
• Naïve matrix product
• Alternative implementations
• Blocking

3
Matrix

• Field of m rows and n columns
• If mn called square matrix
• Otherwise rectangular

4
Submatrices

• Vector can be considered as m?1 matrix
• Matlab does this (everything is a matrix in
  matlab)
• Scalar values are 1?1 matrices
• However, vectors are special cases of matrices
  concerning representation
• Mathematically usually different meaning
• Submatrices gained by deleting rows and columns
  of matrix

5
Partitioning

• Splitting matrix into submatrices
• Partitioned matrix sometimes called block matrix
• In most cases, submatrices of block matrix have
  same size

6
Special Matrices

• An n?n matrix is called
• Diagonal if aij 0 for all i ? j
• Upper triangular if aij 0 for all i gt j
• Lower triangular if aij 0 for all i lt j
• A partitioned matrix is called
• Upper block triangular matrix if Aij 0 for all
  i gt j

7
More Special Matrices

• The transposed AT of an m?n matrix is
• n?m matrix where aij and aji are swapped
• Correspondingly the transposed conjugate matrix
  AH
• A is called symmetric if AAT
• A is called hermitian if AAH

8
Matrix Addition

• Matrices are added element-wise
• Dimensions must be equal
• Is this operation faster or slower then vector
  addition?
• Is it better to iterate first over the rows or
  first over the columns?

9
Run Time of Matrix Addition

• Adding matrices of dimensions m?n takes
• 2 m?n loads
• m?n stores
• m?n additions
• Recall that the additions can be in 1 or 2
  processor cycles
• Loading and storing is more expansive, especially
  from main memory
• But also from L2 (or L3 if exist)

10
Is Matrix or Vector Addition Faster?

• Adding vectors of length l m?n takes
• 2 l loads
• l stores
• l additions
• Thus same arithmetic and memory operations are
  performed
• Should run at the same speed, unless
• See next slides

11
Row Major and Column Major

• In C and C 2D arrays are stored row major
• All elements of one row are stored together
• In Fortran 2D arrays are stored column major
• All elements of one column are stored together
• If you define your own matrix type with your own
  indexing you use either one
• What is done in some libraries, like ???

12
Loop Nesting

• Loop nesting should fit storing scheme
• For row major matrices iterating first over rows
• For column major first over columns
• Stride-1 memory access
• Why is this important?

13
Matrix Norms

• Matrices are operators on vectors
• Matrix vector product Ax can be considered as
  applying the operator A on x
• Therefore most matrix norms are operator norms
• Definitions of vector norms be applied to the
  image of the operation

14
Special Matrix Norms

• .p for p 1, 2, and ? are the most
  interesting
• p 2 difficult to compute
• p 1 and ? much simpler
• The first is the maximum of column sums
• The latter is the maximal row sum

15
Matrix Vector Product

• Number of columns must be equal to vector
  dimension
• Result is a vector
• Dimension is number of rows
• i-th element of result is dot product of i-th row
  and multiplied vector

16
Alternative Computation

• Dot product is loop over matrix rows
• What corresponds to a loop over columns of the
  matrix?
• Scaled vector addition
• Is this faster or slower than dot product?

17
Run Time Analysis

• Operand vector n loads
• Result vector m stores
• Matrix mn loads
• Multiplications mn
• Additions m(n-1)
• Which is the most expensive?

18
Run Time Analysis (ii)

• Number of additions, multiplications and loads
  for the matrix are almost equal
• Loading from memory is much slower then basic
  arithmetic
• Assume we have a row-major matrix and a cache
  line of sc double
• Loads are composes of
• Vector addition mn memory loads
• Dot product 1/sc mn memory loads and (sc-1)/sc
  mn loads from cache
• Safe to assume tl,m tl, ? Sp ? sc
• For column-major matrices it is the opposite