# Generating equations in ModellingToolkit from infinite series

Hi all!

I was advised to open this issue as a follow-up of this thread: https://discourse.julialang.org/t/speeding-up-solution-to-the-large-system-of-non-linear-complex-valued-odes/39075

My issue: I am generating nonlinear complex-valued ODES from the infinite series, by essentially summing out the series for a given N (number of terms). Series have linear, quadratic and cubic terms in it (single, double and triple sums), making equations (summations) very long for large N as the number of terms in the entire system grows approximately as a quadratic function of N. Therefore when I call function containing my equations in the solver, code allocates large amounts of memory, making it impossible for me to solve for realistic N (I would be looking at N=100 as a minimum).
Examples of generated equations, which I give to the solver:

``````using DifferentialEquations, JLD2, StaticArrays

const tspan=(0.0,100.0)
const t_name=tspan[2]
function Nodes(u,p,t) #equations to be solved. For large N, they get very long and ugly
A1=((9*((p[3])^2/2)-p[2]/2*abs2(p[1]))*u[1]-p[2]/2*(p[1])^2*conj(u[7]) - (p[2]/2)*(2*p[1]*u[1]*conj(u[4])+conj(p[1])*u[1]*u[4]+2*p[1]*u[2]*conj(u[5])+conj(p[1])*u[2]*u[3]+2*p[1]*u[3]*conj(u[6])+conj(p[1])*u[3]*u[2]+2*p[1]*u[4]*conj(u[7])+conj(p[1])*u[4]*u[1] + u[1]*conj(u[1])*u[1]+u[1]*conj(u[2])*u[2]+u[2]*conj(u[2])*u[1]+u[1]*conj(u[3])*u[3]+u[2]*conj(u[3])*u[2]+u[3]*conj(u[3])*u[1]+u[1]*conj(u[4])*u[4]+u[2]*conj(u[4])*u[3]+u[3]*conj(u[4])*u[2]+u[4]*conj(u[4])*u[1]+u[1]*conj(u[5])*u[5]+u[2]*conj(u[5])*u[4]+u[3]*conj(u[5])*u[3]+u[4]*conj(u[5])*u[2]+u[5]*conj(u[5])*u[1]+u[1]*conj(u[6])*u[6]+u[2]*conj(u[6])*u[5]+u[3]*conj(u[6])*u[4]+u[4]*conj(u[6])*u[3]+u[5]*conj(u[6])*u[2]+u[6]*conj(u[6])*u[1]+u[1]*conj(u[7])*u[7]+u[2]*conj(u[7])*u[6]+u[3]*conj(u[7])*u[5]+u[4]*conj(u[7])*u[4]+u[5]*conj(u[7])*u[3]+u[6]*conj(u[7])*u[2]+u[7]*conj(u[7])*u[1]))/im
A2=((4*((p[3])^2/2)-p[2]/2*abs2(p[1]))*u[2]-p[2]/2*(p[1])^2*conj(u[6]) - (p[2]/2)*(2*p[1]*u[1]*conj(u[3])+conj(p[1])*u[1]*u[5]+2*p[1]*u[2]*conj(u[4])+conj(p[1])*u[2]*u[4]+2*p[1]*u[3]*conj(u[5])+conj(p[1])*u[3]*u[3]+2*p[1]*u[4]*conj(u[6])+conj(p[1])*u[4]*u[2]+2*p[1]*u[5]*conj(u[7])+conj(p[1])*u[5]*u[1] + u[1]*conj(u[1])*u[2]+u[2]*conj(u[1])*u[1]+u[1]*conj(u[2])*u[3]+u[2]*conj(u[2])*u[2]+u[3]*conj(u[2])*u[1]+u[1]*conj(u[3])*u[4]+u[2]*conj(u[3])*u[3]+u[3]*conj(u[3])*u[2]+u[4]*conj(u[3])*u[1]+u[1]*conj(u[4])*u[5]+u[2]*conj(u[4])*u[4]+u[3]*conj(u[4])*u[3]+u[4]*conj(u[4])*u[2]+u[5]*conj(u[4])*u[1]+u[1]*conj(u[5])*u[6]+u[2]*conj(u[5])*u[5]+u[3]*conj(u[5])*u[4]+u[4]*conj(u[5])*u[3]+u[5]*conj(u[5])*u[2]+u[6]*conj(u[5])*u[1]+u[1]*conj(u[6])*u[7]+u[2]*conj(u[6])*u[6]+u[3]*conj(u[6])*u[5]+u[4]*conj(u[6])*u[4]+u[5]*conj(u[6])*u[3]+u[6]*conj(u[6])*u[2]+u[7]*conj(u[6])*u[1]+u[2]*conj(u[7])*u[7]+u[3]*conj(u[7])*u[6]+u[4]*conj(u[7])*u[5]+u[5]*conj(u[7])*u[4]+u[6]*conj(u[7])*u[3]+u[7]*conj(u[7])*u[2]))/im
A3=((1*((p[3])^2/2)-p[2]/2*abs2(p[1]))*u[3]-p[2]/2*(p[1])^2*conj(u[5]) - (p[2]/2)*(2*p[1]*u[1]*conj(u[2])+conj(p[1])*u[1]*u[6]+2*p[1]*u[2]*conj(u[3])+conj(p[1])*u[2]*u[5]+2*p[1]*u[3]*conj(u[4])+conj(p[1])*u[3]*u[4]+2*p[1]*u[4]*conj(u[5])+conj(p[1])*u[4]*u[3]+2*p[1]*u[5]*conj(u[6])+conj(p[1])*u[5]*u[2]+2*p[1]*u[6]*conj(u[7])+conj(p[1])*u[6]*u[1] + u[1]*conj(u[1])*u[3]+u[2]*conj(u[1])*u[2]+u[3]*conj(u[1])*u[1]+u[1]*conj(u[2])*u[4]+u[2]*conj(u[2])*u[3]+u[3]*conj(u[2])*u[2]+u[4]*conj(u[2])*u[1]+u[1]*conj(u[3])*u[5]+u[2]*conj(u[3])*u[4]+u[3]*conj(u[3])*u[3]+u[4]*conj(u[3])*u[2]+u[5]*conj(u[3])*u[1]+u[1]*conj(u[4])*u[6]+u[2]*conj(u[4])*u[5]+u[3]*conj(u[4])*u[4]+u[4]*conj(u[4])*u[3]+u[5]*conj(u[4])*u[2]+u[6]*conj(u[4])*u[1]+u[1]*conj(u[5])*u[7]+u[2]*conj(u[5])*u[6]+u[3]*conj(u[5])*u[5]+u[4]*conj(u[5])*u[4]+u[5]*conj(u[5])*u[3]+u[6]*conj(u[5])*u[2]+u[7]*conj(u[5])*u[1]+u[2]*conj(u[6])*u[7]+u[3]*conj(u[6])*u[6]+u[4]*conj(u[6])*u[5]+u[5]*conj(u[6])*u[4]+u[6]*conj(u[6])*u[3]+u[7]*conj(u[6])*u[2]+u[3]*conj(u[7])*u[7]+u[4]*conj(u[7])*u[6]+u[5]*conj(u[7])*u[5]+u[6]*conj(u[7])*u[4]+u[7]*conj(u[7])*u[3]))/im
A4=((-p[2]/2*abs2(p[1]))*u[4]-p[2]/2*(p[1])^2*conj(u[4]) - (p[2]/2)*(2*p[1]*u[1]*conj(u[1])+conj(p[1])*u[1]*u[7]+2*p[1]*u[2]*conj(u[2])+conj(p[1])*u[2]*u[6]+2*p[1]*u[3]*conj(u[3])+conj(p[1])*u[3]*u[5]+2*p[1]*u[4]*conj(u[4])+conj(p[1])*u[4]*u[4]+2*p[1]*u[5]*conj(u[5])+conj(p[1])*u[5]*u[3]+2*p[1]*u[6]*conj(u[6])+conj(p[1])*u[6]*u[2]+2*p[1]*u[7]*conj(u[7])+conj(p[1])*u[7]*u[1] + u[1]*conj(u[1])*u[4]+u[2]*conj(u[1])*u[3]+u[3]*conj(u[1])*u[2]+u[4]*conj(u[1])*u[1]+u[1]*conj(u[2])*u[5]+u[2]*conj(u[2])*u[4]+u[3]*conj(u[2])*u[3]+u[4]*conj(u[2])*u[2]+u[5]*conj(u[2])*u[1]+u[1]*conj(u[3])*u[6]+u[2]*conj(u[3])*u[5]+u[3]*conj(u[3])*u[4]+u[4]*conj(u[3])*u[3]+u[5]*conj(u[3])*u[2]+u[6]*conj(u[3])*u[1]+u[1]*conj(u[4])*u[7]+u[2]*conj(u[4])*u[6]+u[3]*conj(u[4])*u[5]+u[4]*conj(u[4])*u[4]+u[5]*conj(u[4])*u[3]+u[6]*conj(u[4])*u[2]+u[7]*conj(u[4])*u[1]+u[2]*conj(u[5])*u[7]+u[3]*conj(u[5])*u[6]+u[4]*conj(u[5])*u[5]+u[5]*conj(u[5])*u[4]+u[6]*conj(u[5])*u[3]+u[7]*conj(u[5])*u[2]+u[3]*conj(u[6])*u[7]+u[4]*conj(u[6])*u[6]+u[5]*conj(u[6])*u[5]+u[6]*conj(u[6])*u[4]+u[7]*conj(u[6])*u[3]+u[4]*conj(u[7])*u[7]+u[5]*conj(u[7])*u[6]+u[6]*conj(u[7])*u[5]+u[7]*conj(u[7])*u[4]))/im
A5=((1*((p[3])^2/2)-p[2]/2*abs2(p[1]))*u[5]-p[2]/2*(p[1])^2*conj(u[3]) - (p[2]/2)*(2*p[1]*u[2]*conj(u[1])+conj(p[1])*u[2]*u[7]+2*p[1]*u[3]*conj(u[2])+conj(p[1])*u[3]*u[6]+2*p[1]*u[4]*conj(u[3])+conj(p[1])*u[4]*u[5]+2*p[1]*u[5]*conj(u[4])+conj(p[1])*u[5]*u[4]+2*p[1]*u[6]*conj(u[5])+conj(p[1])*u[6]*u[3]+2*p[1]*u[7]*conj(u[6])+conj(p[1])*u[7]*u[2] + u[1]*conj(u[1])*u[5]+u[2]*conj(u[1])*u[4]+u[3]*conj(u[1])*u[3]+u[4]*conj(u[1])*u[2]+u[5]*conj(u[1])*u[1]+u[1]*conj(u[2])*u[6]+u[2]*conj(u[2])*u[5]+u[3]*conj(u[2])*u[4]+u[4]*conj(u[2])*u[3]+u[5]*conj(u[2])*u[2]+u[6]*conj(u[2])*u[1]+u[1]*conj(u[3])*u[7]+u[2]*conj(u[3])*u[6]+u[3]*conj(u[3])*u[5]+u[4]*conj(u[3])*u[4]+u[5]*conj(u[3])*u[3]+u[6]*conj(u[3])*u[2]+u[7]*conj(u[3])*u[1]+u[2]*conj(u[4])*u[7]+u[3]*conj(u[4])*u[6]+u[4]*conj(u[4])*u[5]+u[5]*conj(u[4])*u[4]+u[6]*conj(u[4])*u[3]+u[7]*conj(u[4])*u[2]+u[3]*conj(u[5])*u[7]+u[4]*conj(u[5])*u[6]+u[5]*conj(u[5])*u[5]+u[6]*conj(u[5])*u[4]+u[7]*conj(u[5])*u[3]+u[4]*conj(u[6])*u[7]+u[5]*conj(u[6])*u[6]+u[6]*conj(u[6])*u[5]+u[7]*conj(u[6])*u[4]+u[5]*conj(u[7])*u[7]+u[6]*conj(u[7])*u[6]+u[7]*conj(u[7])*u[5]))/im
A6=((4*((p[3])^2/2)-p[2]/2*abs2(p[1]))*u[6]-p[2]/2*(p[1])^2*conj(u[2]) - (p[2]/2)*(2*p[1]*u[3]*conj(u[1])+conj(p[1])*u[3]*u[7]+2*p[1]*u[4]*conj(u[2])+conj(p[1])*u[4]*u[6]+2*p[1]*u[5]*conj(u[3])+conj(p[1])*u[5]*u[5]+2*p[1]*u[6]*conj(u[4])+conj(p[1])*u[6]*u[4]+2*p[1]*u[7]*conj(u[5])+conj(p[1])*u[7]*u[3] + u[1]*conj(u[1])*u[6]+u[2]*conj(u[1])*u[5]+u[3]*conj(u[1])*u[4]+u[4]*conj(u[1])*u[3]+u[5]*conj(u[1])*u[2]+u[6]*conj(u[1])*u[1]+u[1]*conj(u[2])*u[7]+u[2]*conj(u[2])*u[6]+u[3]*conj(u[2])*u[5]+u[4]*conj(u[2])*u[4]+u[5]*conj(u[2])*u[3]+u[6]*conj(u[2])*u[2]+u[7]*conj(u[2])*u[1]+u[2]*conj(u[3])*u[7]+u[3]*conj(u[3])*u[6]+u[4]*conj(u[3])*u[5]+u[5]*conj(u[3])*u[4]+u[6]*conj(u[3])*u[3]+u[7]*conj(u[3])*u[2]+u[3]*conj(u[4])*u[7]+u[4]*conj(u[4])*u[6]+u[5]*conj(u[4])*u[5]+u[6]*conj(u[4])*u[4]+u[7]*conj(u[4])*u[3]+u[4]*conj(u[5])*u[7]+u[5]*conj(u[5])*u[6]+u[6]*conj(u[5])*u[5]+u[7]*conj(u[5])*u[4]+u[5]*conj(u[6])*u[7]+u[6]*conj(u[6])*u[6]+u[7]*conj(u[6])*u[5]+u[6]*conj(u[7])*u[7]+u[7]*conj(u[7])*u[6]))/im
A7=((9*((p[3])^2/2)-p[2]/2*abs2(p[1]))*u[7]-p[2]/2*(p[1])^2*conj(u[1]) - (p[2]/2)*(2*p[1]*u[4]*conj(u[1])+conj(p[1])*u[4]*u[7]+2*p[1]*u[5]*conj(u[2])+conj(p[1])*u[5]*u[6]+2*p[1]*u[6]*conj(u[3])+conj(p[1])*u[6]*u[5]+2*p[1]*u[7]*conj(u[4])+conj(p[1])*u[7]*u[4] + u[1]*conj(u[1])*u[7]+u[2]*conj(u[1])*u[6]+u[3]*conj(u[1])*u[5]+u[4]*conj(u[1])*u[4]+u[5]*conj(u[1])*u[3]+u[6]*conj(u[1])*u[2]+u[7]*conj(u[1])*u[1]+u[2]*conj(u[2])*u[7]+u[3]*conj(u[2])*u[6]+u[4]*conj(u[2])*u[5]+u[5]*conj(u[2])*u[4]+u[6]*conj(u[2])*u[3]+u[7]*conj(u[2])*u[2]+u[3]*conj(u[3])*u[7]+u[4]*conj(u[3])*u[6]+u[5]*conj(u[3])*u[5]+u[6]*conj(u[3])*u[4]+u[7]*conj(u[3])*u[3]+u[4]*conj(u[4])*u[7]+u[5]*conj(u[4])*u[6]+u[6]*conj(u[4])*u[5]+u[7]*conj(u[4])*u[4]+u[5]*conj(u[5])*u[7]+u[6]*conj(u[5])*u[6]+u[7]*conj(u[5])*u[5]+u[6]*conj(u[6])*u[7]+u[7]*conj(u[6])*u[6]+u[7]*conj(u[7])*u[7]))/im

