1
$\begingroup$

I have a circuit with a pulsed voltage source connected across a resistor and inductor that are in series:

circuit

Potential difference of V1: V1

The solution by LTSpice gives:

Current through the inductor: L1 current

Potential difference across the inductor: L1 voltage

Now I solve this problem for current and potential difference of L1 in Mathematica 11.3:

ClearAll["Global`*"]

(*inductance*) l1 = 0.001;
(*resistance*) r1 = 3;
(*pulsed voltage 1*) vUp = 24;
(*pulsed voltage 0*) vDown = 0;
(*voltage source definition*) 
v1[t_] := 
  Piecewise[{{vUp, 0.001 <= t <= 0.011}, {vUp, 0.013 <= t <= 0.023}}, 
   vDown];
(*ImageSize in plots*) imgSize = 350;
(*Model time*) time = 0.025;

eq1 = v1[t] == i[t]*r1 + l1*i'[t];
ic1 = i[0] == 0;

sol = NDSolveValue[{eq1, ic1}, i, {t, 0, time}];

p1 = Plot[v1[t], {t, 0, time},
   PlotRange -> All,
   PlotPoints -> 200,
   AxesOrigin -> {0, 0},
   Frame -> True,
   GridLines -> Automatic,
   GridLinesStyle -> LightGray,
   PlotLabel -> "Applied potential difference",
   PlotStyle -> Thick,
   Exclusions -> 
    None (*for connection of piecewise function v1 in step up/down*),
   ImageSize -> imgSize
   ];

p2 = Plot[sol[t], {t, 0, time},
   PlotRange -> All,
   PlotPoints -> 200,
   AxesOrigin -> {0, 0},
   Frame -> True,
   GridLines -> Automatic,
   GridLinesStyle -> LightGray,
   PlotLabel -> "Current through resistor r1 and inductor l1",
   PlotStyle -> Thick,
   ImageSize -> imgSize
   ];

p3 = Plot[sol'[t]*l1, {t, 0, time},
   PlotRange -> All,
   PlotPoints -> 300,
   AxesOrigin -> {0, 0},
   Frame -> True,
   GridLines -> Automatic,
   GridLinesStyle -> LightGray,
   PlotLabel -> "Potential difference across inductor l1",
   PlotStyle -> Thick,
   Exclusions -> None,
   ImageSize -> imgSize
   ];

Grid[{{p2, p1}, {p3}}, Frame -> True]

The solution by Mathematica corresponds to what I have in LTSpice (so I can assume that my model functions well): MathematicaSolution

Now let's switch resistance of the resistor R1 to 1.0 Ohm

(*resistance*) r1 = 1;

This results in current ramping up to 24 amps: 24amps

However, my real power supply can provide only up to 10 amps. How to do I properly setup a current limit in this model?

$\endgroup$
1
  • 1
    $\begingroup$ This strikes me as more of a modeling/mathematical problem than a Mathematica one? $\endgroup$
    – march
    Commented Nov 19, 2018 at 9:23

2 Answers 2

3
$\begingroup$

Another approach

The simulation starts normally with the equation : v1[t] == i[t] r1 + l1 i'[t].

But when the current reaches the limit imax=10A, the simulation is stopped and restarted with a new equation (and a new initial condition). The new equation corresponds to i[t] = constant. It is simply : i'[t] == 0. The equation v1[t] == i[t] r1 + l1 i'[t] which is not valid anymore disappears.

Of course, this idea is only interesting if the system returns to the first equation when necessary, like the real system. It turns out that this is possible with the same kind of switching.

(*inductance*) l1 = 0.001;
(*resistance*) r1 = 1;
(*pulsed voltage 1*) vUp = 24;
(*pulsed voltage 0*) vDown = 0;
(*voltage source definition*) 
v1[t_] := 
    Piecewise[{{vUp, 0.001 <= t <= 0.011}, {vUp, 
     0.013 <= t <= 0.023}}, 
      vDown];
(*ImageSize in plots*) imgSize = 350;
(*Model time*) time = 0.025;


myNDSolveOption = 
  "ExtrapolationHandler" -> {Indeterminate &, 
    "WarningMessage" -> False};
imax = 10; (* maximum current *)

