sparse Tutorial 03: Linear algebra support.#

Contributors:
    Mauricio Aristizabal
        University of Texas at San Antonio
        Universidad EAFIT, Colombia.

Date of creation:  Jun 19 2024
Last modification: Jun 19 2024

Introduction#

The main purpose of the Order Truncated Imaginary (OTI) numbers was to provide an efficient method for computing arbitrary-order derivatives in Finite Element applications. Finite Element analysis is highly dependent on linear algbera operations, thus it was fundamental to develop an efficient interface to linear algebra, including vector and matrix operations.

This document shows applications of pyoti for the computation of linear-algebra related operations.

For this, first import the libraries required.

[1]:
import pyoti.sparse as oti  # Import the library.
import numpy as np

Array class.#

The support for linear algebra operations is provided by the definition of an array class. The array class is named ‘MATrix of Sparse Otis’ matso. As its name suggest, its fundamental mode of operation is as a matrix element. Thus this class always consider the behavior as a rank-2 array.

Although this class has its own constructor, pyoti provides a numpy-like array creation functions that allow an easy migration to the new class.

An array that can hold an OTI number can be created in multiple ways. Here will be discussed a few of them.

Filled constructors#

An array can be created using equivalent zeros or ones function, like with numpy. For instance, creating a 2x2 array of zeroes is created using the function zeros(shape,**kwargs), with shape being an integer or a list with up-to 2 elements defining the shape of the array.

[2]:
X = oti.zeros((2,2))
print(X)
matso< shape: (2, 2),
 - Column 0
(0,0) 0
(1,0) 0
 - Column 1
(0,1) 0
(1,1) 0
>

Remarks:

  • The number of dimensions of a matso is always 2.

  • Printing a matso array gives its shape and the contents are listed by columns.

  • It may be more convenient to print only the real part of the matso array. The real part can be extracted with the “attribute” real. This returns a numpy array with the real coefficients of the matso array.

[3]:
print(X.real)
print(type(X.real))
[[0. 0.]
 [0. 0.]]
<class 'numpy.ndarray'>

Other constructors are shown below:

Ones constructor:#

[4]:
print('"ones" constructor:')
print( oti.ones( ( 2, 3 ) ).real )
"ones" constructor:
[[1. 1. 1.]
 [1. 1. 1.]]

Identity matrix (eye)#

[5]:
print('"eye" constructor:')
print( oti.eye( 2 ).real  )
"eye" constructor:
[[1. 0.]
 [0. 1.]]

Array constructor#

Similar to numpy, pyoti has an array function that allows to create arrays from lists, numpy arrays, and other oti arrays. The array( val, ) is used for this matter, as follows:

Array from list#

[6]:
Y = oti.array([[1,2,3],
               [4,5,6]])
Y
[6]:
matso< shape: (2, 3),
 - Column 0
(0,0) 1
(1,0) 4
 - Column 1
(0,1) 2
(1,1) 5
 - Column 2
(0,2) 3
(1,2) 6
>

List items can also contain OTI scalars:

[7]:
Z = oti.array([[1+oti.e(1),2,3+oti.e(3)],
               [4+oti.e(2),5,6         ]])
Z
[7]:
matso< shape: (2, 3),
 - Column 0
(0,0) 1 + 1 * e([1])
(1,0) 4 + 1 * e([2])
 - Column 1
(0,1) 2
(1,1) 5
 - Column 2
(0,2) 3 + 1 * e([3])
(1,2) 6
>

Numpy arrays can also be used to construct the OTI array.

[8]:
npArray =  np.arange(0,20,2).reshape((5,-1))
print(npArray)

W = oti.array( npArray  )
W
[[ 0  2]
 [ 4  6]
 [ 8 10]
 [12 14]
 [16 18]]
[8]:
matso< shape: (5, 2),
 - Column 0
(0,0) 0
(1,0) 4
(2,0) 8
(3,0) 12
(4,0) 16
 - Column 1
(0,1) 2
(1,1) 6
(2,1) 10
(3,1) 14
(4,1) 18
>

Getter, setters and slicing operations#

Getting and setting elements in OTI arrays work similar to numpy arrays and lists. use the getter setter as follows

Note: All OTI arrays are rank-2 arrays: Thus, accessing individual elements must be done using two-indices. Indices are zero-based.

Setting an individual item for an OTI array:#

An Considering an OTI matrix \(X^*\in\mathbb{OTI}^{ n_r \times n_c }\), where \(n_r\) and \(n_c\) correspond to number of rows and number of columns of the matrix, an item \(x_{ij}^*\) is set as X[i,j]=x_ij.

For example, consider the matrix \(X^*\in\mathbb{OTI}^{ 2 \times 2 }\), defined as follows:

\[\begin{split}X^* = \begin{bmatrix} 0 & 1.5+\epsilon_2 \\ 2.5+\epsilon_1 & 0 \end{bmatrix}\end{split}\]

This matrix can be created by

[9]:
# Define zeros.
X = oti.zeros((2,2))

# Assign non-zeros.
X[0,1] = 1.5+oti.e(1)
X[1,0] = 2.5+oti.e(1)

X
[9]:
matso< shape: (2, 2),
 - Column 0
(0,0) 0
(1,0) 2.5 + 1 * e([1])
 - Column 1
(0,1) 1.5 + 1 * e([1])
(1,1) 0
>

Slicing#

Slicing uses the same notation as in numpy arrays.

[10]:
X[0,:]
[10]:
matso< shape: (1, 2),
 - Column 0
(0,0) 0
 - Column 1
(0,1) 1.5 + 1 * e([1])
>

Scalar operations#

Scalar operations with OTI arrays are supported. This means that OTI arrays support addition/subtraction/multiplication and division with OTI scalars. The opration is simply propagated elementwise.

[11]:
2*X
[11]:
matso< shape: (2, 2),
 - Column 0
(0,0) 0
(1,0) 5 + 2 * e([1])
 - Column 1
(0,1) 3 + 2 * e([1])
(1,1) 0
>
[12]:
X + oti.e(2)
[12]:
matso< shape: (2, 2),
 - Column 0
(0,0)  + 1 * e([2])
(1,0) 2.5 + 1 * e([1]) + 1 * e([2])
 - Column 1
(0,1) 1.5 + 1 * e([1]) + 1 * e([2])
(1,1)  + 1 * e([2])
>

Vectorized evaluations#

Similar to numpy arrays, OTI arrays allow elementwise operations. The way this works is that the operation, say a function evaluation, is performed elementwise.

[13]:
oti.cos(X)
[13]:
matso< shape: (2, 2),
 - Column 0
(0,0) 1
(1,0) -0.801144 - 0.598472 * e([1])
 - Column 1
(0,1) 0.0707372 - 0.997495 * e([1])
(1,1) 1
>

Linear algebra operations#

Matrix matrix multiplication#

[14]:
X = oti.array([[ 4.5, 2.8],
               [-2.8, 3.0]])

Y = oti.array([[ 0.5, 4.8, 3.9 ],
               [-0.8, 0.0, 4.8 ] ])

oti.dot(X,Y)
[14]:
matso< shape: (2, 3),
 - Column 0
(0,0) 0.01
(1,0) -3.8
 - Column 1
(0,1) 21.6
(1,1) -13.44
 - Column 2
(0,2) 30.99
(1,2) 3.48
>