@SVector [A1, A2, A3, A4, A5, A6, A7]
end
``````

I tried putting parenthesis around every 4-th term or so to avoid splatting of +() but this did not solve the issue.

Implementation:

Currently, I generate equations as string arrays, write them to a text file, and then load them to the solver. Example code I use for cubic term:

``````using OffsetArrays, ..dnDot_String

function C(N)
Cubic=Vector{String}(undef,6N+1)
oCubic=OffsetArray(Cubic,-3N:3N)
dummy=Array{String}(undef,2N+1,2N+1,6N+1)
odummy=OffsetArray(dummy,-N:N,-N:N,-3N:3N)
fun=Vector{String}(undef,3)
cTuple=NTuple{3,Int64}

for q = -3N:1:3N     #q, n and m are indices appearing in triple summation
for n = -N:1:N
for m = -N:1:N
if abs(q-n+m) <= abs(N)   #condition on indices to select only the terms I need
cTuple = (n,m,q-n+m)
fun = dnDot_String.boolMatch(cTuple,N)
@inbounds odummy[n,m,q] = fun[1]*"*conj("*fun[2]*")*"*fun[3]*"+" #storing term into the array.
elseif abs(q-n+m) > abs(N)
@inbounds odummy[n,m,q] = ""
end
end
end
@inbounds oCubic[q] = chop(string(odummy1[:,:,q]...))::SubString # concantenating the array for each q. This gives vector where q-th component contains all the cubic terms for q-th equation
end
``````

