Combines one or more LinearOperators in to a Block Diagonal matrix.

Inherits From: LinearOperator, Module

Used in the notebooks

Used in the tutorials

This operator combines one or more linear operators [op1,...,opJ], building a new LinearOperator, whose underlying matrix representation has each operator opi on the main diagonal, and zero's elsewhere.

Shape compatibility

If opj acts like a [batch] matrix Aj, then op_combined acts like the [batch] matrix formed by having each matrix Aj on the main diagonal.

Each opj is required to represent a matrix, and hence will have shape batch_shape_j + [M_j, N_j].

If opj has shape batch_shape_j + [M_j, N_j], then the combined operator has shape broadcast_batch_shape + [sum M_j, sum N_j], where broadcast_batch_shape is the mutual broadcast of batch_shape_j, j = 1,...,J, assuming the intermediate batch shapes broadcast.

Arguments to matmul, matvec, solve, and solvevec may either be single Tensors or lists of Tensors that are interpreted as blocks. The jth element of a blockwise list of Tensors must have dimensions that match opj for the given method. If a list of blocks is input, then a list of blocks is returned as well.

When the opj are not guaranteed to be square, this operator's methods might fail due to the combined operator not being square and/or lack of efficient methods.

# Create a 4 x 4 linear operator combined of two 2 x 2 operators.
operator_1 = LinearOperatorFullMatrix([[1., 2.], [3., 4.]])
operator_2 = LinearOperatorFullMatrix([[1., 0.], [0., 1.]])
operator = LinearOperatorBlockDiag([operator_1, operator_2])

==> [[1., 2., 0., 0.],
     [3., 4., 0., 0.],
     [0., 0., 1., 0.],
     [0., 0., 0., 1.]]

==> [4, 4]

==> scalar Tensor

x1 = ... # Shape [2, 2] Tensor
x2 = ... # Shape [2, 2] Tensor
x = tf.concat([x1, x2], 0)  # Shape [2, 4] Tensor
==> tf.concat([operator_1.matmul(x1), operator_2.matmul(x2)])

# Create a 5 x 4 linear operator combining three blocks.
operator_1 = LinearOperatorFullMatrix([[1.], [3.]])
operator_2 = LinearOperatorFullMatrix([[1., 6.]])
operator_3 = LinearOperatorFullMatrix([[2.], [7.]])
operator = LinearOperatorBlockDiag([operator_1, operator_2, operator_3])

==> [[1., 0., 0., 0.],
     [3., 0., 0., 0.],
     [0., 1., 6., 0.],
     [0., 0., 0., 2.]]
     [0., 0., 0., 7.]]

==> [5, 4]

# Create a [2, 3] batch of 4 x 4 linear operators.
matrix_44 = tf.random.normal(shape=[2, 3, 4, 4])
operator_44 = LinearOperatorFullMatrix(matrix)

# Create a [1, 3] batch of 5 x 5 linear operators.
matrix_55 = tf.random.normal(shape=[1, 3, 5, 5])
operator_55 = LinearOperatorFullMatrix(matrix_55)

# Combine to create a [2, 3] batch of 9 x 9 operators.
operator_99 = LinearOperatorBlockDiag([operator_44, operator_55])

# Create a shape [2, 3, 9] vector.
x = tf.random.normal(shape=[2, 3, 9])
==> Shape [2, 3, 9] Tensor

# Create a blockwise list of vectors.
x = [tf.random.normal(shape=[2, 3, 4]), tf.random.normal(shape=[2, 3, 5])]
==> [Shape [2, 3, 4] Tensor, Shape [2, 3, 5] Tensor]


The performance of LinearOperatorBlockDiag on any operation is equal to the sum of the individual operators' operations.

Matrix property hints

This LinearOperator is initialized with boolean flags of the form is_X, for X = non_singular, self_adjoint, positive_definite, square. These have the following meaning:

  • If is_X == True, callers should expect the operator to have the property X. This is a promise that should be fulfilled, but is not a runtime assert. For example, finite floating point precision may result in these promises being violated.
  • If is_X == False, callers should expect the operator to not have X.
  • If is_X == None (the default), callers should have no expectation either way.

operators Iterable of LinearOperator objects, each with the same dtype and composable shape.
is_non_singular Expect that this operator is non-singular.
is_self_adjoint Expect that this operator is equal to its hermitian transpose.
is_positive_definite Expect that this operator is positive definite, meaning the quadratic form x^H A x has positive real part for all nonzero x. Note that we do not require the operator to be self-adjoint to be positive-definite. See:
is_square Expect that this operator acts like square [batch] matrices. This is true by default, and will raise a ValueError otherwise.
name A name for this LinearOperator. Default is the individual operators names joined with _o_.

TypeError If all operators do not have the same dtype.
ValueError If operators is empty or are non-square.

H Returns the adjoint of the current LinearOperator.

Given A representing this LinearOperator, return A*. Note that calling self.adjoint() and self.H are equivalent.

batch_shape TensorShape of batch dimensions of this LinearOperator.

If this operator acts like the batch matrix A with A.shape = [B1,...,Bb, M, N], then this retu