second order differential equations & Laplace

[MrAl]

Hello again,

Oh ok, i think you meant this:

v1(t)=1-e^(-t)
v2(t)=(1-e^(-(t-a)))*(t>a)
v(t)=v1(t)-v2(t)

where a is the delay, and note that the factor (t>a) is a logical expression that evaluates to zero for t<=a and evaluates to 1 for t>a. That's just to make the math work for the unit step function u(t-a). The other part missing in your expression (but present previously) is the exponent t-a in v2(t).

Given the above expression for v(t), that's not periodic, but it does work for the first cycle of the square wave. With a period of 1 and 50 percent duty cycle then we have a=0.5
a=0.5
v1(t)=1-e^(-t)
v2(t)=(1-e^(-(t-a)))*(t>a)
v(t)=v1(t)-v2(t)

and now when we evaluate v(1) for example we get the right result, which is v=0.23865 to five significant figures. This agrees with other ways of calculating this.

But again that is not periodic. If we want it to be periodic, we have to do (for example) the following...

Define another interval:
v3(t)=(1-e^(-(t-2*a)))*(t>2*a)

and note all we did here was multiply 'a' by 2, and that brings us to 1.5 seconds because 2*a=1 and so this v3(t) starts to take effect after t=1 up to t=1.5 seconds. So using that next, we get:

a=0.5
v1(t)=1-e^(-t)
v2(t)=(1-e^(-(t-a)))*(t>a)
v3(t)=(1-e^(-(t-2*a)))*(t>2*a)
v(t)=v1(t)-v2(t)+v3(t)

and now we evaluate v(1.5) and we get:
v=0.53822 to five significant digits.

If we want the solution at t=2, we have to do that again. We get:
v4(t)=(1-e^(-(t-3*a)))*(t>3*a)

and note all we did here was multiply 'a' by 3 this time which brings us to 3*a+a or 2 seconds. So we now have in total:
a=0.5
v1(t)=1-e^(-t)
v2(t)=(1-e^(-(t-a)))*(t>a)
v3(t)=(1-e^(-(t-2*a)))*(t>2*a)
v4(t)=(1-e^(-(t-3*a)))*(t>3*a)
v(t)=v1(t)-v2(t)+v3(t)-v4(t)

and now we can evaluate v(2) and we get:
v(2)=0.32645

and that agrees with other methods too.

So you see to get the periodic nature in the result we have to apply that technique repeatedly. For each 'a' interval we have to use the next integer multiple of 'a' in the equation, and depending on if we are charging or discharging we either make the sign positive or negative.

If we do this enough times (roughly until more than 5 time constants) we start to see a repeat in the waveform which means we've reached the periodic portion of the wave. We can then look for the min and max to calculate the ripple.

An alternate method relies on calculating the initial value for each interval, then using the same equation over and over with the previous final value becoming the new initial value.


Here are the results of doing it the way explained above up to t=50 seconds. We can see that after 5 time constants we start to see repetition, and after ten times that (50 time constants) we see very close to the same result for every pulse. We can subtract the min from the max to get the ripple in this case, but note we cant always do that because the peaks dont always occur (for other circuits) at the transition period of the input wave (at the rise and fall times). We are lucky to have this here because it's only a first order equation.
Note also that if we take the average of the two last transition peaks, we get the average DC value. Again we are lucky we can do this here because it's only a first order equation. If it was second order we might have to hunt for the min and max a bit more carefully.

Another thing we can do is form this procedure into an equation as a series, and find the simplification of this series and it will result in a closed form equation for the N^th interval. It's not always easy to do though so it's up to you if you want to try this or not. It may be easy for this circuit but not for the second order circuit.

   t     V(t)
-----  --------
  0.5  0.393469
  1.0  0.238651
  1.5  0.538219
  2.0  0.326446
  2.5  0.591469
  3.0  0.358744
  3.5  0.611059
  4.0  0.370626
  4.5  0.618265
  5.0  0.374997
  5.5  0.620916
  6.0  0.376605
  6.5  0.621892
  7.0  0.377196
  7.5  0.622251
  8.0  0.377414
  8.5  0.622383
  9.0  0.377494
  9.5  0.622431
 10.0  0.377524
 10.5  0.622449
 11.0  0.377534
 11.5  0.622456
 12.0  0.377538
 12.5  0.622458
 13.0  0.377540
 13.5  0.622459
 14.0  0.377540
 14.5  0.622459
 15.0  0.377541
 15.5  0.622459
 16.0  0.377541
 16.5  0.622459
 17.0  0.377541
 17.5  0.622459
 18.0  0.377541
 18.5  0.622459
 19.0  0.377541
 19.5  0.622459
 20.0  0.377541
 20.5  0.622459
 21.0  0.377541
 21.5  0.622459
 22.0  0.377541
 22.5  0.622459
 23.0  0.377541
 23.5  0.622459
 24.0  0.377541
 24.5  0.622459
 25.0  0.377541
 25.5  0.622459

 