and function `dnDot_String.boolMatch`

``````function boolMatch(n::NTuple{T},N::Int64) where T<:Any
x=lastindex(n)
dummy=Vector{String}(undef,x)
for i=1:x
dummy[i]="u[\$(n[i]+N+1)]"
end
return dummy
end
``````

This essentially generates `u[i]` terms recognised by DifferentialEquations package base on the values of indices in `cTuple`.

So my question is, is there a more efficient way of generating equations in ModellingToolkit or otherwise that would avoid memory overhead I get when trying to solve the system?

Thank you!

My computer is currently occupied, but this should send you in the right direction:

``````using ModelingToolkit, OffsetArrays
N = 5

@variables t u[1:N](t)
dummyeqs = Array{Operation}(undef,2N+1,2N+1,6N+1)
odummyeqs= OffsetArray(dummy,-N:N,-N:N,-3N:3N)
eqs = Array{Equation}(undef,6N+1)
oeqs= OffsetArray(dummy,-3N:3N)

for q = -3N:1:3N     #q, n and m are indices appearing in triple summation
for n = -N:1:N
for m = -N:1:N
if abs(q-n+m) <= abs(N)   #condition on indices to select only the terms I need
cTuple = (n,m,q-n+m)
fun = boolMatch(cTuple,N)
@inbounds odummyeqs[n,m,q] += fun[1]*conj(fun[2])*fun[3] #storing term into the array.
#elseif abs(q-n+m) > abs(N)
#    @inbounds odummy[n,m,q] = ""
end
end
end
@inbounds oeqs[q] = reduce(+,odummyeqs[:,:,q]) # concantenating the array for each q. This gives vector where q-th component contains all the cubic terms for q-th equation
end

