[MATLAB Grader] Comparing Transfer Functions

Question

BW 2022-6-29

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1750415-matlab-grader-comparing-transfer-functions

评论： Agustín 2023-9-15

Just wondering if there is a good way to compare transfer functions in MATLAB Grader?

I am testing a pseudo student response against the reference solution, purposely adding some slight deviation to check the test function.

Both the pseudo and reference solution output the same transfer function (at least to 4 dp), however from my testing, I believe the the condition:

abs(expected-actual)

cant be met, as it is the incorrect argument data type. Further to this,

(expected-actual)

creates an error function with higher order terms, I believe

minreal(expected-actual)

would make more sense.

Any help appreciated.

10 个评论
显示 8更早的评论隐藏 8更早的评论

Paul 2022-7-14

在 MATLAB Online 中打开

@Matt Rich

I do like the L-inf approach, but it all really depends on what @BW considers equivalent.

For example:

sys1 = zpk([],-1,1)
sys1 =
 
    1
  -----
  (s+1)
 
Continuous-time zero/pole/gain model.
sys2 = zpk(1,[-1 1],1)
sys2 =
 
     (s-1)
  -----------
  (s+1) (s-1)
 
Continuous-time zero/pole/gain model.

Of course the difference from an L-inf norm perspective is going to be very small

norm(sys1-sys2,Inf)
ans = 3.2149e-16

But from a grading perspective, whether or not sys1 is equivalent to sys2 depends on whether or not the problem requires the student to perform the pole zero cancellation.

Perhaps a more concrete example would be

sys1 = tf(1,[1 1]);
sys2 = tf([1 2],[1 3]);
sysfb1 = feedback(sys1,sys2)
sysfb1 =
 
      s + 3
  -------------
  s^2 + 5 s + 5
 
Continuous-time transfer function.
sysfb2 = sys1/(1 + sys1*sys2)
sysfb2 =
 
      s^2 + 4 s + 3
  ----------------------
  s^3 + 6 s^2 + 10 s + 5
 
Continuous-time transfer function.

Here, sysfb1 and sysfb2 are identical by the L-inf test

norm(sysfb1 - sysfb2,Inf)
ans = 0

But, if the point of exercise is learn how to use feedback() instead of doing transfer function algebra in Matlab, then perhaps sysfb2 is to be considered incorrect.

And then of course we have to make sure the tolerance used in the L-inf test is smaller than the gain of the systems in question. Are sys1 and sys2 suppsed to be equal.

sys1 = 1e-10*tf(1,[1 1]);
sys2 = 1.1e-10*tf(1,[1 1]);
norm(sys1-sys2,inf)
ans = 1.0000e-11

Finally, as to your excellent point about delays, I think it really extends even to systems with input/ouput delays, because even those can become internal delays after the subtraction

sys1 = tf(1,[1 1],'InputDelay',1);
sys2 = sys1;
norm(parallel(sys1,-sys2),Inf)  % works ok
ans = 0
sys2 = tf(1,[1 1],'InputDelay',2);  % doesn't work ok
norm(parallel(sys1,-sys2),Inf)
Error using DynamicSystem/norm
The "norm" command cannot be used for continuous-time models with internal delays. Use the "pade" command to approximate such delays.

Matt Rich 2022-7-19

编辑：Matt Rich 2022-7-19

在 MATLAB Online 中打开

Yes, @BW @Paul excellent point about systems with poles on the imaginary axis! In those cases, one potential modification to the

approach I can think of right now is to evaluate

for some nonzero, conveniently chosen U which cancels the poles of the error system

that are on the imaginary axis.

For example, form two "nearly identical" systems, each with 3 poles on the imaginary axis:

s = tf('s'); 
num = -7.0554*(s+2.919)*(s^2 + 0.04626*s + 0.3162);
den = s*(s+2.924)*(s+0.01636)*(s^2 + 0.5739*s + 6.944)*(s^2+1);
G1 = num/den
G1 =
 
                   -7.055 s^3 - 20.92 s^2 - 3.184 s - 6.512
  --------------------------------------------------------------------------
  s^7 + 3.514 s^6 + 9.679 s^5 + 23.96 s^4 + 9.011 s^3 + 20.45 s^2 + 0.3322 s
 
