Help solving coupled differential equations involving complex numbers

Specific Domains Modelling & Simulations
1

Original title: “Some help solving set of differential equations”

Hi. I am stuck in doing some homework, and I was hoping I could get some help.
The task is to analytically or numerically solve some differential equations. I will provide the full task for completeness:

image
image1796×594 174 KB

These are the DE’s to solve, in the following task:

image
image1960×334 77.8 KB

“Real valued” Rabi frequencies means that \Omega is simply some real number, which reduces the complexity quite a bit. Below is my attempt at defining this, using ModelingToolkit.jl:

using ModelingToolkit
using DifferentialEquations

@variables t c̃₁(t) c̃₂(t)    # independent and dependent variables
@parameters Ω Δ             # parameters 

∂ₜ = Differential(t) # define an operator for the differentiation w.r.t. time

Inspired by the docs, I tried to define a vector of equations:

@named eqs = ODESystem([∂ₜ(c̃₁) ~ Ω*im*c̃₂*cis(t*Δ),
                        ∂ₜ(c̃₂) ~ Ω*im*c̃₁*cis(-t*Δ)
])

Resulting in

ERROR: type Array has no field lhs
Stacktrace:

So instead, I tried what I found inntuitive - defining two equations, and giving those as a vector:

julia> @named eq1 = ODESystem(∂ₜ(c̃₁) ~ Ω*im*c̃₂*cis(t*Δ))
Model eq1 with 2 equations
States (2):
  c̃₁(t)
  c̃₂(t)
Parameters (2):
  Ω
  Δ

julia> @named eq2 = ODESystem(∂ₜ(c̃₂) ~ Ω*im*c̃₁*cis(-t*Δ))
Model eq2 with 2 equations
States (2):
  c̃₂(t)
  c̃₁(t)
Parameters (2):
  Ω
  Δ

So far so good I think. I do not get why both equations say “with 2 equations”.
I then attempt:

julia> using DifferentialEquations: solve

julia> initial_conds = [c̃₁ => 1.0, c̃₂ => 0.0];

julia> tspan = (0.0, 10.0);

julia> params = [Ω => 1.0, Δ => 1.0];

julia> prob = ODEProblem([eq1, eq2], initial_conds, tspan, params)
ERROR: No methods were found for the model function passed to the equation solver.
The function `f` needs to have dispatches, for example, for an ODEProblem
`f` must define either `f(u,p,t)` or `f(du,u,p,t)`. For more information
on how the model function `f` should be defined, consult the docstring for
the appropriate `AbstractSciMLFunction`.

Offending function: f
Stacktrace:
 [1] isinplace(f::Vector{ODESystem}, inplace_param_number::Int64, fname::String, iip_preferred::Bool; has_two_dispatches::Bool, isoptimization::Bool)
   @ SciMLBase ~/.julia/packages/SciMLBase/RAGXU/src/utils.jl:236
 [2] isinplace (repeats 2 times)
   @ ~/.julia/packages/SciMLBase/RAGXU/src/utils.jl:229 [inlined]
 [3] ODEProblem(f::Vector{ODESystem}, u0::Vector{Pair{Num, Float64}}, tspan::Tuple{Float64, Float64}, p::Vector{Pair{Num, Float64}}; kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
   @ SciMLBase ~/.julia/packages/SciMLBase/RAGXU/src/problems/ode_problems.jl:168
 [4] ODEProblem(f::Vector{ODESystem}, u0::Vector{Pair{Num, Float64}}, tspan::Tuple{Float64, Float64}, p::Vector{Pair{Num, Float64}})
   @ SciMLBase ~/.julia/packages/SciMLBase/RAGXU/src/problems/ode_problems.jl:167
 [5] top-level scope
   @ REPL[33]:1

Now I realize that I could be going about it the wrong way. Is it even an ODE system, or something else? I do however not know how to proceed, and would much appreciate some pointers.

Do anyone have the know-how and time to help me out?

2

ModelingToolkit variables are all real, as written you get

julia> eqs = [∂ₜ(c̃₁) ~ Ω*im*c̃₂*cis(t*Δ),
              ∂ₜ(c̃₂) ~ Ω*im*c̃₁*cis(-t*Δ)
       ]
2-element Vector{Vector{Equation}}:
 [Differential(t)(c̃₁(t)) ~ -Ω*c̃₂(t)*sin(t*Δ), 0 ~ Ω*c̃₂(t)*cos(t*Δ)]
 [Differential(t)(c̃₂(t)) ~ -Ω*c̃₁(t)*sin(-t*Δ), 0 ~ Ω*c̃₁(t)*cos(-t*Δ)]

Which probably isn’t what you mean. Try defining c̃₁ and c̃₂ with explicit real and imaginary parts.

@variables c̃r₁ c̃i₁ c̃₁ c̃r₂ c̃i₂ c̃₂
c̃₁ = c̃r₁ + im*c̃i₁
c̃₂ = c̃r₂ + im*c̃i₂
1 Like
3

Thank you for your suggestion.
With the following definitions:

@variables t c̃r₁ c̃i₁ c̃₁ c̃r₂ c̃i₂ c̃₂  # independent and dependent variables
@parameters Ω Δ                     # parameters 

c̃₁ = c̃r₁ + im*c̃i₁
c̃₂ = c̃r₂ + im*c̃i₂

∂ₜ = Differential(t) # define an operator for the differentiation w.r.t. time

I get the following error:

julia> @named eqs = ODESystem([∂ₜ(c̃₁) ~ Ω*im*c̃₂*cis(t*Δ),
                               ∂ₜ(c̃₂) ~ Ω*im*c̃₁*cis(-t*Δ)
       ])
ERROR: type ComplexTerm has no field metadata
Stacktrace:
  [1] getproperty(x::Symbolics.ComplexTerm{Real}, f::Symbol)
 ...

I get the same error when defining a single equation:

julia> @named eq1 = ODESystem(∂ₜ(c̃₁) ~ Ω*im*c̃₂*cis(t*Δ))
ERROR: type ComplexTerm has no field metadata

Should I define seperate equations for the real and imaginary parts, giving me 4 equations? If yes, I am unsure how to even go about that. Should I use the real and imag functions in defining the equations?

4

Open an issue. This is a very specific question to ModelingToolkit with complex numbers and it would be helpful if you title it as such.

5

Done. Are you saying this could be solved without ModelingToolkit? Any quick pointers on how?

6

Just use complex numbers in DiffEq directly.

2 Likes
7

Thanks guys. This ended up doing the trick:

using DifferentialEquations
using ParameterizedFunctions

my_ode = @ode_def begin
    dc̃₁ = Ω*im*c̃₂*exp( im*t*Δ)
    dc̃₂ = Ω*im*c̃₁*exp(-im*t*Δ)
end Ω Δ

u₀ = Complex{Float64}[1, 0]
tspan = (0.0, 10.0)
p = [1.0, 1.0]

prob = ODEProblem(my_ode, u₀, tspan, p)
sol = solve(prob)
3 Likes