function boolMatch(n::NTuple{T},N::Int64) where T<:Any
x=lastindex(n)
[u[(n[i]+N+1)] for i=1:x]
end

sys = ODESystem(eqs,t,u,[])
...
``````

Essentially, its symbolic language works like Julia, so if `x` and `y` are symbols, then `x+y` makes the expression `x+y`.

Hi Chris!

Once again thank you for your input. I have slightly mended your code:

``````using ModelingToolkit, OffsetArrays
N = 3

@variables t u[1:2N+1](t) #2N+1 because sum run from -N to N. So for N normal modes I get 2N+1 equations
dummyeqs = Array{Operation}(undef,2N+1,2N+1,6N+1)
odummyeqs= OffsetArray(dummyeqs,-N:N,-N:N,-3N:3N)
eqs = Array{Equation}(undef,6N+1)
oeqs= OffsetArray(eqs,-3N:3N)
function boolMatch(n::NTuple{T},N::Int64) where T<:Any
x=lastindex(n)
[u[(n[i]+N+1)] for i=1:x]
end

for q = -3N:1:3N     #q, n and m are indices appearing in triple summation
for n = -N:1:N
for m = -N:1:N
if abs(q-n+m) <= abs(N)   #condition on indices to select only the terms I need
cTuple = (n,m,q-n+m)
fun = boolMatch(cTuple,N)
@inbounds odummyeqs[n,m,q] += fun[1]*conj(fun[2])*fun[3] #storing term into the array.
elseif abs(q-n+m) > abs(N)
@inbounds odummyeqs[n,m,q] = 0 # I have brought this back, because I assume reduce() will throw error for undefined reference
end
end
end
@inbounds oeqs[q] = reduce(+,odummyeqs[:,:,q]) # concantenating the array for each q. This gives vector where q-th component contains all the cubic terms for q-th equation
end

