# Exponential matrix calculation with the exp() function vs procedure with the use of the eigen() function

**URL:** <https://discourse.julialang.org/t/exponential-matrix-calculation-with-the-exp-function-vs-procedure-with-the-use-of-the-eigen-function/48798>\
**Category:** New to Julia\
**Tags:** question, linearalgebra\
**Created:** [October 22, 2020, 5:28am UTC](https://discourse.julialang.org/t/exponential-matrix-calculation-with-the-exp-function-vs-procedure-with-the-use-of-the-eigen-function/48798 "2020-10-22T05:28:49Z")\
**Posts on this page:** 8\
**Page:** 1

<div class="post-metadata">

**Author:** ![HerAdri](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/heradri/32/5816_2.png) [@HerAdri](https://discourse.julialang.org/u/HerAdri)\
**Post date:** [October 22, 2020, 5:28am UTC](https://discourse.julialang.org/t/exponential-matrix-calculation-with-the-exp-function-vs-procedure-with-the-use-of-the-eigen-function/48798/1 "2020-10-22T05:28:49Z")

</div>

Using the procedure shown in [Matrix Exponentiation](https://brilliant.org/wiki/matrix-exponentiation/?), I have performed the calculation of the exponential matrix:  
for the case of the first matrix, the result of the procedure matches the result returned by the function exp ():

```julia
08:22:59->>2×2 Array{Int64,2}:
 1 -1
 2 4
Julia>exp(A)
08:22:59->>2×2 Array{Float64,2}:
 -5.30742 -12.6965
 25.393 32.782
Julia>F=eigen(A)
08:22:59->>Eigen{Float64,Float64,Array{Float64,2},Array{Float64,1}}
values:
2-element Array{Float64,1}:
 2.0
 3.0
vectors:
2×2 Array{Float64,2}:
 -0.707107 0.447214
  0.707107 -0.894427
Julia>v=F.vectors
08:22:59->>2×2 Array{Float64,2}:
 -0.707107 0.447214
  0.707107 -0.894427
Julia>Λ=exp.(F.values)
08:22:59->>2-element Array{Float64,1}:
  7.38905609893065
 20.085536923187668
Julia>D=diagm(Λ)
08:22:59->>2×2 Array{Float64,2}:
 7.38906 0.0
 0.0 20.0855
Julia>v*D*inv(v)
08:22:59->>2×2 Array{Float64,2}:
 -5.30742 -12.6965
 25.393 32.782

```

But for the second matrix, the results do not match !!  
I would like to know why and what I should do to obtain the expected result

```julia
Julia>C=[-0.4 -1;1 0.45]
08:27:43->>2×2 Array{Float64,2}:
 -0.4 -1.0
  1.0 0.45
Julia>exp(C)
08:27:43->>2×2 Array{Float64,2}:
 0.254525 -0.890921
 0.890921 1.01181
Julia>F=eigen(C)
08:27:43->>Eigen{Complex{Float64},Complex{Float64},Array{Complex{Float64},2},Array{Complex{Float64},1}}
values:
2-element Array{Complex{Float64},1}:
 0.02499999999999991 - 0.9051933495115836im
 0.02499999999999991 + 0.9051933495115836im
vectors:
2×2 Array{Complex{Float64},2}:
 -0.30052-0.640068im -0.30052+0.640068im
 0.707107-0.0im 0.707107+0.0im
Julia>v=F.vectors
08:27:43->>2×2 Array{Complex{Float64},2}:
 -0.30052-0.640068im -0.30052+0.640068im
 0.707107-0.0im 0.707107+0.0im
Julia>Λ=exp.(F.values)
08:27:43->>2-element Array{Complex{Float64},1}:
 0.6331664487902933 - 0.8064560400308953im
 0.6331664487902933 + 0.8064560400308953im
Julia>D=diagm(Λ)
08:27:43->>2×2 Array{Complex{Float64},2}:
 0.633166-0.806456im 0.0+0.0im
      0.0+0.0im 0.633166+0.806456im
Julia>v*D*inv(v)
08:27:43->>2×2 Array{Complex{Float64},2}:
 0.254525+0.0im -0.890921+5.55112e-17im
 0.890921-5.55112e-17im 1.01181-1.11022e-16im

```

---

<div class="post-metadata">

**Author:** ![liuyxpp](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/liuyxpp/32/9870_2.png) [@liuyxpp](https://discourse.julialang.org/u/liuyxpp)\
**Post date:** [October 22, 2020, 5:34am UTC](https://discourse.julialang.org/t/exponential-matrix-calculation-with-the-exp-function-vs-procedure-with-the-use-of-the-eigen-function/48798/2 "2020-10-22T05:34:51Z")

</div>

It seems both are correct. Only the later one has imaginary parts which is about rounding error.

---

<div class="post-metadata">

**Author:** ![HerAdri](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/heradri/32/5816_2.png) [@HerAdri](https://discourse.julialang.org/u/HerAdri)\
**Post date:** [October 22, 2020, 6:26am UTC](https://discourse.julialang.org/t/exponential-matrix-calculation-with-the-exp-function-vs-procedure-with-the-use-of-the-eigen-function/48798/3 "2020-10-22T06:26:32Z")

</div>

```julia
Julia>exp(C)
09:22:30->>2×2 Array{Float64,2}:
 0.254525 -0.890921
 0.890921 1.01181

Julia>r=v*D*inv(v)
09:22:38->>2×2 Array{Complex{Float64},2}:
 0.254525+0.0im -0.890921+5.55112e-17im
 0.890921-5.55112e-17im 1.01181-1.11022e-16im

Julia>abs.(r)
09:22:38->>2×2 Array{Float64,2}:
 0.254525 0.890921
 0.890921 1.01181

Julia>abs.(r)==exp(C)
09:22:38->>false

```

The difference is in the sign of the element X [1,2]

---

<div class="post-metadata">

**Author:** ![liuyxpp](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/liuyxpp/32/9870_2.png) [@liuyxpp](https://discourse.julialang.org/u/liuyxpp)\
**Post date:** [October 22, 2020, 6:34am UTC](https://discourse.julialang.org/t/exponential-matrix-calculation-with-the-exp-function-vs-procedure-with-the-use-of-the-eigen-function/48798/4 "2020-10-22T06:34:32Z")

</div>

I think you really need is

```julia
real.(r)

```

not

```julia
abs.(r)

```

---

<div class="post-metadata">

**Author:** ![HerAdri](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/heradri/32/5816_2.png) [@HerAdri](https://discourse.julialang.org/u/HerAdri)\
**Post date:** [October 22, 2020, 7:40am UTC](https://discourse.julialang.org/t/exponential-matrix-calculation-with-the-exp-function-vs-procedure-with-the-use-of-the-eigen-function/48798/5 "2020-10-22T07:40:49Z")

</div>

> [@liuyxpp](#):
>
> `real.(r)`

using r = real. (…) I can plot the result:

```julia
C=[-0.4 -1;1 0.45]
x0=[0 0.22]
∆t = .25
F=eigen(C*∆t)
v=F.vectors
Λ=exp.(F.values)
D=diagm(Λ)
r=v*D*inv(v)
T=real.(r)
x = x0
x1 =x0 
for i = 1:100
x = x*T #repeatedly multiply by T
x1=vcat(x1, x) # & store current x1(t) in the array x1
end
p1=plot( x1)
p2=plot(x1[:,1],x1[:,2])
plot(p1, p2,layout = (2, 1), legend = false)

```

![plot_Real](https://global.discourse-cdn.com/julialang/original/3X/0/2/02cd1e58e0075a798b50cb85d08193df444933e2.jpeg)

how to plot the results using it using imaginary numbers?:

```julia
C=[-0.4 -1;1 0.45]
x0=[0 0.22]
∆t = .25
F=eigen(C*∆t)
v=F.vectors
Λ=exp.(F.values)
D=diagm(Λ)
r=v*D*inv(v)
x = x0
x1 =x0 
for i = 1:100
x = x*r #repeatedly multiply by T
x1=vcat(x1, x) # & store current x1(t) in the array x1
end
plot(x1)
plot(x1[:,1],x1[:,2])

```

![plot_Real Imag](https://global.discourse-cdn.com/julialang/original/3X/1/9/1925a81ec8e720d532592e52b3c391d535f79ed8.jpeg)

```julia
plot(x1[:,1],x1[:,2])
ERROR: MethodError: no method matching isless(::Float64, ::Complex{Float64})
Closest candidates are:
  isless(::Float64, ::Float64) at float.jl:465
  isless(::Missing, ::Any) at missing.jl:87
  isless(::AbstractFloat, ::AbstractFloat) at operators.jl:165
  ...
Stacktrace:
 [1] min(::Complex{Float64}, ::Float64) at .\operators.jl:431
 [2] expand_extrema! at C:\Users\hermesr\.julia\packages\Plots\qZHsp\src\axes.jl:335 [inlined]
 [3] expand_extrema!(::Plots.Axis, ::Array{Complex{Float64},1}) at C:\Users\hermesr\.julia\packages\Plots\qZHsp\src\axes.jl:358
 [4] expand_extrema!(::Plots.Subplot{Plots.GRBackend}, ::Dict{Symbol,Any}) at C:\Users\hermesr\.julia\packages\Plots\qZHsp\src\axes.jl:388
 [5] _expand_subplot_extrema(::Plots.Subplot{Plots.GRBackend}, ::Dict{Symbol,Any}, ::Symbol) at C:\Users\hermesr\.julia\packages\Plots\qZHsp\src\pipeline.jl:366
 [6] _process_seriesrecipe(::Plots.Plot{Plots.GRBackend}, ::Dict{Symbol,Any}) at C:\Users\hermesr\.julia\packages\Plots\qZHsp\src\pipeline.jl:402
 [7] _plot!(::Plots.Plot{Plots.GRBackend}, ::Dict{Symbol,Any}, ::Tuple{Array{Complex{Float64},1},Array{Complex{Float64},1}}) at C:\Users\hermesr\.julia\packages\Plots\qZHsp\src\plot.jl:234
 [8] plot(::Array{Complex{Float64},1}, ::Vararg{Array{Complex{Float64},1},N} where N; kw::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}) at C:\Users\hermesr\.julia\packages\Plots\qZHsp\src\plot.jl:57
 [9] plot(::Array{Complex{Float64},1}, ::Array{Complex{Float64},1}) at C:\Users\hermesr\.julia\packages\Plots\qZHsp\src\plot.jl:51
 [10] top-level scope at REPL[154]:1

```

---

<div class="post-metadata">

**Author:** ![mzaffalon](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mzaffalon/32/214168_2.png) [@mzaffalon](https://discourse.julialang.org/u/mzaffalon)\
**Post date:** [October 22, 2020, 8:44am UTC](https://discourse.julialang.org/t/exponential-matrix-calculation-with-the-exp-function-vs-procedure-with-the-use-of-the-eigen-function/48798/6 "2020-10-22T08:44:41Z")

</div>

Maybe `plot(real.(x1), imag.(x1))`?

---

<div class="post-metadata">

**Author:** ![HerAdri](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/heradri/32/5816_2.png) [@HerAdri](https://discourse.julialang.org/u/HerAdri)\
**Post date:** [October 22, 2020, 9:16am UTC](https://discourse.julialang.org/t/exponential-matrix-calculation-with-the-exp-function-vs-procedure-with-the-use-of-the-eigen-function/48798/7 "2020-10-22T09:16:45Z")

</div>

> [@mzaffalon](#):
>
> plot(real.(x1), imag.(x1))

```julia
p3=plot(real.(x1))
p4=plot(real.(x1[:,1]),real.(x1[:,2]))
p5=plot(real.(x1), imag.(x1))
plot(p3, p4,p5,layout = (3, 1), legend = false)

```

 ![plot_Real 2](https://global.discourse-cdn.com/julialang/original/3X/7/8/782d6e10f43648faa51f1a5f3cf2155ce0e89e76.jpeg)

---

<div class="post-metadata">

**Author:** ![Tamas\_Papp](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamas_papp/32/25949_2.png) [@Tamas\_Papp](https://discourse.julialang.org/u/Tamas_Papp)\
**Post date:** [October 22, 2020, 12:11pm UTC](https://discourse.julialang.org/t/exponential-matrix-calculation-with-the-exp-function-vs-procedure-with-the-use-of-the-eigen-function/48798/8 "2020-10-22T12:11:02Z")

</div>

> [@HerAdri](#):
>
> Using the procedure shown in [Matrix Exponentiation](https://brilliant.org/wiki/matrix-exponentiation/?), I have performed the calculation of the exponential matrix

This is well-known to be a difficult problem numerically. Please read

[https://www.researchgate.net/publication/228590312\_Nineteen\_Dubious\_Ways\_to\_Compute\_the\_Exponential\_of\_a\_Matrix\_Twenty-Five\_Years\_Later](https://www.researchgate.net/publication/228590312_Nineteen_Dubious_Ways_to_Compute_the_Exponential_of_a_Matrix_Twenty-Five_Years_Later)

where Section 6 talks about the method you linked. TL;DR: just use `exp`.