# How to solve a boost circuit (or any switched RLC circuit) in time domain?

I need to model a boost circuit with the following topology:

simulate this circuit – Schematic created using CircuitLab

I am interested in the resistor voltage, so I found the differential equation of each state and solved them by Euler's implicit method, alternating the equation to be solved as the state was changed.

I compared the results with a circuit simulator. Each differential equation is right alone, but when I put them together, one after the other, forming the boost switched dynamic, the result is wrong.

These are my steps:

'ON state' differential equation (just the capacitor discharge)

$$\frac{dV_{out}}{dt} = -\frac{V_{out}}{RC}$$

'OFF state' differential equation

$$\frac{d^2V_{out}}{dt^2} + \frac{1}{RC}\frac{dV_{out}}{dt} + \frac{1}{LC}V_{out} = \frac{V_{in}}{LC}$$

This is my (MATLAB) code:

R = 10;
L = 1e-5;
C = 1e-4;
Vin = 10;

h = 1e-6;   % time step
sizeVector = 10000; % number of steps
t = 0:h:sizeVector*h-h;  % time vector

Vout = zeros(1,sizeVector); % output voltage vector
dVout = zeros(1,sizeVector); % first derivative vector

freqSw = 5e3; % switching frequency
D = zeros(1,sizeVector); % State of switch for each simulation step (duty cycle = 0.5)
D = sign(cos(2*pi*freqSw.*t + pi/2 + 0.001));
D(D<0) = 0;

for k=2:sizeVector
if D(k) == 0
% Off state
Vout(k) = (Vout(k-1) + h*((dVout(k-1) + (h/(L*C))*Vin)/(1 + h/(R*C))))/(1 + (h^2/(L*C))/(1 + h/(R*C)));
dVout(k) = (dVout(k-1) + h*(Vin/(L*C) - Vout(k)/(L*C)))/(1 + h/(R*C));
end
if D(k) == 1
% On state
Vout(k) = Vout(k-1)/(1+h/(R*C));
dVout(k) = (Vout(k)-Vout(k-1))/h; % just for convenience
end

end

figure(1)
plot(t,Vout)
xlabel('Time (s)')
ylabel('Voltage (V)')
hold on


These are the comparations between PSIM results (red) and my code in MATLAB (blue). (For the capacitor discharge, I considered in both cases a initial voltage of 10V)

Open state:

Closed state:

Boost converter as in the code above:

I can never get steady state voltage superior than my input voltage (10V, in this case), but the waveform seems to consistent with each one of the states waveforms separately..

There must be an error in my code or am I missing some important point? Is this kind of solving possible?

Following the sugestions, I changed the second order differential equation into a system of two first order equations, in order to calculate the inductor's current, which is necessary because when the 'off state' begins, the inductors current must be a initial condition. In the previous model, I wasn't considering the current change on the 'on state'.

The 'off state' system:

$$\begin{bmatrix} \frac{dV_{out}}{dt}\\ \frac{dI_{ind}}{dt} \end{bmatrix} = \begin{bmatrix} -\frac{1}{RC} & \frac{1}{C}\\ -\frac{1}{L} & 0 \end{bmatrix} \times \begin{bmatrix} V_{out}\\ I_{ind} \end{bmatrix} + \begin{bmatrix} 0\\ \frac{Vin}{L} \end{bmatrix}$$