sys = ODESystem(eqs,t,u,[])
``````

However, when running this I get

``````ERROR: LoadError: UndefRefError: access to undefined reference
Stacktrace:
[1] getindex at .\array.jl:745 [inlined]
[2] getindex(::OffsetArrays.OffsetArray{ModelingToolkit.Operation,3,Array{ModelingToolkit.Operation,3}}, ::Int64, ::Int64, ::Int64) at C:\Users\visenkon\.julia\packages\OffsetArrays\WmVaB\src\OffsetArrays.jl:150
[3] top-level scope at C:\Users\visenkon\Desktop\chris\chris_cubic.jl:21
[4] include at .\boot.jl:328 [inlined]
[6] include(::Module, ::String) at .\Base.jl:31
[7] exec_options(::Base.JLOptions) at .\client.jl:287
[8] _start() at .\client.jl:460
in expression starting at C:\Users\visenkon\Desktop\chris\chris_cubic.jl:15
``````

which occurs on this line `@inbounds odummyeqs[n,m,q] += fun[1]*conj(fun[2])*fun[3]`
I dont get this when using strings despite updating array in the same way, so I am not sure whats going on?

Thank you!

It might be better to make that

`dummyeqs = zeros(Operation,2N+1,2N+1,6N+1)`

Hi Chris!

