[help] Implementing a basic NPN BJT model

[RagibHasin]

While fiddling with Modia, I tried to implement a basic NPN BJT model referring to the MATLAB code from (https://gist.github.com/bencholmes/d3ce40d1e9450f692320f68fd0a239f7) like this

BJT_NPN = Model(
    b = Pin,
    c = Pin,
    e = Pin,
    
    Ids = 1.e-6u"A",
    Vt = 0.04u"V",
    gainForward = 300.0,
    gainReverse = 120.0,
    equations = :[
        alphaR = (gainReverse + 1)/gainReverse
        Vbe = b.v - e.v
        Vbc = b.v - c.v
        If = Ids*(exp(Vbe/Vt) - 1)
        Ir = Ids*(exp(Vbc/Vt) - 1)
        c.i = If - Ir*alphaR
        b.i = If/gainForward + Ir/gainReverse
        e.i = c.i + b.i
    ]
)

Using it like this:

CommonEmitter= Model(
    Vbb = ConstantVoltage | Map(V=5.0u"V"),
    Vcc = ConstantVoltage | Map(V=15.0u"V"),
    RC = Resistor | Map(R=10e3u"Ω"),
    RB = Resistor | Map(R= 1e3u"Ω"),
    X = BJT_NPN,
    ground = Ground,

    connect = :[
        (Vcc.p, RC.p)
        (Vbb.p, RB.p)
        (X.c, RC.n)
        (X.b, RB.n)
        (ground.p, Vcc.n, Vbb.n, X.e)
    ]
)

commonEmitter1 = @instantiateModel(CommonEmitter)
simulate!(commonEmitter1, stopTime = 1.0u"s")
plot(commonEmitter1, [("X.b.v", "X.c.v", "X.e.v"]), ("X.b.i", "X.c.i", "X.e.i")], figure=1)

gives the following message

Information message from getSortedAndSolvedAST for model CommonEmitter:
The following variables are iteration variables but have no start/init values defined.
If units are used in the model, start/init values with correct units should be defined
to avoid unit errors during compilation.
Involved variables:
    X.Vbe
    RC.n.v

ERROR: LoadError: DimensionError: 0.0 and 0.0 V are not dimensionally compatible.
Stacktrace:
  [1] +(x::Quantity{Float64, NoDims, Unitful.FreeUnits{(), NoDims, nothing}}, y::Quantity{Float64, 𝐋^2 𝐌 𝐈^-1 𝐓^-3, Unitful.FreeUnits{(V,), 𝐋^2 𝐌 𝐈^-
1 𝐓^-3, nothing}} )

Some help pointing in the right direction would be very helpful.

[MartinOtter]

As a quick fix: use @instantiateModel(CommonEmitter, unitless=true) to switch off unit handling for code generation.

There was a warning message from @instantiateModel:

Information message from getSortedAndSolvedAST for model CommonEmitter:
The following variables are iteration variables but have no start/init values defined.
If units are used in the model, start/init values with correct units should be defined
to avoid unit errors during compilation.
Involved variables:
    X.Vbe
    RC.n.v

This means that the mentioned variables need to be explicitly declared with units, such as:

