bugGNU Octave - Bugs: bug #58538, [octave forge] (control) step and...

 
 

bug #58538: [octave forge] (control) step and lsim fail on high order transfer functions

Submitted by:  None
Submitted on:  Wed 10 Jun 2020 05:16:56 PM UTC  
 
Category:  Octave Forge Package Severity:  3 - Normal
Priority:  3 - Low Item Group:  Incorrect Result
Status:  In Progress Assigned to:  None
Originator Name:  gs Originator Email:  -email is unavailable-
Open/Closed:  Open Release:  dev
Operating System:  GNU/Linux

Add a New Comment (Rich Markup)
   

You are not logged in

Please log in, so followups can be emailed to you.

 

( Jump to the original submission)

Fri 19 Jun 2020 07:49:27 PM UTC, comment #8: 

Changeset 1618:9bfa416611fd solves the issue for lsim ().

Torsten Lilge <ttl>
Project Member
Fri 19 Jun 2020 07:19:10 PM UTC, comment #7: 

I have pushed changeset 1617:104e9c0fc974 which solves the problem for step () by converting a continuous-time system into a state space representation before using c2d for discretization.

The issue is still present in lsim() and impulse(). In the latter case, the discretization is done using the "impulse" option of c2d, which obviously leads to an inaccurate discrete-time model even for a state space representation.

Torsten Lilge <ttl>
Project Member
Wed 10 Jun 2020 08:52:04 PM UTC, comment #6: 

Nicholas, thanks for posting.
Just to avoid confusions: I changed the order of the two cases before posting here. So the results posted by in comment #5 are in different order compared to my comment #2.

Torsten Lilge <ttl>
Project Member
Wed 10 Jun 2020 08:31:20 PM UTC, comment #5: 

(adding details from the mailing list)
running step and lsim commands in comment #0  result in the following warning being displayed repeatedly (for step, it appears every time the figure is given focus until a ctrl-c)

warning: opengl_renderer: data values greater than float capacity.  (1) Scale data, or (2) Use gnuplot

switching to gnuplot results in the output attached as file OctaveControlGnuplot.png

looking at the output from Torsten's tests:

output on Octave:

>> pkg load control
>> q=tf([  -1.46668928632582  38.5234944088261  -1023.42433722006  18055.8729050165  -275471.128598611  3680685.32154243  -40413629.6869944  433109318.542883  -3699840188.09672  32827566143.7738  -225798332168.488  1689921426466.71  -9494747382069.7  60613579740811.5  -278626623858487  1.52747193482069e+15  -5.6883783399098e+15  2.69067309068308e+16  -7.92226372994905e+16  3.25360727881732e+17  -7.2059792475876e+17  2.61175754638012e+18  -3.8909506684981e+18  1.31637705830625e+19  -9.12640045084024e+18  3.78290213329238e+19  1.41630594163597e+19  5.02932079881115e+19  1.29547036349465e+20  1.05983336007477e+19  2.43633391646665e+20  -1.31208497581345e+19  1.04713601706961e+20  7.83465823853832e+18  ],[  1  17.1698120499521  542.934121038049  7294.98595152311  126155.971925398  1377306.37948524  16829405.6386297  152789570.438515  1448656865.6342  11095921096.5808  85365123037.7488  556372911296.259  3555399791580.02  19799716986834.5  106353561602035  506292672690482  2.29624277698132e+15  9.31618125768562e+15  3.56349262092232e+16  1.22398996473335e+17  3.92373264437545e+17  1.12860046096021e+18  2.99632129507514e+18  7.09934073973815e+18  1.53176948230893e+19  2.91760427412351e+19  4.97010301791734e+19  7.33491370510952e+19  9.42777562297368e+19  1.01604047339533e+20  9.12510568398965e+19  6.44576522525933e+19  3.46367980244933e+19  1.23309424846899e+19  2.4145011614912e+18  ]);
>> sys_ss_ct = ss (q);
>> sys_dt2   = c2d (sys_ss_ct, 0.01);
>> eig (sys_dt2.a)
ans =

   0.99422 + 0.08037i
   0.99422 - 0.08037i
   0.99346 + 0.07615i
   0.99346 - 0.07615i
   0.99333 + 0.07155i
   0.99333 - 0.07155i
   0.99368 + 0.06688i
   0.99368 - 0.06688i
   0.99488 + 0.00426i
   0.99488 - 0.00426i
   0.99458 + 0.00842i
   0.99458 - 0.00842i
   0.99435 + 0.01262i
   0.99435 - 0.01262i
   0.99443 + 0.01718i
   0.99443 - 0.01718i
   0.99379 + 0.02253i
   0.99379 - 0.02253i
   0.99340 + 0.02768i
   0.99340 - 0.02768i
   0.99335 + 0.03249i
   0.99335 - 0.03249i
   0.99329 + 0.03711i
   0.99329 - 0.03711i
   0.99324 + 0.04174i
   0.99324 - 0.04174i
   0.99327 + 0.04646i
   0.99327 - 0.04646i
   0.99345 + 0.05142i
   0.99345 - 0.05142i
   0.99368 + 0.05637i
   0.99368 - 0.05637i
   0.99406 + 0.06152i
   0.99406 - 0.06152i

