In [1]:
## The goal is to look at powers of Jordan blocks. 
# (I need a better name than Jordan blocks!)


function jordan_block(lam, k)
  n = k
  A = zeros(n,n)
  for i = 1:n
    A[i,i] = lam
    if i < n
      A[i,i+1] = 1
    end
  end
  return A
end 

A = jordan_block(0.99, 30) # This is a 5x5 matrix with 1s on the superdiagonal.
Out[1]:
30×30 Matrix{Float64}:
 0.99  1.0   0.0   0.0   0.0   0.0   …  0.0   0.0   0.0   0.0   0.0   0.0
 0.0   0.99  1.0   0.0   0.0   0.0      0.0   0.0   0.0   0.0   0.0   0.0
 0.0   0.0   0.99  1.0   0.0   0.0      0.0   0.0   0.0   0.0   0.0   0.0
 0.0   0.0   0.0   0.99  1.0   0.0      0.0   0.0   0.0   0.0   0.0   0.0
 0.0   0.0   0.0   0.0   0.99  1.0      0.0   0.0   0.0   0.0   0.0   0.0
 0.0   0.0   0.0   0.0   0.0   0.99  …  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.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   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.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
 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.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   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.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      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.0   0.0
 0.0   0.0   0.0   0.0   0.0   0.0      1.0   0.0   0.0   0.0   0.0   0.0
 0.0   0.0   0.0   0.0   0.0   0.0      0.99  1.0   0.0   0.0   0.0   0.0
 0.0   0.0   0.0   0.0   0.0   0.0   …  0.0   0.99  1.0   0.0   0.0   0.0
 0.0   0.0   0.0   0.0   0.0   0.0      0.0   0.0   0.99  1.0   0.0   0.0
 0.0   0.0   0.0   0.0   0.0   0.0      0.0   0.0   0.0   0.99  1.0   0.0
 0.0   0.0   0.0   0.0   0.0   0.0      0.0   0.0   0.0   0.0   0.99  1.0
 0.0   0.0   0.0   0.0   0.0   0.0      0.0   0.0   0.0   0.0   0.0   0.99
In [2]:
##
using LinearAlgebra
norms_jordan = map(k -> opnorm(A^k, 1), 1:100000)
norms_diag = map(k->opnorm(Diagonal([0.99 for i = 1:5])^k, 1), 1:1000)
Out[2]:
1000-element Vector{Float64}:
 0.99
 0.9801
 0.9702989999999999
 0.96059601
 0.9509900498999999
 0.941480149401
 0.9320653479069899
 0.9227446944279201
 0.9135172474836408
 0.9043820750088044
 0.8953382542587164
 0.8863848717161292
 0.8775210229989678
 ⋮
 4.8217807298316476e-5
 4.7735629225333305e-5
 4.725827293307997e-5
 4.6785690203749175e-5
 4.631783330171168e-5
 4.585465496869457e-5
 4.539610841900762e-5
 4.494214733481754e-5
 4.449272586146937e-5
 4.404779860285467e-5
 4.360732061682612e-5
 4.317124741065786e-5
In [3]:
##
lines(1:1000, norms_jordan, label="Jordan block", axis=(;yscale=log10))
#lines!(f.axis, norms_diag, label="Diagonal matrix")
#f
UndefVarError: `lines` not defined

Stacktrace:
 [1] top-level scope
   @ In[3]:2
In [4]:
##
using Plots
Plots.plot(1:1000, norms_jordan, yscale=:log10)
Plots.plot!(norms_diag)
┌ Warning: No strict ticks found
└ @ PlotUtils ~/.julia/packages/PlotUtils/jEGKP/src/ticks.jl:191
┌ Warning: No strict ticks found
└ @ PlotUtils ~/.julia/packages/PlotUtils/jEGKP/src/ticks.jl:191
┌ Warning: No strict ticks found
└ @ PlotUtils ~/.julia/packages/PlotUtils/jEGKP/src/ticks.jl:191
┌ Warning: No strict ticks found
└ @ PlotUtils ~/.julia/packages/PlotUtils/jEGKP/src/ticks.jl:191
Out[4]:
┌ Warning: No strict ticks found
└ @ PlotUtils ~/.julia/packages/PlotUtils/jEGKP/src/ticks.jl:191
┌ Warning: No strict ticks found
└ @ PlotUtils ~/.julia/packages/PlotUtils/jEGKP/src/ticks.jl:191
┌ Warning: No strict ticks found
└ @ PlotUtils ~/.julia/packages/PlotUtils/jEGKP/src/ticks.jl:191