# Broadcast raising a matrix to a specific power
**URL:** https://discourse.julialang.org/t/broadcast-raising-a-matrix-to-a-specific-power/96159
**Category:** General Usage
**Tags:** question
**Created:** [March 16, 2023, 12:42am UTC](https://discourse.julialang.org/t/broadcast-raising-a-matrix-to-a-specific-power/96159 "2023-03-16T00:42:03Z")
**Posts on this page:** 3
**Page:** 1
### Author:  [@Justin\_Booth](https://discourse.julialang.org/u/Justin_Booth)
#### Post date: [March 16, 2023, 12:42am UTC](https://discourse.julialang.org/t/broadcast-raising-a-matrix-to-a-specific-power/96159/1 "2023-03-16T00:42:03Z")
Hi, new to Julia here. I’m trying to model the accumulation of mutations in a population using a Markov chain, I’ve got a transition matrix and an initial state vector, basically what I want to do is raise the transition matrix to the power of the number of generations and then multiply it by the initial state vector and plot the values on a graph.
My code is as follows:
```julia
using Plots
μ=0.00005
T=[μ^2-2μ+1 μ-μ^2 μ-μ^2 μ^2; 0 1-μ 0 μ; 0 0 1-μ μ; 0 0 0 1]
i=[1000 0 0 0]
x = range(0, 10, length=100)
p_every = i*(Ref(T).^x)
y = p_every[1]
plot(x, y;dpi=460)
```
When I run this though, I get the following error:
DimensionMismatch(“matrix A has dimensions (1,4), vector B has length 100”)
1. **generic\_matvecmul!** (::Array{Array{Float64,2},1}, ::Char, ::Array{Int64,2}, ::Array{Array{Float64,2},1}, ::LinearAlgebra.MulAddMul{true,true,Bool,Bool})@_matmul.jl:639_
2. **mul!** @_matmul.jl:81_ [inlined]
3. **mul!** @_matmul.jl:208_ [inlined]
4. \*\*\*\*\*(::Array{Int64,2}, ::Array{Array{Float64,2},1})@_matmul.jl:51_
5. **top-level scope** @\*[Local: 3]
I’ve been having a hard time finding information about the broadcast function online, so I was hoping someone might be able to clarify what I’m doing wrong.
---
### Author:  [@Krastanov](https://discourse.julialang.org/u/Krastanov)
#### Post date: [March 16, 2023, 1:13am UTC](https://discourse.julialang.org/t/broadcast-raising-a-matrix-to-a-specific-power/96159/2 "2023-03-16T01:13:35Z")
Broadcasting “peels” one layer of containers and broadcasts over it (kinda like implicit list/array/set/etc comprehension). You use it corectly to broadcast over the `Ref` container and the `x` vector, i.e. for each element in the Ref and for each element in the vector, perform the `^` operation. The result of the `Ref * vector` will be a `vector`. Then you ask for the `*` multiplication between the `i 1xn matrix` and the result of the previous operation, which happens to be a vector. There is no multiplication defined between `1xn matrix` and `vector`. You probably want to broadcast put `i` in another container and broadcast over that as well as `Ref(i) .* Ref(T).^x`.
One source of mistakes here is that people usually think of `.*` as elementwise multiplication. You want to actually broadcast matrix multiplication over containers that themselves contain matrices.
On a side note, matrix-vector multiplication is much faster than matrix-matrix multiplication (quadratic vs qubic time complexity). You probably will get better results with something that repeatedly does matrix-vector multiplication:
```julia
p_every = [i] # maybe add a type to this
for iter in 1:iterations
push!(p_every, T*p_every[end])
end
```
or if you need only the first component
```julia
y = zeros(iterations)
buffer = copy(i)
for iter in 1:iterations
buffer = T*buffer # this can be done more efficiently with in-place ops
y[iter] = buffer[1]
end
```
or you can use the ~~`reduce`~~ `foldr/foldl` function or a few other ways, depending on what is important for your use case.
Two more notes:
`y=p_every[1]` would just give you the initial i, not the first scalar in each vector (which is what I assume you want)
`range(0,10,length=100)` would include fractions, not only integers. Fractional power of a matrix would be quite a bit more expensive. Make sure you actually need it (and if you actually need it, it would be more efficient to do some factorization of the matrix first, to make the computation cheaper)
---
### Author:  [@Justin\_Booth](https://discourse.julialang.org/u/Justin_Booth)
#### Post date: [March 16, 2023, 1:17am UTC](https://discourse.julialang.org/t/broadcast-raising-a-matrix-to-a-specific-power/96159/3 "2023-03-16T01:17:31Z")
Thank you so much, this is super helpful!
Yeah the efficiency of this is not something I’m at all considering because 1. this is basically my first time using julia, and 2. I just want to get a quick sketch of this problem, but I really do appreciate the insight and it’s certainly something I’ll look more into in the future.
Again thanks!