>> sys_tf_dt = c2d (q, 0.01);
>> sys_dt1   = ss (sys_tf_dt);
>> eig (sys_dt1.a)
ans =

   1.42735 + 0.92459i
   1.42735 - 0.92459i
   1.17401 + 0.91528i
   1.17401 - 0.91528i
   0.97044 + 0.84467i
   0.97044 - 0.84467i
   0.81966 + 0.74285i
   0.81966 - 0.74285i
   0.70737 + 0.63075i
   0.70737 - 0.63075i
   0.61757 + 0.53082i
   0.61757 - 0.53082i
   0.55208 + 0.43837i
   0.55208 - 0.43837i
   0.50582 + 0.34596i
   0.50582 - 0.34596i
   0.47599 + 0.25501i
   0.47599 - 0.25501i
   0.45893 + 0.16733i
   0.45893 - 0.16733i
   0.45075 + 0.08282i
   0.45075 - 0.08282i
   0.44843 + 0.00000i
   0.79468 + 0.52856i
   0.79468 - 0.52856i
   0.64503 + 0.00000i
   0.93407 + 0.00000i
   1.70651 + 0.84064i
   1.70651 - 0.84064i
   1.96533 + 0.64519i
   1.96533 - 0.64519i
   2.14777 + 0.35049i
   2.14777 - 0.35049i
   2.21291 + 0.00000i

>> diary off

output on Matlab 2020a:

q=tf([  -1.46668928632582  38.5234944088261  -1023.42433722006  18055.8729050165  -275471.128598611  3680685.32154243  -40413629.6869944  433109318.542883  -3699840188.09672  32827566143.7738  -225798332168.488  1689921426466.71  -9494747382069.7  60613579740811.5  -278626623858487  1.52747193482069e+15  -5.6883783399098e+15  2.69067309068308e+16  -7.92226372994905e+16  3.25360727881732e+17  -7.2059792475876e+17  2.61175754638012e+18  -3.8909506684981e+18  1.31637705830625e+19  -9.12640045084024e+18  3.78290213329238e+19  1.41630594163597e+19  5.02932079881115e+19  1.29547036349465e+20  1.05983336007477e+19  2.43633391646665e+20  -1.31208497581345e+19  1.04713601706961e+20  7.83465823853832e+18  ],[  1  17.1698120499521  542.934121038049  7294.98595152311  126155.971925398  1377306.37948524  16829405.6386297  152789570.438515  1448656865.6342  11095921096.5808  85365123037.7488  556372911296.259  3555399791580.02  19799716986834.5  106353561602035  506292672690482  2.29624277698132e+15  9.31618125768562e+15  3.56349262092232e+16  1.22398996473335e+17  3.92373264437545e+17  1.12860046096021e+18  2.99632129507514e+18  7.09934073973815e+18  1.53176948230893e+19  2.91760427412351e+19  4.97010301791734e+19  7.33491370510952e+19  9.42777562297368e+19  1.01604047339533e+20  9.12510568398965e+19  6.44576522525933e+19  3.46367980244933e+19  1.23309424846899e+19  2.4145011614912e+18  ]);
sys_ss_ct = ss (q);
sys_dt2   = c2d (sys_ss_ct, 0.01);
eig (sys_dt2.a)

