Issue using ModelingToolkit to simulate electrical circuits

Hello everyone,

Lately, I have been working quite extensively with ModelingToolkit, particularly migrating some dynamic control models that I had previously implemented in various closed-source programs (often achieving significantly better performance with MTK by the way!). While I have verified that ModelingToolkit can indeed be used for this purpose without issue, I had never tried linking those control systems with electrical systems.

Yesterday, I began modeling a few simple electrical circuits with ModelingToolkit. At first, everything seemed to work correctly; however, I noticed a behavior that turns out to be rather inconvenient. I will describe it below through an example, as I believe that will make it easier to understand.

Let us take the RC circuit as an example, whose operation and construction are explained in RC Circuit · ModelingToolkitStandardLibrary.jl.

If we implement the circuit, the result is as expected — no issues arise, and the following equations are generated (as shown below).

using PlutoUI
using ModelingToolkit, Plots, OrdinaryDiffEq
using ModelingToolkit: t_nounits as t, D_nounits as D
using ModelingToolkitStandardLibrary
using ModelingToolkitStandardLibrary.Electrical
using ModelingToolkitStandardLibrary.Blocks


@mtkmodel RC begin

	@parameters begin

		R = 1.0
		C = 1.0
		V = 1.0

	end

	@components begin

		resistor = Resistor(R = R)
		capacitor = Capacitor(C = C, v = 0.0)
		source = Voltage()
		constant = Constant(k = V)
		ground = Ground()
		
		
	end

    @equations begin
		
        connect(constant.output, source.V)
        connect(source.p, resistor.p)
        connect(resistor.n, capacitor.p)
        connect(capacitor.n, source.n, ground.g)
		
    end
end

@mtkbuild sys = RC()

begin
	prob = ODEProblem(sys, Pair[], (0, 10.0))
	sol = solve(prob)
end
\frac{\mathrm{dcapacitor.} v(t)}{\mathrm{d}t} = \frac{\mathrm{capacitor.} i(t)}{\mathrm{capacitor.} C}

But what happens if we express the resistor as two resistors in series whose total value equals that of the original resistor? In principle, one would expect the system to behave the same, and therefore the resulting system of equations to be identical and to behave in the same way. However, when we define the model and generate the equations we can see that a new equation is created. Moreover, when solving the system the following error appears:

@mtkmodel RCAC2 begin

	@parameters begin

		R1 = 0.7
		R2 = 0.3
		C = 1.0
		V = 1.0
		

	end

	@components begin

		resistor1 = Resistor(R = R1)
		resistor2 = Resistor(R = R2)		
		capacitor = Capacitor(C = C, v = 0.0)
		source = Voltage()
		constant = Constant(k = V)
		ground = Ground()
		
		
	end

    @equations begin
		
        connect(constant.output, source.V)
        connect(source.p, resistor1.p)
        connect(resistor1.n, resistor2.p)
		connect(resistor2.n,capacitor.p)
        connect(capacitor.n, source.n, ground.g)
		
    end
end

@mtkbuild sysAC2 = RCAC2()


probAC2 = ODEProblem(sysAC2, Pair[], (0, 10.0))
solAC2 = solve(probAC2)

\begin{align*} \frac{\mathrm{dcapacitor.} v(t)}{\mathrm{d}t} &= \frac{\mathrm{capacitor.} i(t)}{\mathrm{capacitor.} C} \\[1em] 0 &= \mathrm{resistor1.} i(t) - \mathrm{capacitor.} i(t) \end{align*}

Cyclic guesses detected in the system. Symbolic values were found for the following variables/parameters in the map: resistor1₊v(t) => -0.3capacitor₊i(t) In order to resolve this, please provide additional numeric guesses so that the chain can be resolved to assign numeric values to each variable. "

This error can be easily solved by assigning initial values to one of the resistors. In this case, that is straightforward, but of course, in more complex circuits this becomes considerably more complicated. Moreover, I have verified through other examples that if the initial value does not match the actual value at the start, the system does not solve correctly. Given this situation, I wonder whether there is something wrong with my code (most likely due to my being a beginner in this area), or if there is something I am overlooking when modeling electrical circuits in ModelingToolkit. Obviously, in this particular case the issue is not serious, since the two resistors can simply be combined into a single one without any problem. However, in cases such as two transistors connected in series, this issue is not so easily avoided.

If possible, I would appreciate if someone could provide an explanation or offer some advice on how to prevent this behavior.

Thanks in advance.

@cryptic.ax isn’t this the case where you should have a linear SCC so the guess isn’t required?

No initial conditions are provided. I’d hazard a guess that the initialization is underdetermined. Unfortunately, that appears as a warning and then the absence of guesses shows up as an error because it doesn’t use SCCNonlinearProblem.

So, if I understand correctly, the issue arises because the system is underdetermined during initialization when no initial conditions are provided, and the solver does not employ SCCNonlinearProblem in this case.

While assigning manually initial conditions to each element solves the issue, this approach can become quite overwhelming when dealing with larger circuits. Is there any recommended workaround or suggested approach to handle this more efficiently?

Any non-empty ODE/DAE requires some initial conditions to solve. They can be provided during model construction as defaults, or to the second argument of ODEProblem as a Dict or list of pairs mapping variables to values. You should only need to provide as many initial conditions as differential equations/variables.

With electrical circuits, you could assume all currents and voltages are zero at t=0. This is easy to do automatically, even for very large circuits.

But then you need to provide a ramp for the supply voltage of your circuit. Also, not difficult, and it scales well (worst case, you need, e.g., five supply voltages, and perhaps turn them on in a specific order).

This will not work for every circuit. But ramping the supplies & inputs from 0 is indeed a convenient way to get easy initialization. I mean to open an issue to do this homotopy automatically.

Another way to get a consistent initialization is to simulate the circuit in SPICE, which is typically more lenient and use the starting state from there as an initial guess in MTK.

Assuming all currents are zero at t = 0 seemed like a convenient way to obtain a good initialization (as @ufechner7 pointed out). However, after doing some small experiments, this method does not always work and can lead to poor initializations as well. I will try ramping up all inputs, as both @ufechner7 and @GeorgeGkountouras suggested. Nevertheless, I am quite surprised by this behavior. One thing I really like about ModelingToolkit is that everything usually works as expected quite easily (and extremely efficiently), so I was surprised that ModelingToolkit does not solve this initialization problem automatically.