The 'on state' system: (it's actually not a system, but I keep it matricial for simplicity)

$$\begin{bmatrix} \frac{dV_{out}}{dt}\\ \frac{dI_{ind}}{dt} \end{bmatrix} = \begin{bmatrix} -\frac{1}{RC} & 0\\ 0 & 0 \end{bmatrix} \times \begin{bmatrix} V_{out}\\ I_{ind} \end{bmatrix} + \begin{bmatrix} 0\\ \frac{Vin}{L} \end{bmatrix}$$

There is still one problem. In this way, I can't prevent the inductor's current from being negative in the 'off state', so I inserted an if statement in every 'off state' iteration. If $I_{ind} < 0$, solve the system again imposing $I_{ind} = 0$. This can be made just setting to zero the elements related to $I_{ind}$ in the matrix.

Below are the result comparison (PSIM, first, MATLAB second) and the code.

R = 10;
L = 1e-3;
C = 1e-4;
Vin = 10; % DC input voltage

h = 1e-7;                % step
sizeVector = 120000;     % number of steps
t = 0:h:sizeVector*h-h;  % time vector

X = zeros(2,sizeVector); % state vector

A = [-1/(R*C) 1/C ; -1/L 0]; % off state matrix
b = [0; Vin/L];              % off state vector

Aaux = [-1/(R*C) 0 ; 0 0];   % off state auxiliar matrix
baux = [0 ; 0];              % off state auxiliar vector

A2 = [-1/(R*C) 0; 0 0];      % on state matrix
b2 = [0; Vin/L];             % on state vector

freqSw = 5e3;               % switching frequency
D = zeros(1,sizeVector); % State of switch for each simulation step (duty cycle = 0.5)
D = sign(cos(2*pi*freqSw.*t + pi/2 + 0.001));
D(D<0) = 0;

for k=2:sizeVector

if D(k) == 0 % off state
X(:,k) = (inv(eye(length(A)) - h*A))*(X(:,k-1) + h*b);
if X(2,k) < 0 % avoid Iind < 0
X(:,k) = (inv(eye(length(Aaux)) - h*Aaux))*(X(:,k-1) + h*baux);
end

end
if D(k) == 1 % on state
X(:,k) = (inv(eye(length(A2)) - h*A2))*(X(:,k-1) + h*b2);
end
end

figure(1)
plot(t,X(1,:))
hold on
plot(t,X(2,:))
xlabel('Time (s)')
ylabel('Voltage (V)\Current(A)')
legend('V_{out}','I_{in}')
ylim([-6 36])
grid
hold on

• Welcome to EE.SE. What is the value of 'L' ? Did you forget that with the switch closed L charges up as a magnetic field, and dumps its charge as a 'boost' when the switch is open. Normally a mosfet is the switch running at tens of KHZ, and L is 10uH or greater. 'C' filters out the ripple so the load sees a smooth DC voltage. The diode is a high-speed type. – Sparky256 Jan 16 '18 at 3:23
• You need to have an equation for what happens to the inductor current in the "on state". – The Photon Jan 16 '18 at 3:36
• In this example, L = 0.01 H. As I am just interessed in output voltage, I completely missed the inductor charging in the 'on state'. I got the point, now I just need to think how to take the inductor energy into account every time the 'off state' is reached. Thank you very much, Sparky256 and The Photon! – Felipe Dicler Jan 16 '18 at 3:50
• Yeah when you are running the "Off state" part of the simulation you need to start with the inductor current not zero but equal to where the "On state" left it. The challenge is how to express that IC (initial condition) to your (Matlab) simulation. – Brian Drummond Jan 16 '18 at 11:34

I can see at least one thing wrong with your equations (I haven't reviewed your code):

• There should be an equation describing the evolution of the inductor current during the "on state".

I would probably just keep a system of two first-order equations rather than try to write a 2nd order equation. It won't take Matlab any longer to solve it that way.

One way to analyze boost conversion is to compute the energy stored in the inductor based on Switch time, Ton, current, I inductance,L when the switch is closed then equate that to the energy dissipated in the load in a regulated loop.

P_in = 0.5 * L I²/T ( Energy/time when switch is ON)
Pout = 0.5 * I * Vin * Toff/T (Power Out)
Assuming lossless switching
output voltage $Vo=Vin\sqrt{\dfrac{kR_LT_{on}}{2L}}$ for 0< k< 1 with a regulated loop

• Perharps I didn't emphasized, but I need a time domain representation of voltage, including its transient. I am working on an algorithm to handle real time interstep delays which occur in switched circuit simulations. This model is for evaluating the algorithm and is supposed to run on an embedded environment, aiming real time simulation. – Felipe Dicler Jan 16 '18 at 11:25
• Then your motives are certainly unclear. Reinvent the wheel? or improve existing models and simulators? – Tony Stewart Sunnyskyguy EE75 Jan 16 '18 at 17:01
• My purpose is to study algorithms that handle switching delays in real time simulations. In a first approach, I will try to reproduce some papers' results regarding this subject, like this: link. The boost model will be used as a case study for evaluating the algorithms performance. – Felipe Dicler Jan 16 '18 at 20:08
• So you want to compare EMT, ATP, PSCAD, etc... – Tony Stewart Sunnyskyguy EE75 Jan 16 '18 at 20:38
• Perhaps what you need instead is PLECS software plexim.com/sites/default/files/flyers/flyer_smps_letter.pdf – Tony Stewart Sunnyskyguy EE75 Jan 16 '18 at 20:51