ans =

   0.9942 + 0.0804i
   0.9942 - 0.0804i
   0.9935 + 0.0762i
   0.9935 - 0.0762i
   0.9933 + 0.0716i
   0.9933 - 0.0716i
   0.9937 + 0.0669i
   0.9937 - 0.0669i
   0.9941 + 0.0615i
   0.9941 - 0.0615i
   0.9937 + 0.0564i
   0.9937 - 0.0564i
   0.9935 + 0.0514i
   0.9935 - 0.0514i
   0.9933 + 0.0465i
   0.9933 - 0.0465i
   0.9932 + 0.0417i
   0.9932 - 0.0417i
   0.9933 + 0.0371i
   0.9933 - 0.0371i
   0.9934 + 0.0325i
   0.9934 - 0.0325i
   0.9934 + 0.0277i
   0.9934 - 0.0277i
   0.9938 + 0.0225i
   0.9938 - 0.0225i
   0.9944 + 0.0172i
   0.9944 - 0.0172i
   0.9943 + 0.0126i
   0.9943 - 0.0126i
   0.9946 + 0.0084i
   0.9946 - 0.0084i
   0.9949 + 0.0043i
   0.9949 - 0.0043i

sys_tf_dt = c2d (q, 0.01);
sys_dt1   = ss (sys_tf_dt);
eig (sys_dt1.a)

ans =

   2.0691 + 0.0000i
   2.0125 + 0.2961i
   2.0125 - 0.2961i
   1.8581 + 0.5393i
   1.8581 - 0.5393i
   1.6475 + 0.6911i
   1.6475 - 0.6911i
   1.4552 + 0.7483i
   1.4552 - 0.7483i
   1.2923 + 0.8011i
   1.2923 - 0.8011i
   1.1007 + 0.8116i
   1.1007 - 0.8116i
   0.9311 + 0.7639i
   0.9311 - 0.7639i
   0.8008 + 0.6860i
   0.8008 - 0.6860i
   0.7056 + 0.5990i
   0.7056 - 0.5990i
   0.6352 + 0.5112i
   0.6352 - 0.5112i
   0.5848 + 0.4239i
   0.5848 - 0.4239i
   0.5529 + 0.3406i
   0.5529 - 0.3406i
   0.5346 + 0.2680i
   0.5346 - 0.2680i
   0.5163 + 0.2043i
   0.5163 - 0.2043i
   0.4994 + 0.1373i
   0.4994 - 0.1373i
   0.4897 + 0.0687i
   0.4897 - 0.0687i
   0.4866 + 0.0000i

diary off


Nicholas Jankowski <nrjank>
Wed 10 Jun 2020 08:12:55 PM UTC, comment #4: 

comment #1:

> So, to clarify because it's not entirely clear from your bug description, the bug you are reporting is that 'step(sys)' tries to plot too many points and overloads the Octave OpenGL plotting system, while 'step(ss(sys))' works as a workaround and does plot the correct figure. Is that accurate?


step() and lsim() create simulations that are completely wrong, the simulations seem to be unstable which create very large values and the plot window errors.

so 2 problems:
1: simulation is incorrect, gives incorrect values, not even close
2. plot window can't handle very large (probably inf) values

Anonymous
Wed 10 Jun 2020 07:50:44 PM UTC, comment #3: 

Changing from "Matlab Compatibility" to "Incorrect Result" and raising severity, since this is not just a problem in plotting the result.

Torsten Lilge <ttl>
Project Member
Wed 10 Jun 2020 07:47:25 PM UTC, comment #2: 

No, there is a problem with the internal system representation used by step (). As posted in the help list, this seems to be related to a problem when discretizing a continuous-time transfer function. To test:

# Get q from the original post

# discretization of the transfer function
# unstable eigenvalues (|lambda| > 1)
sys_tf_dt = c2d (q, 0.01);
sys_dt1   = ss (sys_tf_dt);
lambda = eig (sys_dt1.a)
stable = isstable (sys_dt1)

# discretization of the state space model
# => stable eigenvalues
sys_ss_ct = ss (q);
sys_dt2   = c2d (sys_ss_ct, 0.01);
lambda = eig (sys_dt2.a)
stable = isstable (sys_dt2)

I have attached a script that reads the transfer function q and then does the test given above.