Thank you for your reply. This is what I get if I run the code:

``````ERROR: LoadError: MethodError: Cannot `convert` an object of type ModelingToolkit.Operation to an object of type ModelingToolkit.Equation
Closest candidates are:
convert(::Type{T}, !Matched::T) where T at essentials.jl:168
ModelingToolkit.Equation(::ModelingToolkit.Expression, !Matched::ModelingToolkit.Expression) at C:\Users\visenkon\.julia\packages\ModelingToolkit\MohAp\src\equations.jl:10
ModelingToolkit.Equation(::Any, !Matched::Any) at C:\Users\visenkon\.julia\packages\ModelingToolkit\MohAp\src\equations.jl:10
Stacktrace:
[1] setindex!(::Array{ModelingToolkit.Equation,1}, ::ModelingToolkit.Operation, ::Int64) at .\array.jl:782
[2] setindex!(::OffsetArrays.OffsetArray{ModelingToolkit.Equation,1,Array{ModelingToolkit.Equation,1}}, ::ModelingToolkit.Operation, ::Int64) at C:\Users\visenkon\.julia\packages\OffsetArrays\WmVaB\src\OffsetArrays.jl:165
[3] top-level scope at C:\Users\visenkon\Desktop\chris\chris_cubic.jl:27
[4] include at .\boot.jl:328 [inlined]
[6] include(::Module, ::String) at .\Base.jl:31
[7] exec_options(::Base.JLOptions) at .\client.jl:287
[8] _start() at .\client.jl:460
in expression starting at C:\Users\visenkon\Desktop\chris\chris_cubic.jl:15
``````