processSimulation[equsAndIcs_, resultList_, tfin_] := 
 Module[{ndssdata},
  ndssdata = 
   First @ NDSolve`ProcessEquations[equsAndIcs, {i[t]}, {t, 0, time}, 
     "MaxStepSize" -> 0.001, 
     "ExtrapolationHandler" -> {Indeterminate &, 
       "WarningMessage" -> False}];
  NDSolve`Iterate[ndssdata, tfin];
  {NDSolve`SolutionDataComponent[
    ndssdata @ "SolutionData" ["Forward"], "Time"],
   Cases[Equal @@@ NDSolve`ProcessSolutions[ndssdata, "Forward"], 
     i[_] == _][[1]],
   Append[resultList, i[t] /. NDSolve`ProcessSolutions[ndssdata]]}
  ]


res00 = {0, i[0] == 0, {}};
While[res00[[1]] < time,
 res00 = processSimulation[{v1[t] == i[t] r1 + l1  i'[t], res00[[2]], 
    WhenEvent[i[t] > imax, "StopIntegration" ]}, res00[[3]], time];
 res00 = processSimulation[{ i'[t] == 0, res00[[2]], 
    WhenEvent[i[t] r1 + l1 i'[t] > v1[t], "StopIntegration"]}, 
   res00[[3]], time]
 ]

Plot[Evaluate[res00[[3]]], {t, 0, time}, 
 PlotStyle -> (Directive[#, Thickness[0.01]] & /@ {Black, Red})]
Plot[Evaluate[D[#, t] & /@ (res00[[3]])], {t, 0, time}, 
 PlotStyle -> (Directive[#, Thickness[0.01]] & /@ {Black, Red})]

enter image description here

enter image description here

The black segments correspond to the first equation, the red ones to the second equation.

$\endgroup$
1
  • $\begingroup$ Andre, thank you very much for putting your time into this question even after I accepted your previous answer. The new one is more advanced because it accurately represents unrestrained inductive kick-back. $\endgroup$ Commented Dec 3, 2018 at 13:48
2
$\begingroup$

The results depend of the way the voltage source limits the current (slowly ? quickly ? is there a decoupling capacitor ?).

The following solution supposes that the source switches brutally from to 0V when I>10A and returns brutally to the normal value when I<9.9A.

The relevant code is :

  • eq1 = k[t] v1[t] == (i[t]*r1 + l1*i'[t]); (* k[t] = 0 or 1 *)
  • NDSolve[... eq1 ... WhenEvent[i[t]>10,k[t]->0],WhenEvent[i[t]<9.9,k[t]->1]]

The whole code :

  ClearAll["Global`*"]

(*inductance*) l1 = 0.001;
(*resistance*) r1 = 1;
(*pulsed voltage 1*) vUp = 24;
(*pulsed voltage 0*) vDown = 0;
(*voltage source definition*) 
v1[t_] := 
  Piecewise[{{vUp, 0.001 <= t <= 0.011}, {vUp, 0.013 <= t <= 0.023}}, 
   vDown];
(*ImageSize in plots*) imgSize = 350;
(*Model time*) time = 0.025;

eq1 = k[t]  v1[t] ==  (i[t]*r1 + l1*i'[t]); (* k[t] = 0 or 1 *)
eq2 = vfil[t] + 5. 10^-5 vfil'[t] == l1*i'[t]; (* just a filter to see the 
filtered voltage across the self *)

ic1 = i[0] == 0;

sol = NDSolveValue[{eq1, ic1, eq2, vfil[0]==0,k'[t]==1,k[0]==1,WhenEvent[i[t]>10,k[t]->  0 ],WhenEvent[i[t]<9.9,k[t]->  1 ]}, {i,vfil}, {t, 0, time}]

p1 = Plot[v1[t], {t, 0, time},
   PlotRange -> All,
   PlotPoints -> 200,
   AxesOrigin -> {0, 0},
   Frame -> True,
   GridLines -> Automatic,
   GridLinesStyle -> LightGray,
   PlotLabel -> "Applied potential difference",
   PlotStyle -> Thick,
   Exclusions -> 
    None (*for connection of piecewise function v1 in step up/down*),
   ImageSize -> imgSize
   ];

p2 = Plot[sol[[1]][t], {t, 0, time},
   PlotRange -> All,
   PlotPoints -> 200,
   AxesOrigin -> {0, 0},
   Frame -> True,
   GridLines -> Automatic,
   GridLinesStyle -> LightGray,
   PlotLabel -> "Current through resistor r1 and inductor l1",
   PlotStyle -> Thick,
   ImageSize -> imgSize
   ];

p3 = Plot[sol[[1]]'[t]*l1, {t, 0, time},
   PlotRange -> All,
   PlotPoints -> 300,
   AxesOrigin -> {0, 0},
   Frame -> True,
   GridLines -> Automatic,
   GridLinesStyle -> LightGray,
   PlotLabel -> "Potential difference across inductor l1",
   PlotStyle -> Thick,
   Exclusions -> None,
   ImageSize -> imgSize
   ];

 p4 = Plot[sol[[2]][t] , {t, 0, time},
   PlotRange -> All,
   PlotPoints -> 300,
   AxesOrigin -> {0, 0},
   Frame -> True,
   GridLines -> Automatic,
   GridLinesStyle -> LightGray,
   PlotLabel -> "Potential difference across inductor l1 filtered",
   PlotStyle -> Thick,
   Exclusions -> None,
   ImageSize -> imgSize
   ];

Grid[{{p2, p1}, {p3,p4}}, Frame -> True]

enter image description here

Note : A voltage source that switches to 0 Volts is not the same as a voltage source that switches to 0 Amperes (ie I=0, ie open-circuit).

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.