Complex number error

[TrumeAAA]

Hi,

I tried to solve a second order differential equation with complex number in the equation and failed. However, when I set the imaginary part to be zero, I get the correct result exactly the same as that I get using Matlab.
The equation is in the following picture:

The "-" in the second term on the right of the equation should be "+".

And here is my Julia code:
using DifferentialEquations
using Plots

mₑ = 1
mᵢ = 1836
c = 3.0e8
λ₀ = 1.0
τ₀ = λ₀*(1.0e-6)/c
T₀ = 0.6e-12 # pulse duration in ps
ω₀ = 2pi/τ₀
I₀ = 7e16
I₁ = 1.37e18/(λ₀^2)
a₀₀ = sqrt(I₀/I₁)
E₀ = 2.75e3sqrt(I₀)
nc = 1.1e27λ₀^(-2)
nₑ = 0.45nc
Tᵢ = 500000 # units in eV
Tₑ = Tᵢ
λ_De = sqrt( Tᵢ*(8.854e-12)/(nₑ*(1.602e-19)) )
N_D = 4pinₑ*(λ_De^3)/3
Λ = 9N_D
νₑᵢ = (2.91e-6)nₑlog(Λ)(Tₑ^(-1.5))
ωₚ₀ = (5.641e4)sqrt(nₑ10^(-6))
A₁ = (9.1e-31)(3e8)^2/((1.602e-19)Tᵢ) # mc²/Tᵢ
A₂ = mₑ/mᵢ # mₑ/mᵢ
A₃ = ωₚ₀^2/(ω₀(ω₀)) # ωₚ₀²/(ω₀(ω₀+νₑᵢim))
A₄ = a₀₀^2
τ₁ = T₀ω₀

function f(du,u,p,t)

# calculate ϵ₀
χ₁ = A₃/sqrt(1+A₄/u)
χ₂ = 1 + (A₁*A₄)/( 12*((1+A₄/u)^2)*A₂*u )
χ₃ = -1-6*sqrt(1+A₄/u)*A₂
ϵ₀ = 1-χ₁*(χ₂^(χ₃))

# calculate ϵ₂
B₂ = 1 + A₄/u
B₁ = -1-6*A₂*sqrt(B₂)
C₁ = 12*u*A₂*((B₂)^2)
C₂ = 2*(u^3)*((B₂)^1.5)
C₃ = (u^3)*sqrt(B₂)
C₄ = 6*u*A₂*((B₂)^3)
C₅ = 12*A₂*((B₂)^2)
𝜂₁ = (A₃*A₄/C₂)*((1+A₁*A₄/C₁)^B₁)
𝜂₂ = (A₃*A₄/C₃)*((1+A₁*A₄/C₁)^B₁)
𝜂₃ = (-A₁*A₄/C₄-A₁/C₅)*B₁
𝜂₄ = (3*A₂/sqrt(B₂))*log(1+A₁*A₄/C₁)
ϵ₂ = 𝜂₁-𝜂₂*(𝜂₃-𝜂₄)
3.1415926

ddu = 1/(ϵ₀*((τ₁)^4)*(u^3)) + ϵ₂*u/(ϵ₀*((τ₁)^2))

end

u0 = 1.0
du0 = 0.0
tspan = (0.0,5000.0)
prob = SecondOrderODEProblem(f,du0,u0,tspan,callback=CallbackSet())
sol = solve(prob,RK4(),reltol=1e-8,abstol=1e-8,saveat=1) # Tsit5()
d = sol[2,:]

theme(:dao)

plot(d,
linewidth=5,
label="g",
xaxis="ξ",
yaxis="g",
grid="off",
title="g evolution with dimensionless distance ξ",
# xtickfont = font(15),
# ytickfont = font(15),
border=:box
# xscale = :log,
# size = (1200,800),
)