This happens on this line inside the loop:

``````@inbounds oeqs[q] = reduce(+,odummyeqs[:,:,q])
``````

After reading the documentation on MT I have tried changing the type of `eqs` in the code you provided to Operation. This worked inside the loop. I then went on to define system in the following way:

``````@variables t u[1:2N+1](t)
@derivatives D'~t
@parameters p[1:3]
equations =  [D(u[i]) ~ eqs[i] for i in 1:2N+1]
sys = ODESystem(eqs,t,u,[])
``````

where `eqs` are comprised of linear, quadratic, and cubic terms. When I print the `equations`, I can see that they are correct, i.e the same as I got when using strings. However, when I try and solve equations using the following code, I get a type error:

``````N=3
u1=[0.0+0.0*im for i in 1:2N+1]  # simple initial conditions
u1[Int64((N-1)/2)]=1.0+0.0*im
u1[Int64((N-1)/2)+1]=0.0+0.0*im
u1[Int64((N-1)/2)+2]=1.0+0.0*im
u0=SVector{2N+1}(u1)

function Diff(N::Int64,psi::Float64,g::Float64)
println(tspan)
k=sqrt(g*psi) #defining parameters
p = [psi,g,k]

Analytic_N = N

path=pwd()
path1="C:/Users/visenkon/Desktop/JLD2 for g=\$g, psi=\$psi, N=\$Analytic_N,t=\$t_name"
rm(path1,force=true,recursive=true)
mkdir(path1)
cd(path1)

prob = ODEProblem(sys,u0,tspan,p)
@save "sol_tf=\$(tspan[2]).jld2" sol
cd(path)
return "success"
end
``````

And the error is as follows:

``````ERROR: MethodError: Cannot `convert` an object of type Complex{Float64} to
an object of type Variable
Closest candidates are:
convert(::Type{Variable}, ::Operation) at C:\Users\visenkon\.julia\packages\ModelingToolkit\MohAp\src\ModelingToolkit.jl:72
convert(::Type{T}, ::T) where T at essentials.jl:168
Variable(::Any) at C:\Users\visenkon\.julia\packages\ModelingToolkit\MohAp\src\variables.jl:50
...
Stacktrace:
[1] varmap_to_vars(::SArray{Tuple{7},Complex{Float64},1,7}, ::Array{Variable,1}) at C:\Users\visenkon\.julia\packages\ModelingToolkit\MohAp\src\variables.jl:251
[2] #_#86(::Nothing, ::Bool, ::Bool, ::Bool, ::Bool, ::Bool, ::Bool, ::ModelingToolkit.SerialForm, ::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}, ::Type{ODEProblem{true,tType,isinplace,P,F,K,PT} where PT where K where F where P where isinplace where tType}, ::ODESystem, ::SArray{Tuple{7},Complex{Float64},1,7}, ::Tuple{Float64,Float64}, ::Array{Float64,1}) at C:\Users\visenkon\.julia\packages\ModelingToolkit\MohAp\src\systems\diffeqs\abstractodesystem.jl:218
[3] ODEProblem{true,tType,isinplace,P,F,K,PT} where PT where K where F where P where isinplace where tType(::ODESystem, ::SArray{Tuple{7},Complex{Float64},1,7}, ::Tuple{Float64,Float64}, ::Array{Float64,1}) at C:\Users\visenkon\.julia\packages\ModelingToolkit\MohAp\src\systems\diffeqs\abstractodesystem.jl:215
[4] #ODEProblem#85(::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}, ::Type{ODEProblem}, ::ODESystem, ::SArray{Tuple{7},Complex{Float64},1,7}, ::Vararg{Any,N} where N) at C:\Users\visenkon\.julia\packages\ModelingToolkit\MohAp\src\systems\diffeqs\abstractodesystem.jl:191
[5] ODEProblem(::ODESystem, ::SArray{Tuple{7},Complex{Float64},1,7}, ::Vararg{Any,N} where N) at C:\Users\visenkon\.julia\packages\ModelingToolkit\MohAp\src\systems\diffeqs\abstractodesystem.jl:191
[6] Diff(::Int64, ::Float64, ::Float64) at C:\Users\visenkon\Desktop\chris\chris_solver.jl:30
[7] top-level scope at REPL[2]:1
``````

