Ideal Diode and Ideal Switch instability problem

Greetings,

I am new to Julia community and, in particular, to the use of ModelingToolkit package. I would like to implement an IdealDiode component to be replaced by a short circuit when the controlled variable (i.e. integral of the current flowing into the Diode) has reached one of the two thresholds (being either 0 or Vin_th).

Initially, i had thought of implementing a model for a simple ideal diode and a model of an ideal controlled switch. Such devices would be connected in parallel in the final circuit.

I had faced some difficulties:

  • Ideal Diode returns instability (dt < dtmin) probably dued to the rising edge. This is why I have decided to go for a Shockley Diode.
  • Ideal Switch works this way: if controlled variable is <= 0 or >= Vin_th, the switch is closed, so that the diode is excluded. The switch is otherwise open. I am getting the very same error. I decided to introduce a tanh function in the current equation of such device in order to reduce this behavior.

These workarounds seem not to have worked when simulating, considering the full circuit. In some Threads on Discourse and in the Documentation, I’ve read about Callbacks but I haven’t understood how to write them properly and where to insert them in the DSL code and if there exist some macro helping in this process. Should I drop the usage of new macros (i.e. @mtkmodel, @mtkbuild 
) for the old format in order to add Callbacks?

Another possibility could be to use conditionals, in order to include/exclude the IdealDiode in the circuit but i don’t think this will solve instability, in the end. I write here to ask some suggestion to avoid workarounds and for some help facing these problems directly.

This is the code i used for modeling. Forgive my English, it’s not my first language. Thanks in advance for any feedback.

Used libraries

using ModelingToolkit
using ModelingToolkitStandardLibrary.Blocks
using ModelingToolkitStandardLibrary.Electrical
using OrdinaryDiffEq
using Plots

Diode

@mtkmodel Diode begin
    @extend v, i = oneport = OnePort()
    @parameters begin
        threshold = 0, [description = "Voltage Threshold"]
    end
    @components begin
        input = RealInput()
    end
    @equations begin
        0 ~ ifelse(input.u > 0, v, i)
    end
end

ver. 2

exlin(x, max_x) = ifelse(x > max_x,
                         exp(max_x)*(1 + x - max_x),
                         exp(x))
@mtkmodel Diode begin
    @extend v, i = oneport = OnePort()
    @parameters begin
        Ids     = 1e-6,
        [description = "Reverse-bias current"]
        max_exp = 15,
        [description = "Value after which linearization is applied"]
        R       = 1e8,
        [description = "Diode Resistance"]
        Vth     = 1e-3,
        [description = "Threshold voltage"]
        k       = 1e3,
        [description = "Speed of exponential"]
    end
    @equations begin
        i ~ Ids * (exlin(k * (v - Vth) / Vth, max_exp) - 1) + (v / R)
    end
end

Switch:

@mtkmodel Switch begin
    @extend v, i = oneport = OnePort()
    @parameters begin
        V_FRC, [description = "Airway Volume at FRC"]
    end
    @variables begin
        ∫i(t) = 0, [description = "Current integral"]
    end
    @equations begin
        0 ~ ifelse(∫i / V_FRC <= 0, v,
                   ifelse(∫i / V_FRC < 0.9, i, v))
    end
end

Switch ver. 2

@mtkmodel Switch begin
    @extend v, i = oneport = OnePort()
    @parameters begin
        V_FRC, [description = "Airway Volume at FRC"]
        Rclosed = 1e-12, [description = "Switch Resistance when Closed"]
        Ropen = 2.5e5, [description = "Switch Resistance when Open"]
        k = 1e3
    end
    @variables begin
        ∫i(t) = 0, [description = "Current integral"]
        R(t) = 0,  [description = "Switch Resistance"]
    end
    @equations begin
        R ~ Rclosed + (Ropen - Rclosed) * (1 / 2) * (tanh(k * (∫i / V_FRC) * (1 - (∫i / V_FRC))))
        v ~ R * i
    end
end

For performing the simulation I’ve used:

@mtkbuild system = System()
prob = ODEProblem(system, Pair[], (0, 2));
sol = solve(prob, Rodas4());

You can find one (complex) example how to use callbacks in my tutorial: Examples · Tethers.jl

In priciple you just have to define:
@named sys = ODESystem(eqs, t; continuous_events = cb)

and cb must be an array of equations that passes zero when your switches switch


What I am doing in my code is to create an array of equations, because I have segment elements that can “switch”.

Did you read the documentation: Event Handling and Callback Functions · ModelingToolkit.jl ?

Hi, can you elaborate on this Diode model. Is it specific to your circuit? Because I don’t see in general the need to account for the integral of the current to model the switching transition of a Diode.

That being said, I’m familiar with Modelica’s Diode model (based on switching between large and small resistance values), but have no experience in doing so in ModelingToolkit.

For the reference, these are the simplified, compactified equations of Modelica Standard Libary’s 4.0 IdealDiode:

parameter Real Ron = 1e-5
    "Forward state-on differential resistance (closed resistance)";
parameter Real Goff = 1e-5
    "Backward state-off conductance (opened conductance)";

 Real v "voltage drop";
 Real i "current";
 Boolean off(start=true) "Switching state";
 Real s "Auxiliary variable for actual position on the ideal diode characteristic";
  /* s = 0: knee point
     s < 0: below knee point, blocking
     s > 0: above knee point, conducting */

equation

  v = s * (if off then 1 else Ron);
  i = s * (if off then Goff else 1) ;
  off = s < 0;

Would this approach be implementable in MTK?

References:

Hi! Has there been any progress on this? I would also like to use a diode in ModelingToolkit.jl

The Modelica implementation seems fine: Diode - Simple diode with heating port

If there hasn’t been any progress then I will try to make a PR for it. Shouldn’t be that hard to do.

Do you need discontinuity detection for this?

No, if we use the Shockley diode equation like Modelica then it is continuous.

For reference: Shockley diode equation

Yeah if there’s no derivative discontinuities it should be implementable, I just don’t think anyone has done it yet.

1 Like

I implemented the diode.

There isn’t a place to put common constants yet so I just put them inside the component. If you want I can also make a separate PR for that?

I also implemented the thermal diode which has an additional thermalport and models the temperature dependent behaviour, I’ll make a separate PR for that once I have added tests.

"""
        Diode(; name, Is = 1e-12, n = 1)

Ideal diode based on the Shockley diode equation.

# States

    - See [OnePort](@ref)

# Connectors
    
    - `p` Positive pin
    - `n` Negative pin

# Parameters
     
    - `Is`: [`A`] Saturation current
    - `n`: Ideality factor
"""
@mtkmodel Diode begin
    begin
        k = 1.380649e-23 # Boltzmann constant (J/K)
        q = 1.602176634e-19 # Elementary charge (C)
    end

    @extend v, i = oneport = OnePort(; v = 0.0)
    @parameters begin
        Is = 1e-6, [description = "Saturation current (A)"]
        n = 1, [description = "Ideality factor"]
        T = 300.15, [description = "Ambient temperature"]
    end
    @equations begin
        i ~ Is * (exp(v * q / (n * k * T)) - 1)
    end
end

Made a small demo with an AC source:

1 Like

PR here: Add diode component by langestefan · Pull Request #343 · SciML/ModelingToolkitStandardLibrary.jl · GitHub