dolfinx.la

Linear algebra functionality

Functions

orthonormalize(basis)

Orthogoalise set of PETSc vectors in-place

is_orthonormal(basis[, eps])

Check that list of PETSc vectors are orthonormal

matrix_csr(sp[, block_mode, dtype])

Create a distributed sparse matrix.

vector(map[, bs, dtype])

Create a distributed vector.

create_petsc_vector(map, bs)

Create a distributed PETSc vector.

Classes

MatrixCSR(A)

A distributed sparse matrix that uses compressed sparse row storage.

Vector(x)

A distributed vector object.

class dolfinx.la.InsertMode(self: dolfinx.cpp.la.InsertMode, value: int)

Bases: pybind11_object

Members:

add

insert

add = <InsertMode.add: 0>
insert = <InsertMode.insert: 1>
property name
property value
class dolfinx.la.MatrixCSR(A)[source]

Bases: object

A distributed sparse matrix that uses compressed sparse row storage.

Parameters:

A – The C++/pybind11 matrix object.

Note

Objects of this type should be created using matrix_csr() and not created using the class initialiser.

add(x, rows, cols, bs=1) None[source]

Add a block of values in the matrix.

property block_size

Block sizes for the matrix.

property data: ndarray[Any, dtype[floating]]

Underlying matrix entry data.

index_map(i: int) IndexMap[source]

Index map for row/column.

Parameters:

i – 0 for row map, 1 for column map.

property indices: ndarray[Any, dtype[int32]]

Local column indices.

property indptr: ndarray[Any, dtype[int64]]

Local row pointers.

scatter_reverse() None[source]

Scatter and accumulate ghost values.

set(x, rows, cols, bs=1) None[source]

Set a block of values in the matrix.

set_value(x: floating) None[source]

Set all non-zero entries to a value.

Parameters:

x – The value to set all non-zero entries to.

squared_norm() floating[source]

Compute the squared Frobenius norm.

Note

This operation is collective and requires communication.

to_dense() ndarray[Any, dtype[floating]][source]

Copy to a dense 2D array.

Note

Typically used for debugging.

to_scipy(ghosted=False)[source]

Convert to a SciPy CSR/BSR matrix. Data is shared.

SciPy must be available.

Parameters:

ghosted – If True rows that are ghosted in parallel are included in the return SciPy matrix, otherwise ghost rows are not included.

Returns:

SciPy compressed sparse row (both block sizes equal to one) or a SciPy block compressed sparse row matrix.

class dolfinx.la.Norm(self: dolfinx.cpp.la.Norm, value: int)

Bases: pybind11_object

Members:

l1

l2

linf

frobenius

frobenius = <Norm.frobenius: 3>
l1 = <Norm.l1: 0>
l2 = <Norm.l2: 1>
linf = <Norm.linf: 2>
property name
property value
class dolfinx.la.Vector(x)[source]

Bases: object

A distributed vector object.

Parameters:
  • map – Index map the describes the size and distribution of the vector

  • bs – Block size

Note

Objects of this type should be created using vector() and not created using the class initialiser.

property array: ndarray

Local representation of the vector.

property block_size: int

Block size for the vector.

property index_map: IndexMap

Index map that describes size and parallel distribution.

norm(type=<Norm.l2: 1>) floating[source]

Compute a norm of the vector.

Parameters:

type – Norm type to computed.

Returns:

Computed norm.

scatter_forward() None[source]

Update ghost entries.

scatter_reverse(mode: InsertMode) None[source]

Scatter ghost entries to owner.

Parameters:

mode – Control how scattered values are set/accumulated by owner.

dolfinx.la.create_petsc_vector(map, bs: int)[source]

Create a distributed PETSc vector.

Parameters:
  • map – Index map that describes the size and parallel layout of the vector to create.

  • bs – Block size of the vector.

dolfinx.la.is_orthonormal(basis, eps: float = 1e-12) bool[source]

Check that list of PETSc vectors are orthonormal

dolfinx.la.matrix_csr(sp, block_mode=<BlockMode.compact: 0>, dtype=<class 'numpy.float64'>) MatrixCSR[source]

Create a distributed sparse matrix.

The matrix uses compressed sparse row storage.

Parameters:
  • sp – The sparsity pattern that defines the nonzero structure of the matrix the parallel distribution of the matrix.

  • dtype – The scalar type.

Returns:

A sparse matrix.

dolfinx.la.orthonormalize(basis)[source]

Orthogoalise set of PETSc vectors in-place

dolfinx.la.vector(map, bs=1, dtype=<class 'numpy.float64'>) Vector[source]

Create a distributed vector.

Parameters:
  • map – Index map the describes the size and distribution of the vector.

  • bs – Block size.

  • dtype – The scalar type.

Returns:

A distributed vector.