# Study SH model

**URL:** <https://discourse.julialang.org/t/study-sh-model/61246>\
**Category:** New to Julia\
**Tags:** diffeq, ode\
**Created:** [May 16, 2021, 2:05pm UTC](https://discourse.julialang.org/t/study-sh-model/61246 "2021-05-16T14:05:07Z")\
**Posts on this page:** 7\
**Page:** 1

<div class="post-metadata">

**Author:** ![Frazze](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/frazze/32/25142_2.png) [@Frazze](https://discourse.julialang.org/u/Frazze)\
**Post date:** [May 16, 2021, 2:05pm UTC](https://discourse.julialang.org/t/study-sh-model/61246/1 "2021-05-16T14:05:07Z")

</div>

Hi everyone,

I’m studying the Swift Hohenberg model with only the linear term of the potential. I see that a tutorial is proposed here ([https://rveltz.github.io/BifurcationKit.jl/dev/tutorials2/#Snaking-in-the-2d-Swift-Hohenberg-equation-1](https://rveltz.github.io/BifurcationKit.jl/dev/tutorials2/#Snaking-in-the-2d-Swift-Hohenberg-equation-1)) but for the moment I’m not interested to analyze the bifurcation, I wish only to solve the model. My idea was to use DifferentialEquations.jl to solve the equation, but maybe I’ve made some mistakes…  
Here my definition of the problem:

```julia
function F_sh(u, p)
	@unpack l, L1 = p
	return -L1 * u .+ (l .* u)
end
X = -lx .+ 2lx/(Nx) * collect(0:Nx-1)
Y = -ly .+ 2ly/(Ny) * collect(0:Ny-1)

# define parameters for the PDE
Δ, _ = Laplacian2D(Nx, Ny, lx, ly)
D = 10
#L1 = D*(I + Δ)^2

p = (l = -0.1 , L1 = D*(I + Δ)^2)

u0=rand(100,100)
tspan = (0.0,100.0)
prob = ODEProblem(F_sh,u0,tspan)

sol = solve(prob, saveat =[0.2,25.0,50.0,99.0])

```

(I used the Laplacian construction with sparse matrix used in the SH tutorial).  
Thank you for the help, I’m very new to Julia, so every comment is welcome!

---

<div class="post-metadata">

**Author:** ![rveltz](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/rveltz/32/2707_2.png) [@rveltz](https://discourse.julialang.org/u/rveltz)\
**Post date:** [May 16, 2021, 3:49pm UTC](https://discourse.julialang.org/t/study-sh-model/61246/2 "2021-05-16T15:49:10Z")

</div>

without running it, I’d say you have to pass a flat vector because of the sparse matrix:

```julia
prob = ODEProblem(F_sh,vec(u0),tspan)

```

In passing, you’d better be using a `SplitProblem` I think

---

<div class="post-metadata">

**Author:** ![Frazze](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/frazze/32/25142_2.png) [@Frazze](https://discourse.julialang.org/u/Frazze)\
**Post date:** [May 16, 2021, 8:30pm UTC](https://discourse.julialang.org/t/study-sh-model/61246/3 "2021-05-16T20:30:56Z")

</div>

Thank you so much,  
I tried using the flat vector vec(u0), but it return the same error as before:

> MethodError: no method matching F\_sh(::Matrix{Float64}, ::SciMLBase.NullParameters, ::Float64)  
> Closest candidates are:  
> F\_sh(::Any, ::Any) at In[4]:1

I tried also with SplitODEProblem, splicing F\_sh into two functions sh1 and sh2, but still the error is returned:

> MethodError: no method matching sh1(::Matrix{Float64}, ::SciMLBase.NullParameters, ::Float64)

---

<div class="post-metadata">

**Author:** ![rveltz](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/rveltz/32/2707_2.png) [@rveltz](https://discourse.julialang.org/u/rveltz)\
**Post date:** [May 16, 2021, 8:39pm UTC](https://discourse.julialang.org/t/study-sh-model/61246/4 "2021-05-16T20:39:43Z")

</div>

Ah, try:

```julia
function F_sh(u, p, t = 0)
	@unpack l, L1 = p
	return -L1 * u .+ (l .* u)
end

```

---

<div class="post-metadata">

**Author:** ![Frazze](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/frazze/32/25142_2.png) [@Frazze](https://discourse.julialang.org/u/Frazze)\
**Post date:** [May 17, 2021, 6:50am UTC](https://discourse.julialang.org/t/study-sh-model/61246/5 "2021-05-17T06:50:57Z")

</div>

I’ve already tried this, but the error still remain.  
I put here the complete notebook, maybe it could be useful.  
Best regards

Francesco

---

<div class="post-metadata">

**Author:** ![rveltz](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/rveltz/32/2707_2.png) [@rveltz](https://discourse.julialang.org/u/rveltz)\
**Post date:** [May 17, 2021, 8:35am UTC](https://discourse.julialang.org/t/study-sh-model/61246/6 "2021-05-17T08:35:00Z")

</div>

At least it works but slowly because of bad solver options

```julia

function F_sh(u, p, t)
	@unpack l, L1 = p
	return -L1 * u .+ (l .* u)
end
X = -lx .+ 2lx/(Nx) * collect(0:Nx-1)
Y = -ly .+ 2ly/(Ny) * collect(0:Ny-1)

# define parameters for the PDE
Δ, _ = Laplacian2D(Nx, Ny, lx, ly)
D = 10
#L1 = D*(I + Δ)^2

p = (l = -0.1 , L1 = D*(I + Δ)^2)

u0=rand(Nx,Ny) |> vec
tspan = (0.0,1.0)
prob = ODEProblem(F_sh,u0,tspan, p)

sol = solve(prob, Tsit5(), saveat =[0.2,25.0,50.0,99.0])

```

---

<div class="post-metadata">

**Author:** ![Frazze](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/frazze/32/25142_2.png) [@Frazze](https://discourse.julialang.org/u/Frazze)\
**Post date:** [May 17, 2021, 9:18am UTC](https://discourse.julialang.org/t/study-sh-model/61246/7 "2021-05-17T09:18:41Z")

</div>

I did the same, but still the error is returned… I really don’t understand what I’m missing.

```julia
using DifferentialEquations
using DiffEqOperators, Setfield, Parameters
using BifurcationKit, LinearAlgebra, Plots, SparseArrays
const BK = BifurcationKit

# helper function to plot solution
heatmapsol(x) = heatmap(reshape(x,Nx,Ny)',color=:viridis)

Nx = 151
Ny = 100
lx = 4*2pi
ly = 2*2pi/sqrt(3)

# we use DiffEqOperators to compute the Laplacian operator
function Laplacian2D(Nx, Ny, lx, ly)
	hx = 2lx/Nx
	hy = 2ly/Ny
	D2x = CenteredDifference(2, 2, hx, Nx)
	D2y = CenteredDifference(2, 2, hy, Ny)
	Qx = Neumann0BC(hx)
	Qy = Neumann0BC(hy)
	
	A = kron(sparse(I, Ny, Ny), sparse(D2x * Qx)[1]) + kron(sparse(D2y * Qy)[1], sparse(I, Nx, Nx))
	return A, D2x
end

function F_sh(u, p, t)
	@unpack l, L1 = p
	return -L1 .* u .+ (l .* u)
end

X = -lx .+ 2lx/(Nx) * collect(0:Nx-1)
Y = -ly .+ 2ly/(Ny) * collect(0:Ny-1)

# define parameters for the PDE
Δ, _ = Laplacian2D(Nx, Ny, lx, ly)
D = 10
#L1 = D*(I + Δ)^2

p = (l= -0.1, L1 = D*(I + Δ)^2);

u0=rand(Nx,Ny) |> vec
tspan = (0.0,100.0)
prob = ODEProblem(F_sh,u0,tspan,p)
sol = solve(prob, Tsit5(), saveat =[0.2,50.0,99.0])

```

> MethodError: no method matching Vector{Float64}(::SparseMatrixCSC{Float64, Int64})