(file #49237)

Torsten Lilge <ttl>
Project Member
Wed 10 Jun 2020 07:25:29 PM UTC, comment #1: 

So, to clarify because it's not entirely clear from your bug description, the bug you are reporting is that 'step(sys)' tries to plot too many points and overloads the Octave OpenGL plotting system, while 'step(ss(sys))' works as a workaround and does plot the correct figure. Is that accurate?

Mike Miller <mtmiller>
Project Administrator
Wed 10 Jun 2020 05:16:56 PM UTC, original submission:  

%Octave 5.2.1 and Control 3.2.0

pkg load control
q=tf([  -1.46668928632582  38.5234944088261  -1023.42433722006  18055.8729050165  -275471.128598611  3680685.32154243  -40413629.6869944  433109318.542883  -3699840188.09672  32827566143.7738  -225798332168.488  1689921426466.71  -9494747382069.7  60613579740811.5  -278626623858487  1.52747193482069e+15  -5.6883783399098e+15  2.69067309068308e+16  -7.92226372994905e+16  3.25360727881732e+17  -7.2059792475876e+17  2.61175754638012e+18  -3.8909506684981e+18  1.31637705830625e+19  -9.12640045084024e+18  3.78290213329238e+19  1.41630594163597e+19  5.02932079881115e+19  1.29547036349465e+20  1.05983336007477e+19  2.43633391646665e+20  -1.31208497581345e+19  1.04713601706961e+20  7.83465823853832e+18  ],[  1  17.1698120499521  542.934121038049  7294.98595152311  126155.971925398  1377306.37948524  16829405.6386297  152789570.438515  1448656865.6342  11095921096.5808  85365123037.7488  556372911296.259  3555399791580.02  19799716986834.5  106353561602035  506292672690482  2.29624277698132e+15  9.31618125768562e+15  3.56349262092232e+16  1.22398996473335e+17  3.92373264437545e+17  1.12860046096021e+18  2.99632129507514e+18  7.09934073973815e+18  1.53176948230893e+19  2.91760427412351e+19  4.97010301791734e+19  7.33491370510952e+19  9.42777562297368e+19  1.01604047339533e+20  9.12510568398965e+19  6.44576522525933e+19  3.46367980244933e+19  1.23309424846899e+19  2.4145011614912e+18  ]);

step(q)
%lsim(q,ones(1,1001),[0:.01:10])  %another simulator, same result

*****************************************

these workarounds work which indicate that the internal tf2ss step needs correction:

1:
pkg load signal
[a,b,c,d] = tf2ss(q);
step(ss(a,b,c,d))

2:
only with the control pkg
step(ss(q))

**************************************

picture show result from Matlab 2020a, matlab can handle the tf's. Was only possible to upload one picture.

see this thread for mote info:
https://octave.1599824.n4.nabble.com/control-package-can-t-simulate-complex-systems-td4697481.html

Anonymous

 

(Note: upload size limit is set to 16384 kB, after insertion of the required escape characters.)

Attach Files:
   
   
Comment:
   

Attached Files
file #49238:  OctaveControlGnuplot.png added by nrjank (110KiB - image/png)
file #49237:  bug58538_test.m added by ttl (2KiB - text/x-objcsrc)
file #49236:  MatlabControlOutput.png added by None (103KiB - image/png)

 

Depends on the following items: None found

Items that depend on this one: None found

 

Carbon-Copy List
  • -email is unavailable- added by nrjank (Updated the item)
  • -email is unavailable- added by ttl (Updated the item)
  •  

    Do you think this task is very important?
    If so, you can add your encouragement to it.
    This task has 0 encouragements so far.

    Only project members can vote.

     

     

     

    Follow 12 latest changes.

    Date Changed by Updated Field Previous Value => Replaced by
    2020-06-19 ttl StatusConfirmed => In Progress
    2020-06-10 nrjank Attached File- => Added OctaveControlGnuplot.png, #49238
    2020-06-10 ttl Severity2 - Minor => 3 - Normal
        Item GroupMatlab Compatibility => Incorrect Result
    2020-06-10 ttl Attached File- => Added bug58538_test.m, #49237
    2020-06-10 mtmiller CategoryNone => Octave Forge Package
        Severity3 - Normal => 2 - Minor
        Priority5 - Normal => 3 - Low
        StatusNone => Confirmed
        Release5.2.0 => dev
        Summarycontrol package, step() and lsim() can't simulate complex transfer functions => [octave forge] (control) step and lsim fail on high order transfer functions
    2020-06-10 None Attached File- => Added MatlabControlOutput.png, #49236

    Back to the top


    Powered by Savane 3.5