# Heterogeneous blocked matrices

**URL:** <https://discourse.julialang.org/t/heterogeneous-blocked-matrices/2130>\
**Category:** Numerics\
**Tags:** question\
**Created:** [February 15, 2017, 10:39pm UTC](https://discourse.julialang.org/t/heterogeneous-blocked-matrices/2130 "2017-02-15T22:39:25Z")\
**Posts on this page:** 4\
**Page:** 1

<div class="post-metadata">

**Author:** ![dmbates](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dmbates/32/44_2.png) [@dmbates](https://discourse.julialang.org/u/dmbates)\
**Post date:** [February 15, 2017, 10:39pm UTC](https://discourse.julialang.org/t/heterogeneous-blocked-matrices/2130/1 "2017-02-15T22:39:25Z")

</div>

I am struggling to come up with a sufficiently specific method signature to handle the Cholesky decomposition of heterogeneous blocked matrices, by which I mean that I have a square matrix of `AbstractMatrix{Float64}`, say, but the lower-level matrices can be of different storage types: dense, sparse, diagonal, etc.

I can define the union of all the possible storage types

```nohighlight
typealias Block{T} Union{Diagonal{T}, Diagonal{StaticArrays.MMatrix},
    Hermitian{T, Matrix{T}}, LowerTriangular{T, Matrix{T}}, Matrix{T}}

```

but I don’t think `AbstractMatrix{Block{T}}` is meaningful in a method signature. Or at least I can’t seem to get it to work.

Currently I am patterning my code after the (unreleased) `LinearAlgebra` package by @andreasnoack. I want to create a method for his `cholBlocked!` function that will apply to my blocked matrices or square views of them and is more specific than `AbstractMatrix`. I need to somehow characterize `AbstractMatrix{AbstractMatrix{Float64}}`.

The gory details are in [https://github.com/dmbates/MixedModels.jl/tree/LowerCholesky](https://github.com/dmbates/MixedModels.jl/tree/LowerCholesky) Look for the string `cholBlocked!` in the src tree.

---

<div class="post-metadata">

**Author:** ![mbauman](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mbauman/32/31082_2.png) [@mbauman](https://discourse.julialang.org/u/mbauman)\
**Post date:** [February 15, 2017, 11:10pm UTC](https://discourse.julialang.org/t/heterogeneous-blocked-matrices/2130/2 "2017-02-15T23:10:40Z")

</div>

If I’m understanding your post correctly, I think you’re just butting up against parametric invariance. You just need to do `f{T<:Block}(B::AbstractMatrix{T})`. A concrete example here is `Diagonal{Matrix{Float64}}`:

```julia
julia> D = Diagonal([[1 2; 3 4], [5 6; 7 8]])
2×2 Diagonal{Array{Int64,2}}:
 [1 2; 3 4] ⋅
     ⋅ [5 6; 7 8]

julia> isblock{T<:Block}(x::AbstractArray{T}) = true
       isblock(x) = false
isblock (generic function with 2 methods)

julia> isblock(D)
true

```

---

<div class="post-metadata">

**Author:** ![dmbates](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dmbates/32/44_2.png) [@dmbates](https://discourse.julialang.org/u/dmbates)\
**Post date:** [February 16, 2017, 3:56pm UTC](https://discourse.julialang.org/t/heterogeneous-blocked-matrices/2130/3 "2017-02-16T15:56:11Z")

</div>

Thanks for the reply, Matt, but I don’t think that will work.

I should have been more specific in my explanation. It is the blocks that are heterogeneous, not the wrapper. The array of blocks is dense, square and lower triangular. Blocks on the diagonal are square and, by design, decreasing in size. The largest block, in the [1,1] position, is always diagonal. The trailing diagonal blocks are dense but often their size is negligible compared to the [1,1] block.

If I understand correctly, a signature like

```julia
f{T <: Block}(x::AbstactArray{T}) end

```

applies to a homogeneous array of blocks of a single type, `T`, which is one of the types in the `Block` union.

For the time being I will forgo trying to emulate the elegant recursive `cholBlocked!` function in @andreasnoack package and instead use explicit loops.

---

<div class="post-metadata">

**Author:** ![mbauman](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mbauman/32/31082_2.png) [@mbauman](https://discourse.julialang.org/u/mbauman)\
**Post date:** [February 16, 2017, 5:15pm UTC](https://discourse.julialang.org/t/heterogeneous-blocked-matrices/2130/4 "2017-02-16T17:15:49Z")

</div>

Ah, I think I understand. It’ll depend upon how you construct the blocked matrix. If you just do something like this:

```julia
julia> A = reshape([Diagonal([1.,2.,3.,4.]), zeros(2,4), zeros(4,2), Diagonal([5.,6.])], 2, 2)
2×2 Array{AbstractArray{Float64,2},2}:
 [1.0 0.0 0.0 0.0; 0.0 2.0 0.0 0.0; 0.0 0.0 3.0 0.0; 0.0 0.0 0.0 4.0] [0.0 0.0; 0.0 0.0; 0.0 0.0; 0.0 0.0]
 [0.0 0.0 0.0 0.0; 0.0 0.0 0.0 0.0] [5.0 0.0; 0.0 6.0]

```

Then the eltype is widened beyond your `Union` to simply be a `AbstractArray{Float64,2}`. Depending upon how you do it, the eltype might even get widened the whole way to simply `Any`. So instead of using `Block`, you something like `isblock{T<:AbstractMatrix}(A::AbstractMatrix{T}) = true` should do the trick.