Continuous-time transfer function.
% duplicate the transfer function above introducing small errors 
err = 100*eps; 
num2 = (-7.0554 + err )*(s+2.919) *(s^2 + 0.04626*s + 0.3162);
den2 =   s *(s+2.924) *(s+0.01636 + err) *(s^2 + 0.5739*s + 6.944)*(s^2+1+err);
G2 = num2/den2
G2 =
 
                   -7.055 s^3 - 20.92 s^2 - 3.184 s - 6.512
  --------------------------------------------------------------------------
  s^7 + 3.514 s^6 + 9.679 s^5 + 23.96 s^4 + 9.011 s^3 + 20.45 s^2 + 0.3322 s
 
Continuous-time transfer function.

Check

(expected to be ∞ in this case since both systems have a pole on the imaginary axis)

norm(G1-G2,inf)
ans = Inf

Now, calculate the poles of each system, and extract those on the imaginary axis (allowing for some numerical error):

p1 = pole(G1);
p1imag = p1( abs(real(p1)) < 10*eps )
p1imag = 
   0.0000 + 0.0000i
   0.0000 + 1.0000i
   0.0000 - 1.0000i
p2 = pole(G2);
p2imag = p2( abs(real(p2)) < 10*eps )
p2imag = 
   0.0000 + 0.0000i
   0.0000 + 1.0000i
   0.0000 - 1.0000i

Take the union of the two sets of poles:

pImag = union(p1imag,p2imag)
pImag = 
0000 + 0.0000i
0000 - 1.0000i
0000 + 1.0000i
0000 - 1.0000i
0000 + 1.0000i

Use this set of poles to create a "convenient" (in this case such that

) nonzero U:

U = zpk(pImag,-1*ones(size(pImag)) ,1)
U =
 
  s (s^2 + 1)^2
  -------------
     (s+1)^5
 
Continuous-time zero/pole/gain model.
U = tf(U)
U =
 
  s^5 - 1.11e-15 s^4 + 2 s^3 - 1.11e-15 s^2 + s
  ---------------------------------------------
     s^5 + 5 s^4 + 10 s^3 + 10 s^2 + 5 s + 1
 
Continuous-time transfer function.

Finally, check

norm((G1-G2)*U,inf)
ans = 2.7102e-11

I have not tested this too thoroughly yet with different amounts of error, MIMO systems, and need to think if there's a nice way to make it work for discrete without adding much code for logical statements but it seems promising.

BW 2022-7-19

编辑：BW 2022-7-19

Thanks you once again,

Another option which worked in this case was the ratio of transfer functions, which should be close to unity (using minreal)

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Matt Rich 2022-7-14

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1750415-matlab-grader-comparing-transfer-functions#answer_1006900

编辑：Matt Rich 2022-7-14

The most robust approach to assess LTI systems in MATLAB Grader is probably to use the

norm of the difference between the reference and learner solution:

for some small tolerance.

% check existence

assert( exist('G', 'var') , "The submission must contain a variable named G.")

% check class is a type of LTI system object

assert( isa(G,'tf') | isa(G,'ss') | isa(G,'zpk') , "The variable G must be an LTI system object (tf, ss, or zpk).")

% check IO size

assert( isequal( size(G), size(referenceVariables.G) ) , "G has has the wrong input/output dimensions.")

% check equivalence -- is L_inf norm of difference (error) close enough to zero?

tol = 1e-10;

assert( norm( referenceVariables.G - G , 'inf' ) < tol , "G is not correct.")

If you want to enforce the system is represented using a specific class (e.g., transfer function) then edit the class check:

% check class is tf

assert( isa(G,'tf') , "The variable G must be an transfer function.")

2 个评论
显示无隐藏无

BW 2022-7-18

Excellent! Thanks all

Agustín 2023-9-15

Could we compare both numerator and denominador of the transfer function instead of evaluating the L_inf norm? This might help students to identify where the error is.

请先登录，再进行评论。

[MATLAB Grader] Comparing Transfer Functions

10 个评论
显示 8更早的评论隐藏 8更早的评论

采纳的回答

2 个评论
显示无隐藏无

更多回答（0 个）

另请参阅

类别

标签

产品

Community Treasure Hunt

[MATLAB Grader] Comparing Transfer Functions

10 个评论 显示 8更早的评论隐藏 8更早的评论

采纳的回答

2 个评论 显示 无隐藏 无

更多回答（0 个）

另请参阅

类别

标签

产品

Community Treasure Hunt

10 个评论
显示 8更早的评论隐藏 8更早的评论

2 个评论
显示无隐藏无