CommonEmitter= Model(
   ...
   Vbe = Var(start=0.0u"V")

Unfortunately, RC.n.v from model library Modia/models/Electric.jl has no unit defined. This will be fixed at some time in the future. In the meantime, please use the unitless=true option.

[RagibHasin]

Thanks a lot. Your suggestion seems to work. But then another model using this errors with NaN value.

Code:

Amplifier1 = Model(
    Vcc = ConstantVoltage | Map(V=22.0u"V"),
    R1 = Resistor | Map(R=56e3u"Ω"),
    R2 = Resistor | Map(R=8.2e3u"Ω"),
    RC = Resistor | Map(R=6.8e3u"Ω"),
    RE = Resistor | Map(R=1.5e3u"Ω"),
    Ci = Capacitor | Map(C=10e-6u"F", v=Var(init=0.0u"V")),
    Co = Capacitor | Map(C=10e-6u"F", v=Var(init=0.0u"V")),
    CE = Capacitor | Map(C=20e-6u"F", v=Var(init=0.0u"V")),
    Vi = SineVoltage | Map(V=8e-3u"V", f=1.5e3u"Hz"),
    Vo = VoltageSensor,
    X = BJT_NPN,
    ground = Ground,
    connect = :[
      (Vcc.p , R1.p, RC.p)
      (Vi.p, Ci.p)
      (R1.n, R2.p, Ci.n, X.b)
      (RC.n, Co.p, X.c)
      (RE.p, CE.p, X.e)
      (Vo.p, Co.n)
      (ground.p, R2.n, RE.n, CE.n, Vcc.n, Vi.n, Vo.n)
    ]
)

amplifier1 = @instantiateModel(Amplifier1, unitless=true, log=true)
simulate!(amplifier1, Tsit5(), stopTime = 5e-3u"s", log=true) # runs fine up to 0.00668s

Output:

Instantiating model Amplifier1
  in module: Main
  in file: #$#$#$
Connect equations in Amplifier1:
  Vcc.p.v = R1.p.v
  R1.p.v = RC.p.v
  0 = (Vcc.p.i + R1.p.i) + RC.p.i
  Vi.p.v = Ci.p.v
  0 = Vi.p.i + Ci.p.i
  R1.n.v = R2.p.v
  R2.p.v = Ci.n.v
  Ci.n.v = X.b.v
  0 = ((R1.n.i + R2.p.i) + Ci.n.i) + X.b.i
  RC.n.v = Co.p.v
  Co.p.v = X.c.v
  0 = (RC.n.i + Co.p.i) + X.c.i
  RE.p.v = CE.p.v
  CE.p.v = X.e.v
  0 = (RE.p.i + CE.p.i) + X.e.i
  Vo.p.v = Co.n.v
  0 = Vo.p.i + Co.n.i
  ground.p.v = R2.n.v
  R2.n.v = RE.n.v
  RE.n.v = CE.n.v
  CE.n.v = Vcc.n.v
  Vcc.n.v = Vi.n.v
  Vi.n.v = Vo.n.v
  0 = (((((ground.p.i + R2.n.i) + RE.n.i) + CE.n.i) + Vcc.n.i) + Vi.n.i) + Vo.n.i
Number of states:    3
Number of equations: 70
Parameters:
   1: :(Vcc.V) => 22.0
   2: :(R1.R) => 56000.0
   3: :(R2.R) => 8200.0
   4: :(RC.R) => 6800.0
   5: :(RE.R) => 1500.0
   6: :(Ci.C) => 1.0e-5
   7: :(Co.C) => 1.0e-5
   8: :(CE.C) => 2.0e-5
   9: :(Vi.V) => 0.008
  10: :(Vi.f) => 1500.0
  11: :(X.Ids) => 1.0e-6
  12: :(X.Vt) => 0.04
  13: :(X.gainForward) => 300.0
  14: :(X.gainReverse) => 120.0
Potential states:
   1: Ci.v
   2: Co.v
   3: CE.v
Unknowns:
   1: Vcc.v
   2: Vcc.p.v
   3: Vcc.n.v
   4: Vcc.p.i
   5: Vcc.n.i
   6: Vcc.i
   7: R1.v
   8: R1.p.v
   9: R1.n.v
  10: R1.p.i
  11: R1.n.i
  12: R1.i
  13: R2.v
  14: R2.p.v
  15: R2.n.v
  16: R2.p.i
  17: R2.n.i
  18: R2.i
  19: RC.v
  20: RC.p.v
  21: RC.n.v
  22: RC.p.i
  23: RC.n.i
  24: RC.i
  25: RE.v
  26: RE.p.v
  27: RE.n.v
  28: RE.p.i
  29: RE.n.i
  30: RE.i
  31: Ci.p.v
  32: Ci.n.v
  33: Ci.p.i
  34: Ci.n.i
  35: Ci.i
  36: der(Ci.v)
  37: Co.p.v
  38: Co.n.v
  39: Co.p.i
  40: Co.n.i
  41: Co.i
  42: der(Co.v)
  43: CE.p.v
  44: CE.n.v
  45: CE.p.i
  46: CE.n.i
  47: CE.i
  48: der(CE.v)
  49: Vi.v
  50: Vi.p.v
  51: Vi.n.v
  52: Vi.p.i
  53: Vi.n.i
  54: Vi.i
  55: Vo.p.i
  56: Vo.n.i
  57: Vo.v
  58: Vo.p.v
  59: Vo.n.v
  60: X.alphaR
  61: X.If
  62: X.b.v
  63: X.e.v
  64: X.Ir
  65: X.c.v
  66: X.c.i
  67: X.b.i
  68: X.e.i
  69: ground.p.v
  70: ground.p.i
Equations:
   1: Vcc.v = Vcc.p.v - Vcc.n.v
   2: 0 = Vcc.p.i + Vcc.n.i
   3: Vcc.i = Vcc.p.i
   4: Vcc.v = Vcc.V
   5: R1.v = R1.p.v - R1.n.v
   6: 0 = R1.p.i + R1.n.i
   7: R1.i = R1.p.i
   8: R1.R * R1.i = R1.v
   9: R2.v = R2.p.v - R2.n.v
  10: 0 = R2.p.i + R2.n.i
  11: R2.i = R2.p.i
  12: R2.R * R2.i = R2.v
  13: RC.v = RC.p.v - RC.n.v
  14: 0 = RC.p.i + RC.n.i
  15: RC.i = RC.p.i
  16: RC.R * RC.i = RC.v
  17: RE.v = RE.p.v - RE.n.v
  18: 0 = RE.p.i + RE.n.i
  19: RE.i = RE.p.i
  20: RE.R * RE.i = RE.v
  21: Ci.v = Ci.p.v - Ci.n.v
  22: 0 = Ci.p.i + Ci.n.i
  23: Ci.i = Ci.p.i
  24: Ci.C * der(Ci.v) = Ci.i
  25: Co.v = Co.p.v - Co.n.v
  26: 0 = Co.p.i + Co.n.i
  27: Co.i = Co.p.i
  28: Co.C * der(Co.v) = Co.i
  29: CE.v = CE.p.v - CE.n.v
  30: 0 = CE.p.i + CE.n.i
  31: CE.i = CE.p.i
  32: CE.C * der(CE.v) = CE.i
  33: Vi.v = Vi.p.v - Vi.n.v
  34: 0 = Vi.p.i + Vi.n.i
  35: Vi.i = Vi.p.i
  36: Vi.v = Vi.V * sin(2 * 3.14159 * Vi.f * time)
  37: Vo.p.i = Vo.n.i
  38: 0 = Vo.p.i
  39: Vo.v = Vo.p.v - Vo.n.v
  40: X.alphaR = (X.gainReverse + 1) / X.gainReverse
  41: X.If = X.Ids * (exp((X.b.v - X.e.v) / X.Vt) - 1)
  42: X.Ir = X.Ids * (exp((X.b.v - X.c.v) / X.Vt) - 1)
  43: X.c.i = X.If - X.Ir * X.alphaR
  44: X.b.i = X.If / X.gainForward + X.Ir / X.gainReverse
  45: X.e.i = X.c.i + X.b.i
  46: ground.p.v = 0.0 * 1
  47: Vcc.p.v = R1.p.v
  48: R1.p.v = RC.p.v
  49: 0 = (Vcc.p.i + R1.p.i) + RC.p.i
  50: Vi.p.v = Ci.p.v
  51: 0 = Vi.p.i + Ci.p.i
  52: R1.n.v = R2.p.v
  53: R2.p.v = Ci.n.v
  54: Ci.n.v = X.b.v
  55: 0 = ((R1.n.i + R2.p.i) + Ci.n.i) + X.b.i
  56: RC.n.v = Co.p.v
  57: Co.p.v = X.c.v
  58: 0 = (RC.n.i + Co.p.i) + X.c.i
  59: RE.p.v = CE.p.v
  60: CE.p.v = X.e.v
  61: 0 = (RE.p.i + CE.p.i) + X.e.i
  62: Vo.p.v = Co.n.v
  63: 0 = Vo.p.i + Co.n.i
  64: ground.p.v = R2.n.v
  65: R2.n.v = RE.n.v
  66: RE.n.v = CE.n.v
  67: CE.n.v = Vcc.n.v
  68: Vcc.n.v = Vi.n.v
  69: Vi.n.v = Vo.n.v
  70: 0 = (((((ground.p.i + R2.n.i) + RE.n.i) + CE.n.i) + Vcc.n.i) + Vi.n.i) + Vo.n.i
Perform alias elimination and remove singularities

    Variables that have been eliminated:
      66: X.e.v = CE.v
      24: RC.i = X.c.i
      13: R2.v = X.b.v
      25: RE.v = CE.v
      43: Co.i = 0
      17: R2.n.i = -R2.i
      16: R2.p.i = R2.i
      11: R1.n.i = -R1.i
      61: Vo.p.v = Co.n.v
      23: RC.n.i = -X.c.i
      39: Co.p.v = X.c.v
      21: RC.n.v = X.c.v
      33: Ci.n.v = X.b.v
      14: R2.p.v = X.b.v
       9: R1.n.v = X.b.v
      56: Vi.n.i = Ci.i
      32: Ci.p.v = Vi.v
      20: RC.p.v = Vcc.v
       8: R1.p.v = Vcc.v
       6: Vcc.i = -Vcc.n.i
      60: Vo.v = Co.n.v
       2: Vcc.p.v = Vcc.v
       4: Vcc.p.i = -Vcc.n.i
      57: Vi.i = -Ci.i
      55: Vi.p.i = -Ci.i
      53: Vi.p.v = Vi.v
      49: CE.n.i = -CE.i
      48: CE.p.i = CE.i
      46: CE.p.v = CE.v
      10: R1.p.i = R1.i
      35: Ci.n.i = -Ci.i
      34: Ci.p.i = Ci.i
      29: RE.n.i = -RE.i
      28: RE.p.i = RE.i
      26: RE.p.v = CE.v
      22: RC.p.i = X.c.i
      62: Vo.n.v = 0
      54: Vi.n.v = 0
       3: Vcc.n.v = 0
      47: CE.n.v = 0
      27: RE.n.v = 0
      15: R2.n.v = 0
      41: Co.p.i = 0
      42: Co.n.i = 0
      72: ground.p.v = 0
      59: Vo.n.i = 0
      58: Vo.p.i = 0

    Remaining transformed linear Integer equations:
      25: 0 = Co.n.v + Co.v - X.c.v
      61: 0 = -X.e.i - RE.i - CE.i
      49: 0 = -Vcc.n.i + X.c.i + R1.i
      70: 0 = -ground.p.i + RE.i + CE.i - Ci.i + R2.i - X.c.i - R1.i
      55: 0 = X.b.i - Ci.i - R1.i + R2.i
      45: 0 = X.e.i - X.c.i - X.b.i
      13: 0 = RC.v - Vcc.v + X.c.v
      21: 0 = Vi.v - Ci.v - X.b.v
       5: 0 = -R1.v + Vcc.v - X.b.v

States:
   1: Ci.v
   2: Co.v
   3: CE.v
Unknowns after alias reduction:
   1: Vcc.v
   2: Vcc.n.i
   3: R1.v
   4: R1.i
   5: R2.i
   6: RC.v
   7: RE.i
   8: Ci.i
   9: der(Ci.v)
  10: Co.n.v
  11: der(Co.v)
  12: CE.i
  13: der(CE.v)
  14: Vi.v
  15: X.alphaR
  16: X.If
  17: X.b.v
  18: X.Ir
  19: X.c.v
  20: X.c.i
  21: X.b.i
  22: X.e.i
  23: ground.p.i
Equations after alias reduction:
   1: Vcc.v = Vcc.V
   2: 0 = (-(R1.v) + Vcc.v) + -(X.b.v)
   3: R1.R * R1.i = R1.v
   4: R2.R * R2.i = X.b.v
   5: 0 = (RC.v + -(Vcc.v)) + X.c.v
   6: RC.R * X.c.i = RC.v
   7: RE.R * RE.i = CE.v
   8: 0 = (Vi.v + -(Ci.v)) + -(X.b.v)
   9: Ci.C * der(Ci.v) = Ci.i
  10: 0 = (Co.n.v + Co.v) + -(X.c.v)
  11: Co.C * der(Co.v) = 0.0
  12: CE.C * der(CE.v) = CE.i
  13: Vi.v = Vi.V * sin(2 * 3.14159 * Vi.f * time)
  14: X.alphaR = (X.gainReverse + 1) / X.gainReverse
  15: X.If = X.Ids * (exp((X.b.v - CE.v) / X.Vt) - 1)
  16: X.Ir = X.Ids * (exp((X.b.v - X.c.v) / X.Vt) - 1)
  17: X.c.i = X.If - X.Ir * X.alphaR
  18: X.b.i = X.If / X.gainForward + X.Ir / X.gainReverse
  19: 0 = (X.e.i + -(X.c.i)) + -(X.b.i)
  20: 0 = (-(Vcc.n.i) + X.c.i) + R1.i
  21: 0 = ((X.b.i + -(Ci.i)) + -(R1.i)) + R2.i
  22: 0 = (-(X.e.i) + -(RE.i)) + -(CE.i)
  23: 0 = (((((-(ground.p.i) + RE.i) + CE.i) + -(Ci.i)) + R2.i) + -(X.c.i)) + -(R1.i)
  The DAE is structurally nonsingular.

Sorted highest derivative equations:
[assigned variable]: [differentiation] equation
Strongly connected components are enclosed in []
   1:                                             Vcc.v:  Vcc.v = Vcc.V
  13:                                              Vi.v:  Vi.v = Vi.V * sin(2 * 3.14159 * Vi.f * time)
   8:                                             X.b.v:  0 = (Vi.v + -(Ci.v)) + -(X.b.v)
   2:                                              R1.v:  0 = (-(R1.v) + Vcc.v) + -(X.b.v)
   3:                                              R1.i:  R1.R * R1.i = R1.v
   4:                                              R2.i:  R2.R * R2.i = X.b.v
  15:                                              X.If:  X.If = X.Ids * (exp((X.b.v - CE.v) / X.Vt) - 1)
  14:                                          X.alphaR:  X.alphaR = (X.gainReverse + 1) / X.gainReverse
[
  16:                                              X.Ir:  X.Ir = X.Ids * (exp((X.b.v - X.c.v) / X.Vt) - 1)
  17:                                             X.c.i:  X.c.i = X.If - X.Ir * X.alphaR
   6:                                              RC.v:  RC.R * X.c.i = RC.v
   5:                                             X.c.v:  0 = (RC.v + -(Vcc.v)) + X.c.v
]
   7:                                              RE.i:  RE.R * RE.i = CE.v
  18:                                             X.b.i:  X.b.i = X.If / X.gainForward + X.Ir / X.gainReverse
  21:                                              Ci.i:  0 = ((X.b.i + -(Ci.i)) + -(R1.i)) + R2.i
   9:                                         der(Ci.v):  Ci.C * der(Ci.v) = Ci.i
  10:                                            Co.n.v:  0 = (Co.n.v + Co.v) + -(X.c.v)
  11:                                         der(Co.v):  Co.C * der(Co.v) = 0.0
  19:                                             X.e.i:  0 = (X.e.i + -(X.c.i)) + -(X.b.i)
  22:                                              CE.i:  0 = (-(X.e.i) + -(RE.i)) + -(CE.i)
  12:                                         der(CE.v):  CE.C * der(CE.v) = CE.i
  20:                                           Vcc.n.i:  0 = (-(Vcc.n.i) + X.c.i) + R1.i
  23:                                        ground.p.i:  0 = (((((-(ground.p.i) + RE.i) + CE.i) + -(Ci.i)) + R2.i) + -(X.c.i)) + -(R1.i)

... Simulate model Amplifier1
      Initialization at time = 0.0 s
      Initialization finished within  0.105574 seconds (70.49 k allocations: 3.675 MiB, 98.92% compilation time)
      Termination of Amplifier1 at time = 0.005 s
        nRejectedSteps            = 0
        nRejectedSteps            = 0
        nErrTestFails             = 0
        nTimeEvents               = 0
        nStateEvents              = 0

How can I step-through-debug the model, or failing that find out the exact variable dump of the error state?