The Problem
In its ideal form, the impulse should have an infinite amplitude at t=0
and the area under the curve should be 1:
In other words, the Dirac delta function is a theoretical concept that cannot really be simulated.
Impulse Response
Before seeing how to implement a Dirac delta function in Simulink, let's talk a bit about why would someone need that.
If you are looking into implementing an impulse or Dirac delta function, it's very likely that you want to compute the impulse response of a system. Computing the impulse response is a common way to characterize the dynamic response of a system. This is used in many areas like vibrations and acoustics. As an example, let’s define a system with two poles, one at 5 rad/sec and one at 50 rad/sec.
sys = tf(12,[1 0.5 25]) + tf(100,[1 2 2500]);
[y,t] = impulse(sys,tFinal);
Let's now look at how something similar could be implemented in Simulink.
Fixed-Step Solver
The implementation of an impulse in Simulink will probably be different if implemented fixed versus variable step. Let's begin with the fixed-step case.
With a discrete signal, the Hit Crossing block outputs a value of 1 for one time step when the specified hit crossing offset is crossed. To respect the Dirac definition that the area under the curve must be equal to 1, we can divide this signal by the sample time. mdl = "impulse_ex_discrete";
In this example, I want to trigger the impulse to happen when the output of the Digital Clock block crosses a specific value, let's say 1 second. Using the same transfer function as in our previous example, we can see that the results match.
plot(t+1,y,'LineWidth',2) % Offset the time reference by 1 second
legend('MATLAB','Simulink')
Variable-Step Solvers
With a variable-step solver, there are a few subtle things that you need to take care of.
Impulse at t=0
If what you need is an impulse at t=0, I recommend using the Step block. This will allow you to control the amplitude and duration of the pulse. The important thing to decide is the value to use for the duration and amplitude of the impulse. Since this is impossible to generate an infinitely large amplitude with an infinitely small width, we need to pick a number.
Let's simulate the model with different values for the amplitude and width of the impulse. To respect the assumption that the integral of the pulse is 1, I make the amplitude the inverse of the duration.
in = repmat(Simulink.SimulationInput(mdl),length(p),1);
in(i) = in(i).setVariable('impulseDuration',p(i));
in(i) = in(i).setVariable('impulseAmplitude',1/p(i));
out = sim(in,'ShowProgress','off');
The key thing to realize is that the longer the duration of the impulse, the more it becomes a step instead of an impulse. In the Laplace world, an ideal Diract delta function has an amplitude of 1 for all frequencies, while a step has an amplitude of . The longer the impulse, the less it excites high frequencies.
Looking at the results, we can see that the second mode at 50 rad/sec is significantly less excited when the duration of the step becomes longer.
plot(out(i).yout{1}.Values)
legend(['Ideal' arrayfun(@num2str, p, 'UniformOutput', false)])
Impulse at times other than zero
While the definition of the Dirac delta function says it should return a non-zero value only at t=0, it can be useful to generate an impulse at other times.
As a first attempt, it would be tempting to use a pattern similar to the one used for fixed-step. The challenge here is that with a variable-step solver, we don't have direct control on the duration of the pulse generated by the Hit Crossing block.
Based on what I know from zero-crossing detection, I guessed a gain value of 1e14 for the amplitude of the pulse, but we can see that it's not quite right. in = Simulink.SimulationInput(mdl);
plot(t+1.23456,y,'LineWidth',2)
Scheduling the Impulse Duration
To control the duration of the impulse, it is possible to use the Hit Scheduler block. When the time to trigger the impulse is detected by the Hit Crossing block we do two things: First, the event TrigOn is sent to the Stateflow chart to initiate the impulse. Then the Hit Scheduler schedules a second event TrigOff to be sent to the chart after a specified amount of time to terminate the impulse.
If you prefer, you can slightly simplify the logic by using the Discrete-Event Chart from SimEvents. In that case, you would configure the Hit Crossing block to output a message instead of a signal. This message will wake up the chart, and the chart will keep the impulse on for a duration specified using after. plot(t+1.23456,y,'LineWidth',2)
plot(out1.yout{1}.Values)
plot(out1.yout{1}.Values)
legend("MATLAB", "Hit Scheduler", "Discrete-Event Chart");
Approximation using a normal distribution
To implement that, I can use the same logic as before to initiate and terminate the pulse. Instead of just outputting an on/off signal, I can use an Enabled Subsystem to "play" the normal distribution function.
In this case, we can see that the results match what we got from the impulse function.
plot(out.logsout.get('u').Values)
axis([1.236 1.244 -10 600]);
plot(out.logsout.get('y').Values)
hold on% compare to reference
Now it's your turn
Have fun implementing impulses or Dirac delta functions in your Simulink model. If you have other techniques, share in the comments below.
댓글
댓글을 남기려면 링크 를 클릭하여 MathWorks 계정에 로그인하거나 계정을 새로 만드십시오.