If I instead pass initial conditions with the following syntax:

``````u0 = [u[1] => 1.0+0.0im
u[2] => 0.0 + 0.0]
``````

and so son, I get:

``````ERROR: ArgumentError: invalid index: nothing of type Nothing
Stacktrace:
[1] to_index(::Nothing) at .\indices.jl:273
[2] to_index(::Array{Float64,1}, ::Nothing) at .\indices.jl:250
[3] to_indices at .\indices.jl:301 [inlined]
[4] to_indices at .\indices.jl:298 [inlined]
[5] setindex! at .\abstractarray.jl:1074 [inlined]
[6] varmap_to_vars(::Array{Pair{Operation,Float64},1}, ::Array{Variable,1}) at C:\Users\visenkon\.julia\packages\ModelingToolkit\MohAp\src\variables.jl:253
[7] #_#86(::Nothing, ::Bool, ::Bool, ::Bool, ::Bool, ::Bool, ::Bool, ::ModelingToolkit.SerialForm, ::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}, ::Type{ODEProblem{true,tType,isinplace,P,F,K,PT} where PT where K where F where P where isinplace where tType}, ::ODESystem, ::Array{Pair{Operation,Complex{Float64}},1}, ::Tuple{Float64,Float64}, ::Array{Pair{Operation,Float64},1}) at C:\Users\visenkon\.julia\packages\ModelingToolkit\MohAp\src\systems\diffeqs\abstractodesystem.jl:219
[8] ODEProblem{true,tType,isinplace,P,F,K,PT} where PT where K where F where P where isinplace where tType(::ODESystem, ::Array{Pair{Operation,Complex{Float64}},1}, ::Tuple{Float64,Float64}, ::Array{Pair{Operation,Float64},1}) at C:\Users\visenkon\.julia\packages\ModelingToolkit\MohAp\src\systems\diffeqs\abstractodesystem.jl:215
[9] #ODEProblem#85(::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}, ::Type{ODEProblem}, ::ODESystem, ::Array{Pair{Operation,Complex{Float64}},1}, ::Vararg{Any,N} where N) at C:\Users\visenkon\.julia\packages\ModelingToolkit\MohAp\src\systems\diffeqs\abstractodesystem.jl:191
[10] ODEProblem(::ODESystem, ::Array{Pair{Operation,Complex{Float64}},1},
::Vararg{Any,N} where N) at C:\Users\visenkon\.julia\packages\ModelingToolkit\MohAp\src\systems\diffeqs\abstractodesystem.jl:191
[11] Diff(::Int64, ::Float64, ::Float64) at C:\Users\visenkon\Desktop\chris\chris_solver.jl:32
[12] top-level scope at REPL[2]:1
``````

Am I defining the initial conditions wrong?

Thank you!

EDIT: on a closer look, it appears that I did not input parameters into

``````sys = ODESystem(eqs,t,u,[])
``````

After fixing this, everything appears to be working, however, I need to run some further tests.
Thank you!

If youâ€™re building it directly from the system, it should be like https://mtk.sciml.ai/dev/highlevel/#Example-1:-Symbolically-Building-an-ODEProblem-for-DifferentialEquations.jl-1