2

I have a program to solve 1D Burguer's equation. I have the following values/dimensions, which are giving me a dimensional error. I don't understand what is going wrong. If anyone can help me, I appreciate it.

julia> size(xs) 
(50,)  

julia> size(ts) 
(1000,)  

julia> size(u) 
(50, 1001)  

julia> plot(xs,ts,u[1:50,1:1000],st=:surface, title="Burguer equation", xlabel="X", ylabel="Y", zlabel="U") 
Arrays have incorrect length or dimension.

More on the program itself, in case it helps

Parameters (space)

nx= 50;
delta_x = 15/(nx - 1)
x = range(0, stop=delta_x*(nx-1), length=nx) # Full range of spatial steps for which a solution is desired

Parameters (time)

endTime = 2   # simulation end time
nt = 1000          # nt is the number of timesteps we want to calculate
delta_t = endTime/nt  # Δt is the amount of time each timestep covers
t = range(0, stop=endTime, length=nt) # Full range of time steps for which a solution is desired

#+RESULTS: : 0.0:0.002002002002002002:2.0

Initial conditions (space-time)

# Init array of ones at initial timestep
u_zero = ones(nx) 
  
# Set u₀ = 2 in the interval 0.5 ≤ x ≤ 1 as per our I.C.s
u_zero[0.5 .<= x .<= 3] .= 2  # Note use of . (dot) broadcasting syntax
  
u_zero

#+RESULTS: : [1.0, 1.0, 2.0, 2.0, 2.0, 2.0, 2.0, 2.0, 2.0, 2.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0]

Run the differential equation

# u[:,] = copy(u_zero) # Initialize arbitrary future timestep with initial condition, u_zero
u=zeros((nx,nt+1))
u[:,1]=copy(u_zero)
  
for n in 1:nt       # loop over timesteps, n: nt times
    u[:,n+1] = copy(u[:,n]) # copy the existing values of u^n into u^(n+1)
    for i in 2:nx   
        u[i,n+1] = u[i,n] - u[i,n] * delta_t/delta_x * (u[i,n] - u[i-1,n])
    end
end
using Plots
gr()
plot(xs,ts,u[1:50,1:1000],st=:surface, title="Burguer equation", xlabel="X", ylabel="Y", zlabel="U")

COMMENT:

The solution consisted in changing the back-end for pyplot and to transpose the u matrix, with u'.

pyplot();

Plots.plot(xs,ts,u'[1:1000,:], st=:surface,title="Burguer equation", xlabel="X", ylabel="Time", zlabel="U")

RESULT: enter image description here

(@v1.7) pkg> st Plots
st -      Status `~/.julia/environments/v1.7/Project.toml`
  [91a5bcdd] Plots v1.25.6

(@v1.7) pkg> st -m GR
      Status `~/.julia/environments/v1.7/Manifest.toml`
  [28b8d3ca] GR v0.63.1

As we discussed in https://discourse.julialang.org/t/unmatched-dimensions-for-3d-plot/75180/13 part of the problem seem to be with the GR back-end, currently. PyPlot.jl back-end is working fine.

BuddhiLW
  • 608
  • 3
  • 9
  • Try transposing `u` with `u'`. That gave me a plot, but not sure if that's what you're looking for. – niczky12 Jan 25 '22 at 16:04
  • I still have `julia> plot(xs,ts,u',st=:surface, title="Burguer equation", xlabel="X", ylabel="Y", zlabel="U")` giving `Arrays have incorrect length or dimension.` – BuddhiLW Jan 25 '22 at 16:20
  • Can you post you result, @niczaky12 ? Thank you. – BuddhiLW Jan 25 '22 at 16:40

1 Answers1

2

As the comment suggest, the dimensions of your u matrix which holds the z coordinates are likely flipped. A simple example to illustrate:

julia> x = 1:5
1:5

julia> y = 1:10
1:10

julia> z = [x+y for x in 1:5, y in 1:10]
5×10 Matrix{Int64}:
 2  3  4  5   6   7   8   9  10  11
 3  4  5  6   7   8   9  10  11  12
 4  5  6  7   8   9  10  11  12  13
 5  6  7  8   9  10  11  12  13  14
 6  7  8  9  10  11  12  13  14  15

julia> surface(x, y, z)
Arrays have incorrect length or dimension.

julia> surface(x, y, z')

and the last command gives the correct plot (in this simple example it's easy to verify that z should indeed be 15 for x = 5 and y = 10):

enter image description here

Alternatively, you can pass three vectors to surface which is described in the docs here, but requires building grids out of your x and y vectors first:

f(x, a) = begin
        1 / x + a * x ^ 2
    end
xs = collect(0.1:0.05:2.0)
as = collect(0.2:0.1:2.0)
x_grid = [x for x = xs for y = as]
a_grid = [y for x = xs for y = as]
surface(x_grid, a_grid, f.(x_grid, a_grid))
Nils Gudat
  • 13,222
  • 3
  • 39
  • 60
  • Yes, as we discussed in https://discourse.julialang.org/t/unmatched-dimensions-for-3d-plot/75180/13 also the problem seem to be with the GR back-end currently. PyPlot.jl back-end is working fine. – BuddhiLW Jan 25 '22 at 17:38