RC circuit fails

[muendlein]

I am using a minor modification of RC circuit example by have the resistor and capacitor in parallel instead of in series. However this fails with ERROR: MethodError: no method matching oneunit(::Type{Any}) during the ODE setup.
Any idea what is going on?

using ModelingToolkit, OrdinaryDiffEq, Plots
using ModelingToolkitStandardLibrary.Electrical
using ModelingToolkitStandardLibrary.Blocks: Constant

R = 1.0
C = 1.0
V = 1.0
L = 1.0
@variables t
@named resistor = Resistor(R=R)
@named capacitor = Capacitor(C=C)
@named source = Voltage()
@named constant = Constant(k=V)
@named ground = Ground()

rc_eqs = [
        connect(constant.output, source.V)
        connect(source.p, resistor.p)
        connect(source.p, capacitor.p)
        connect(resistor.n, ground.g)
        connect(capacitor.n, ground.g)
        connect(source.n, ground.g)
        ]

@named rc_model = ODESystem(rc_eqs, t, systems=[resistor, capacitor, constant, source, ground])
sys = structural_simplify(rc_model)
prob = ODAEProblem(sys, Pair[], (0, 10.0))
sol = solve(prob, Tsit5())
plot(sol, vars = [capacitor.v, resistor.i],
     title = "RC Circuit Demonstration",
     labels = ["Capacitor Voltage" "Resistor Current"])

[ValentinKaisermayer]

Connecting a capacitor and a voltage source in parallel will result in numerical difficulties. The only solution is that the capacitor is already charged to the voltage of the source.

I would guess it results in an implicit DAE.
Cc @YingboMa might be a nice test example

using ModelingToolkit, OrdinaryDiffEq, Plots
using ModelingToolkitStandardLibrary.Electrical
using ModelingToolkitStandardLibrary.Blocks: Constant

R = 1.0
C = 1.0
V = 1.0
L = 1.0
@parameters t
D = Differential(t)
@named resistor = Resistor(R=R)
@named capacitor = Capacitor(C=C)
@named source = Voltage()
@named constant = Constant(k=V)
@named ground = Ground()

rc_eqs = [
        connect(constant.output, source.V)
        connect(source.p, resistor.p, capacitor.p)
        connect(resistor.n, capacitor.n, source.n, ground.g)
        ]

@named rc_model = ODESystem(rc_eqs, t, systems=[resistor, capacitor, constant, source, ground])
sys = structural_simplify(rc_model)
prob = DAEProblem(sys, (D.(states(sys)) .=> 0.0), [D(source.V.u) => 0.0], (0, 10.0))
sol = solve(prob, DFBDF())
plot(sol, vars = [capacitor.v, resistor.i],
     title = "RC Circuit Demonstration",
     labels = ["Capacitor Voltage" "Resistor Current"])

[muendlein]

@ValentinKaisermayer I don't really see why connecting a capacitor and voltage source in parallel will cause any stability issues. After all that is just simply charging up a capacitor..
Secondly that does not explain the observed error.

[ValentinKaisermayer]

In reality this works because a capacitor has parasitic resistances. But this ideal model has not.

The error message is unclear to me as well. I guess the reason is that after structual simplify the model is no longer an ODE but rather only algebraic equations. So there is no initial state vector, which Tsit5 can not handle.

[muendlein]

@ValentinKaisermayer Isn't that more of issue of the voltage source rather than the capacitor? Meaning that due to the lack of an internal resistance the current can become infinity?
As far as I know it is not possible to specify an internal resistance of a voltage source, or is that already implemented?

[baggepinnen]

Well, the voltage source is an ideal source, and the capacitor is also ideal, and having both of them being ideal is.. not ideal. If you indeed want to model a parallel connection of two such elements. adding a resistance manually should be trivial. Another option would be to add internal_resistance = 0 as options to all ideal components, which would allow the user to select a non-zero resistance, but it's probably simpler to explicitly model the internal resistance if you need it?

[ValentinKaisermayer]

Making a RealCapacitor with ESR is trivial. If it is part of the library then as separate component.

[ValentinKaisermayer]

Something like this:

function RealCapacitor(;name, C, ESR, ESL, R_leak)
    @named capacitor = Capacitor(C=C)
    @named esr = Resistor(R=ESR)
    @named esl = Inductor(L=ESL)
    @named leakage_res = Resistor(R=R_leak)

    @named p = Pin()
    @named n = Pin()

    connections = [
        connect(p, esr.p)
        connect(esr.n, capacitor.p, leakage_res.p)
        connect(capacitor.n, leakage_res.n, esl.p)
        connect(esl.n, n)
    ]
    return ODESystem(connections, t; systems = [capacitor, esr, esl, leakage_res], name=name)
end

[YingboMa]

This issue seems to be resolved on the latest MTK.