ω₀ should be complex: ω₀+νₑᵢ*im. However, when I add the imaginary part, there's an error, and I can't figure out how to fix it.
The error is:
ERROR: InexactError: Float64(-3.408153847261804e-7 + 1.5016101254502642e-8im)
Stacktrace:
[1] Real
@ .\complex.jl:44 [inlined]
[2] convert
@ .\number.jl:7 [inlined]
[3] #96
@ C:\install\JuliaPkgs\packages\RecursiveArrayTools\IgoeN\src\array_partition.jl:503 [inlined]
[4] ntuple
@ .\ntuple.jl:49 [inlined]
[5] convert(#unused#::Type{ArrayPartition{Float64, Tuple{Float64, Float64}}}, A::ArrayPartition{ComplexF64, Tuple{ComplexF64, Float64}})
@ RecursiveArrayTools C:\install\JuliaPkgs\packages\RecursiveArrayTools\IgoeN\src\array_partition.jl:502
[6] setproperty!(x::OrdinaryDiffEq.ODEIntegrator{RK4, false, ArrayPartition{Float64, Tuple{Float64, Float64}}, Nothing, Float64, SciMLBase.NullParameters, Float64, Float64, Float64, Float64, Vector{ArrayPartition{Float64, Tuple{Float64, Float64}}}, ODESolution{Float64, 2, Vector{ArrayPartition{Float64, Tuple{Float64, Float64}}}, Nothing, Nothing, Vector{Float64}, Vector{Vector{ArrayPartition{Float64, Tuple{Float64, Float64}}}}, ODEProblem{ArrayPartition{Float64, Tuple{Float64, Float64}}, Tuple{Float64, Float64}, false, SciMLBase.NullParameters, DynamicalODEFunction{false, SciMLBase.FullSpecialize, ODEFunction{false, SciMLBase.FullSpecialize, typeof(f), LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, ODEFunction{false, SciMLBase.FullSpecialize, SciMLBase.var"#235#237", LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, Base.Pairs{Symbol, CallbackSet{Tuple{}, Tuple{}}, Tuple{Symbol}, NamedTuple{(:callback,), Tuple{CallbackSet{Tuple{}, Tuple{}}}}}, SecondOrderODEProblem{false}}, RK4, OrdinaryDiffEq.InterpolationData{DynamicalODEFunction{false, SciMLBase.FullSpecialize, ODEFunction{false, SciMLBase.FullSpecialize, typeof(f), LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, ODEFunction{false, SciMLBase.FullSpecialize, SciMLBase.var"#235#237", LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, Vector{ArrayPartition{Float64, Tuple{Float64, Float64}}}, Vector{Float64}, Vector{Vector{ArrayPartition{Float64, Tuple{Float64, Float64}}}}, OrdinaryDiffEq.RK4ConstantCache}, DiffEqBase.DEStats}, DynamicalODEFunction{false, SciMLBase.FullSpecialize, ODEFunction{false, SciMLBase.FullSpecialize, typeof(f), LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, ODEFunction{false, SciMLBase.FullSpecialize, SciMLBase.var"#235#237", LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, OrdinaryDiffEq.RK4ConstantCache, OrdinaryDiffEq.DEOptions{Float64, Float64, Float64, Float64, PIController{Rational{Int64}}, typeof(DiffEqBase.ODE_DEFAULT_NORM), typeof(LinearAlgebra.opnorm), Nothing, CallbackSet{Tuple{}, Tuple{}}, typeof(DiffEqBase.ODE_DEFAULT_ISOUTOFDOMAIN), typeof(DiffEqBase.ODE_DEFAULT_PROG_MESSAGE), typeof(DiffEqBase.ODE_DEFAULT_UNSTABLE_CHECK), DataStructures.BinaryHeap{Float64, DataStructures.FasterForward}, DataStructures.BinaryHeap{Float64, DataStructures.FasterForward}, Nothing, Nothing, Int64, Tuple{}, Int64, Tuple{}}, ArrayPartition{Float64, Tuple{Float64, Float64}}, Float64, Nothing, OrdinaryDiffEq.DefaultInit}, f::Symbol, v::ArrayPartition{ComplexF64, Tuple{ComplexF64, Float64}})
@ Base .\Base.jl:43
[7] initialize!(integrator::OrdinaryDiffEq.ODEIntegrator{RK4, false, ArrayPartition{Float64, Tuple{Float64, Float64}}, Nothing, Float64, SciMLBase.NullParameters, Float64, Float64, Float64, Float64, Vector{ArrayPartition{Float64, Tuple{Float64, Float64}}}, ODESolution{Float64, 2, Vector{ArrayPartition{Float64, Tuple{Float64, Float64}}}, Nothing, Nothing, Vector{Float64}, Vector{Vector{ArrayPartition{Float64, Tuple{Float64, Float64}}}}, ODEProblem{ArrayPartition{Float64, Tuple{Float64, Float64}}, Tuple{Float64, Float64}, false, SciMLBase.NullParameters, DynamicalODEFunction{false, SciMLBase.FullSpecialize, ODEFunction{false, SciMLBase.FullSpecialize, typeof(f), LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, ODEFunction{false, SciMLBase.FullSpecialize, SciMLBase.var"#235#237", LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, Base.Pairs{Symbol, CallbackSet{Tuple{}, Tuple{}}, Tuple{Symbol}, NamedTuple{(:callback,), Tuple{CallbackSet{Tuple{}, Tuple{}}}}}, SecondOrderODEProblem{false}}, RK4, OrdinaryDiffEq.InterpolationData{DynamicalODEFunction{false, SciMLBase.FullSpecialize, ODEFunction{false, SciMLBase.FullSpecialize, typeof(f), LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, ODEFunction{false, SciMLBase.FullSpecialize, SciMLBase.var"#235#237", LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, Vector{ArrayPartition{Float64, Tuple{Float64, Float64}}}, Vector{Float64}, Vector{Vector{ArrayPartition{Float64, Tuple{Float64, Float64}}}}, OrdinaryDiffEq.RK4ConstantCache}, DiffEqBase.DEStats}, DynamicalODEFunction{false, SciMLBase.FullSpecialize, ODEFunction{false, SciMLBase.FullSpecialize, typeof(f), LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, ODEFunction{false, SciMLBase.FullSpecialize, SciMLBase.var"#235#237", LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, OrdinaryDiffEq.RK4ConstantCache, OrdinaryDiffEq.DEOptions{Float64, Float64, Float64, Float64, PIController{Rational{Int64}}, typeof(DiffEqBase.ODE_DEFAULT_NORM), typeof(LinearAlgebra.opnorm), Nothing, CallbackSet{Tuple{}, Tuple{}}, typeof(DiffEqBase.ODE_DEFAULT_ISOUTOFDOMAIN), typeof(DiffEqBase.ODE_DEFAULT_PROG_MESSAGE), typeof(DiffEqBase.ODE_DEFAULT_UNSTABLE_CHECK), DataStructures.BinaryHeap{Float64, DataStructures.FasterForward}, DataStructures.BinaryHeap{Float64, DataStructures.FasterForward}, Nothing, Nothing, Int64, Tuple{}, Int64, Tuple{}}, ArrayPartition{Float64, Tuple{Float64, Float64}}, Float64, Nothing, OrdinaryDiffEq.DefaultInit}, cache::OrdinaryDiffEq.RK4ConstantCache)
@ OrdinaryDiffEq C:\install\JuliaPkgs\packages\OrdinaryDiffEq\4qO6L\src\perform_step\fixed_timestep_perform_step.jl:259
[8] __init(prob::ODEProblem{ArrayPartition{Float64, Tuple{Float64, Float64}}, Tuple{Float64, Float64}, false, SciMLBase.NullParameters, DynamicalODEFunction{false, SciMLBase.FullSpecialize, ODEFunction{false, SciMLBase.FullSpecialize, typeof(f), LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, ODEFunction{false, SciMLBase.FullSpecialize, SciMLBase.var"#235#237", LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, Base.Pairs{Symbol, CallbackSet{Tuple{}, Tuple{}}, Tuple{Symbol}, NamedTuple{(:callback,), Tuple{CallbackSet{Tuple{}, Tuple{}}}}}, SecondOrderODEProblem{false}}, alg::RK4, timeseries_initmbol, Any, NTuple{4, Symbol}, NamedTuple{(:callback, :reltol, :abstol, :saveat), Tuple{CallbackSet{Tuple{}, Tuple{}}, Float64, Float64, Int64}}}) @ OrdinaryDiffEq C:\install\JuliaPkgs\packages\OrdinaryDiffEq\4qO6L\src\solve.jl:5 [10] #solve_call#26 @ C:\install\JuliaPkgs\packages\DiffEqBase\Lq1gG\src\solve.jl:473 [inlined] [11] solve_up(prob::ODEProblem{ArrayPartition{Float64, Tuple{Float64, Float64}}, Tuple{Float64, Float64}, false, SciMLBase.NullParameters, DynamicalODEFunction{false, SciMLBase.FullSpecialize, ODEFunction{false, SciMLBase.FullSpecialize, typeof(f), LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, ODEFunction{false, SciMLBase.FullSpecialize, SciMLBase.var"#235#237", LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, Base.Pairs{Symbol, CallbackSet{Tuple{}, Tuple{}}, Tuple{Symbol}, NamedTuple{(:callback,), Tuple{CallbackSet{Tuple{}, Tuple{}}}}}, SecondOrderODEProblem{false}}, sensealg::Nothing, u0::ArrayPartition{Float64, Tuple{Float64, Float64}}, p::SciMLBase.NullParameters, args::RK4; kwargs::Base.Pairs{Symbol, Real, Tuple{Symbol, Symbol, Symbol}, NamedTuple{(:reltol, :abstol, :saveat), Tuple{Float64, Float64, Int64}}}) @ DiffEqBase C:\install\JuliaPkgs\packages\DiffEqBase\Lq1gG\src\solve.jl:835 [12] #solve#31 @ C:\install\JuliaPkgs\packages\DiffEqBase\Lq1gG\src\solve.jl:802 [inlined] [13] top-level scope
@ d:\Julia\pulsefc.jl:62

I wonder that if I used the wrong type of the numbers. What do I need to do if I want to solve equations with complex?

Any help or suggestion will be helpful.

Thanks.

[ChrisRackauckas]

u0 = 1.0
du0 = 0.0
tspan = (0.0,5000.0)
prob = SecondOrderODEProblem(f,du0,u0,tspan,callback=CallbackSet())

your initial conditions are not complex numbers, so it will not compute in complex numbers. Please change these to complex if you want them to compute the ODE in complex space, i.e.

u0 = 1.0 + 0im
du0 = 0.0 + 0im

[TrumeAAA]

u0 = 1.0
du0 = 0.0
tspan = (0.0,5000.0)
prob = SecondOrderODEProblem(f,du0,u0,tspan,callback=CallbackSet())

your initial conditions are not complex numbers, so it will not compute in complex numbers. Please change these to complex if you want them to compute the ODE in complex space, i.e.

u0 = 1.0 + 0im
du0 = 0.0 + 0im

It works! Thanks a lot!