[derick007]

Hi MrAl,

I have used your suggestions/methods for the first order equation and entered values into the ss and plot on a graph and get a nice waveform of charging, discharging centred around 0.5v. - took an effort to type in all the values.

I have now tried using this on the second order equation and has taken even longer to type in the values. Amplitude = 23.5v, t1=6.43us, freq=43000hz, r= 10, L=100uH and C=22uF.

The ouput values just gradually increase. I have used intervals of 1us upto 100us !! I have ss to prove it - I would like to send it to you.

It looks as if there is a bit more to a second order than first order system ?

Kind Regards, Derek

 

[MrAl]

Hello again,


Yes there is a bit more to a second order system. It takes more care in working the equations.

 

[derick007]

Can I assume my transform for the first period only is correct and that the problem lies with the my interpretation of the equations for every period after the first ?

 

[MrAl]

y(t) = B * [ y1(t) - y2(t) ], were B = amplitude of periodic square wave input = 23.5.

y1(t) = { 1 - e^-t/2RC * (A sinwt + coswt)}

y2(t) = { 1 - e^-(t-t1)/2RC * (A sinw(t-t1) + cosw(t-t1)} * { u(t-t1) }

t1 = Duty Cycle * Period, A = df/sqrt(1-(df*df)) : df = damping factor : w = (1/sqrt(L*C)) - (1/(2*R*C))

Can I assume my transform for the first period only is correct and that the problem lies with the my interpretation of the equations for every period after the first ?

Hi,

What are you using for df ?

 

[derick007]

Hi,

What are you using for df ?

df = sqrt (L/C) / 2 * R, which in this case = 0.1066

 

[MrAl]

Hello again,


I think you meant:

df=sqrt(L/C)/(2*R)

right?

You might have to work on your typing out of expressions a little because this is important to be able to convey information correctly or you could end up with very wrong help replies.

You also should check your expression for 'w' because something is wrong there too.

I think once you get the expressions right you'll find this easy to do.

Basically, for the first cycle, you would use your first equation up to the analysis time point of the first cycle and then subtract your second equation up to the point where you want to analyze as long as the point is beyond where the high part of the cycle ends. If the high part did not end yet you dont use the second equation. You also can not go beyond the first cycle without more additional terms which we can get to later.

For example, say we have a 100us period and 50us high time (50 percent duty cycle) and we want to do the analysis up to 100us only. That includes one complete cycle but no more. We would do this:
y1=ft1(100us)
y2=ft2(100us)
then
ft=y1-y2

That's all there is to it.

But lets say we wanted to know the solution when t=75us, then we have:
y1=ft1(75us)
y2=ft2(75us)
ft=y1=y2

and that's all there is to it because we're just working within the first half cycle.

But remember because the pulse is only 50us wide if we look at times below 50us we dont use y2 because ft2 does not occur until after 50us. We also can not go above 100us yet.

So double check your expression for w as it does not look correct. Your expression for df does look correct.

Just to note a numerical example, same as above looking at 75us i get:
y1=0.9257
y2=0.1336
ft=y1-y2=0.7921

This is using our test set of 100uH, 22uF, and 10 ohms again. This also agrees very closely with a time simulation.

 

[derick007]

Apologies for the mistakes.
I have looked at w again and I think I have made a mistake trying to simplify it. I always had w=(1/sqrt(LC))*(sqrt(1-(L/(4*R*RC))), but I tried to simplify this and believe I have made a mistake in the simplification process. I have reverted to the expression above, hopefully this is correct.

I think I follow what you were saying about when to apply parts of the equations. In my case the square wave has an ON time of 6.43us and then OFF for 16.83us during the first period. Allow the pulse returns to zero after 6.43us, the capacitor (and therefore output) continues to charge via the energy stored in the inductor. However, as the energy stored in the inductor dissipates, the capacitor (and therefore the output voltage) voltage begins to discharge as the inductor cannot supply it any longer. I am not sure at what time the capacitor stops being charge by the inductor and starts to discharge. I assume that because waveforms are sinusoidal and symmetrical, that the charge/discharge point occurs mid way during the OFF time (16.83us) = 8.42us. This would mean the second part of the equation which deals with the capacitor discharging applies after 6.43 + 8.42 us = 14.85us ?

 

[MrAl]

Hello again,


You are so very close with the calculation of w now, the true result is:
w=sqrt(1/(C*L)-1/(4*C^2*R^2))

If you want to try again you should be able to obtain that result too.

The place where the cap starts to discharge may be a little tricky sometimes because it's not as easy to judge how much energy is in the inductor without doing a few calculations, so the point where the cap starts to discharge can be almost anywhere after 6.43us. It could be very close or very far.

We are a little lucky here in that the wave we have to look at to get this time value is the same as the output wave which is what we are solving for anyway. So with that in mind, we have again:

y1=ft1(t)
y2=ft2(t)
y=y1-y2

and since we are only allowing one pulse of 6.43us this last equation is valid for any time after 6.43us.

We also established that ft2(t) is just the same as ft1(t) but with t=t-a, so we have ft2(t)=ft1(t-a) so we can write this a little more straightforward as;

y=f(t)=ft(t)-ft(t-a)

so that makes it a little more evident from a quick inspection what we are doing to get the solution that's all. And because we have 'a' we also have the constraint:

t>=a

so we have:

y=f(t)=ft(t)-ft(t-a), {t>=a}

and this only works because we said we only have one pulse, not a whole train of pulses, and the positive pulse goes from 0 to t=a and then goes to zero.

So now we want to find where the cap starts to discharge. By discharge we mean the voltage starts to decrease from some positive value to some less positive value.

We can start by plotting f(t). We then see the cap charges (because of the pulse) and then starts to discharge at some point. So it goes up, then down, and there may be other points where it goes back up again but we dont need them because the first peak will be the highest. This can be proved if desired. In fact, there are sin and cos terms so the wave probably goes up and down quite a few times before settling out to near zero. But again we dont care because the first hump is the highest unless the system is unstable, and this system is not unstable.

So by plotting f(t) from t=6.43e-6 to say t=400e-6 we can see that the first hump occurs between 50us and 100us. We could narrow this down by plotting a second time from 50us to 100us, and then narrow that down again, and repeat with a more narrow time range, and keep doing that to the desired accuracy. Or we could solve this a bit more analytically.

Knowing that we do have a peak we can solve for the peak by setting the first derivative to zero:

d(f(t))/dt=0

And this means we have to take the time derivative of f(t) which is the subtraction of the two functions as before and set it equal to zero. We then solve for t.
This may not be that easy to do so we may have to do it numerically anyway, but if we substitute the numerical values for the components we can eliminate the exponentials and then solve using trigonometry and we get:

tpk=7.2283705334509e^-5

so the peak where the cap starts to discharge is very close to 72.3us.

You might try doing this yourself too so you can get a feel for this and verify these results.

 

[derick007]

I have checked my derivation for w and it now agrees with yours above.

I have also checked values for when the capacitor starts to discharge (by plugging in values into a ss and NOT differentiation) - again it agees with yours.

But because in my case T = 23.3us and the capacitor discharge time = 72.2us, does this mean the capacitor does not get an opportunity to discharge, with a periodic square input ? It would appear under these circumstances, the capacitor keeps charging past its desired operating point ?

 

[MrAl]

Hi there,


Well if the period is 23.3us and the cap does not start to discharge until we reach a time of 72.3us, then that means the cap gets another boost from the second pulse which will drive it up more than before. But this is only what happens for one pulse where we havent yet looked at what happens more distant into the future with more pulses. The cap voltage being driven up higher is just one part of the picture, the other part is that the increased voltage means the resistor draws more current and thus starts to discharge the cap at a sooner point in time.

If we look at what happens farther into the future with more pulse cycles, we see that as the voltage across the cap increases, the cap tends to start to discharge sooner because the voltage across the resistor is higher and that means more discharge current. So at some point in time with many more pulse cycles we eventually see a sort of dynamic equilibrium being reached where the charge from the next pulse is completely removed by the increased resistor current before the start of the next pulse. This always happens with passive elements because they always eat up at least some tiny amount of energy or more on each cycle.
Even if we increase the pulse width, we still see the same thing except the equilibrium point is reached at a higher output voltage. Higher output voltage means even more discharge current, so the cap discharges as much as it charges, to a high degree of perfection. The only way to get a higher voltage than that is to increase the pulse width even more, but of course we are limited to a 100 percent duty cycle unless we are using a ready made chip that might have an internal limit that is less than 100 percent (very often the case).

 

[derick007]

Hi MrAl,

I have been trying to solve when the capacitor stops charging and starts to discharge and have found a time of 7.00us.

Basically I have found this time by considering how much energy the inductor has stored between 0 - 6.43us and then how fast this energy is dissipated by the resistor and capacitor.

Initially I have only looked at the first period and assumed that the waveform f(t)=f1(t) only applies i.e. the capacitor continues to charge using this equation f1(t) for a time after t=6.43us.

At t = 6.43us I calculate the inductor current IL=1.5A and therefore the Energy stored = 1.12e^-4. At this point I assume the current in the inductor decreases, but the capacitor voltage still increases as it is still being supplied by the inductor. I assume the inductor current decreases at a rate according to the equation :

VL= L * dIL/dt and therefore the current dIL = VL * dt / L . I have calculated values at 6.44us, 6.50us and 7us.
dIL will decrease because VL = 23.5 - Vc and the voltage across the capacitor is increasing (Vc). Vc = f1(t).

This decrease in current means there is a decrease in the amount of energy stored in the inductor (0.5 * L * IL^2). As power is the rate of doing work, I can calculated the power taken out of the inductor by looking at the decrease in energy over a specified time period. Likewise for the resistor and capacitor I can calculate the increase in power dissipated as the voltage continues to rise. I find at t=7us, the amount of energy dissipated by the resistor and capacitor is greater than that stored in the inductor and conclude that the inductor can no longer sustain this increase in capacitor voltage therefore the rate of increase diminishes until the capacitor starts to discharge.

Am I any were near the solution ?

Kind Regards, Derek

 

[MrAl]

Hello again,

Since the inductor is:
v=L*di/dt

if you expect the current to decrease by 1.5 amps in 0.5us then you need a voltage equal to 300 volts. The capacitor voltage only changes a minuscule amount in 0.5us so that's not going to do it

If you need to solve this a different way, then look at the point where the inductor current minus the resistor current goes through zero. That's when the resistor starts to discharge the cap because the additional energy being supplied by the inductor only supplies enough energy for the resistor. None left to charge the cap anymore. You'll find this also occurs around 72us which is 10 times greater than 7us.

It's good that you are looking for an alternate method to figure this out, that shows that you're really thinking about this circuit rather than just throwing formulas at it. I am happy to see you are doing that. It's always good to have a secondary method to solve for something too so that you can double check your original results to make sure they are correct in the first place.

 

[derick007]

Thanks MrAl, for your encouragement and persistence.

Yes, I can see your thinking in that the capacitor changes from a charging mode (requiring current flow into it) to a discharging mode (were the capacitor supplies current). Therefore the current must change from +ve to -ve or vice-versa and therefore pass through zero, which in itself denotes a change of mode.

Many thanks and Kind regards

 

[derick007]

I have been thinking about the equation and proposed to use the following in my ss.

f(t) = f1(t) - f2(t-t1)*u(t-t1) + f3(t-T)*u(t-T) - f4(t-(T+t1))*u(t-(T+t1)) + f5(t-2T)*u(t-2T) -

f6(t-(2T+t1))*u(t-(2T+t1) + ....

I should clarify that all the functions f1, f2 etc. are really the same i.e. exponential, sine, cosine etc. with different delays as noted above built in.

I have already plugged in values for upto f6 and find the output rises from zero and at the present has levelled out at approx. 7.6v. Upto now the graph is a smooth curve without any oscillations.

Kind regards, Derek

 

[MrAl]

Hi again,

You could always check the average output and compare that too. The average is of course:

Vavg=Vin*TH/(TH+TL)

or

Vavg=Vin*TH/TP

where TH is the time high, TL is the time low, TP is the total period.

Note however that the actual wave may sometimes have large bumps and deep valleys so the average is somewhere in between those extremes.

 

[derick007]

I used the standard formula for the unit step response of a second order system and found the Percent Overshoot = 71.4% giving fmax = 11.15v. I also found the peak time to be 149us.

I will continue to enter values into my equation, however I believe the output voltage should not rise above 11.15 and the circuit should not oscillate before 149us ?

Kind Regards, Derek

 

[MrAl]

Hi,

Not sure what you mean by "circuit should not oscillate", do you mean that it wont go up and down just go up before that time?

Also, is that the unit step (times 23.5v) applied indefinitely?

 

[derick007]

I have input quite a few numbers and the output has started to oscillator, around zero ?

The peak overshoot = 11.14 as predicted and the time is around 130 us.

Not sure about the oscillations around zero ? The amplitude of the oscillations decrease exponentially and presumably will level out at zero and not 6.5v ?

I have recorded the first 4 +ve peak amplitudes as 11.14, 5.7, 2.9 and 1.48 ?

Kind regards, Derek

 

[MrAl]

Hello again,

These systems start out going up and down but progressively ratchet upward more and more until they start to taper off at some DC value like the input times the pulse width over the total period. If you dont see that then something is still wrong, unless of course the squareish wave input is at too low of a frequency and then it will go up and down for a long time without ramping up too much. The up and down behavior should definitely not be near zero after a long time period.