bugGNU Octave - Bugs: bug #60384, `mpoles` incorrectly classifies...

 
 

bug #60384: `mpoles` incorrectly classifies distinct poles as unique ones

Submitter:  Luiz Antonio Maccari Junior <luiz>
Submitted:  Mon 12 Apr 2021 05:15:12 PM UTC
   
 
Category:  Octave Function Severity:  3 - Normal
Priority:  5 - Normal Item Group:  Incorrect Result
Status:  Fixed Assigned to:  None
Originator Name:  Open/Closed:  * Closed
Release:  * dev Operating System:  * Any
Fixed Release:  None Planned Release:  None
* Mandatory Fields

Add a New Comment Rich Markup
   

Jump to the original submission

Thu 15 Apr 2021 07:08:36 AM UTC, comment #12: 

Thank you for testing.

I grafted the change to the stable branch here:
https://hg.savannah.gnu.org/hgweb/octave/rev/6dfc06f55cd2

Markus Mützel <mmuetzel>
Group administrator
Wed 14 Apr 2021 08:17:08 PM UTC, comment #11: 

Per Luiz in comment #9, bug is now fixed.  Closing report.

Rik <rik5>
Group administrator
Wed 14 Apr 2021 05:03:49 PM UTC, comment #10: 

@Markus: Since this actually fixes an error in Octave behavior, I think this is appropriate to backport to the stable branch.

Rik <rik5>
Group administrator
Wed 14 Apr 2021 03:48:24 PM UTC, comment #9: 

Thanks. Now it is working correctly.

Luiz Antonio Maccari Junior <luiz>
Wed 14 Apr 2021 03:22:21 PM UTC, comment #8: 

Since these are just changes in .m files, you could download the changed files from the repository:
https://hg.savannah.gnu.org/hgweb/octave/file/tip/scripts/polynomial

Afaict, they should still be compatible with Octave 6.

If you cannot replace the original files for some reason, it should be possible to shadow them by adding the new functions to your load path (see e.g. `addpath`).

Markus Mützel <mmuetzel>
Group administrator
Wed 14 Apr 2021 03:10:43 PM UTC, comment #7: 

Is there a way to apply this patch directly in my installed Octave? Or I have to manually change the indicated lines in mpoles.m and residue.m?

Luiz Antonio Maccari Junior <luiz>
Wed 14 Apr 2021 01:31:13 PM UTC, comment #6: 

I added some BISTs, made sure the test suite still passes, and pushed the patch here:
https://hg.savannah.gnu.org/hgweb/octave/rev/852489d1fcb8

Marking as ready to test.

Markus Mützel <mmuetzel>
Group administrator
Wed 14 Apr 2021 01:09:14 PM UTC, comment #5: 

Turns out the error occurred because one of the poles was orders of magnitude larger than the others and it was negative.

Should be fixed with the attached patch.

That changes the order of poles as returned by `mpoles` to be sorted in descending order of magnitude. That seems to match what Matlab is doing.
Before that change Octave returned the poles in descending order of value.


(file #51264)

Markus Mützel <mmuetzel>
Group administrator
Wed 14 Apr 2021 12:40:14 PM UTC, comment #4: 

I believe, the root of the issue is that `mpoles` incorrectly classifies the two distinct roots in the original example as a single root with multiplicity 2:

A = [1, 1.100000536842105e+04, 1.703789473684211e+03, 0];
poles1 = roots (A);
[e, idx] = mpoles (poles1, 0.001, 1)


In Octave:

>> [e, idx] = mpoles (poles1, 0.001, 1)
e =

   1
   2
   1

idx =

   3
   2
   1


I'd guess that this might be because one of the poles differs in orders of magnitude from the other ones.

Re-titling report.

Fwiw, in Matlab:

>> A = [1, 1.100000536842105e+04, 1.703789473684211e+03, 0];
poles1 = roots (A);
[e, idx] = mpoles (poles1, 0.001, 1)

e =

     1
     1
     1


idx =

     2
     3
     1


Markus Mützel <mmuetzel>
Group administrator
Wed 14 Apr 2021 10:34:47 AM UTC, comment #3: 

I think that the new title is not correct at all. In test 1 example the multiplicity is not greater than one. Only the wrong result is. Could you change title?

Luiz Antonio Maccari Junior <luiz>
Tue 13 Apr 2021 04:49:00 PM UTC, comment #2: 

There are two issues here.  The first is rather simple.  The argument 'e' is an extra Octave-only feature.  It gives information about the multiplicity of the poles.  The Matlab function of the same name does not have this extra feature.  The input is not required and the documentation just needs to be changed to state that this is optional.

The second issue is more substantive.  I've re-written the example and attached it as tst_residue.m.  The code is


B = [1315.789473684211];
A = [1, 1.100000536842105e+04, 1.703789473684211e+03, 0];

poles1 = roots (A)

[r,p,k,e] = residue (B, A)

[B1, A1] = residue (r, p, k, e)


This vividly illustrates the problem because residue() is not able to reproduce the input vectors B and A.  It is not just the values, the dimensions are wrong too.  B is a scalar (1x1) whereas B1 is a vector (2x1).



(file #51257)

Rik <rik5>
Group administrator
Mon 12 Apr 2021 05:26:46 PM UTC, comment #1: 

In matlab it is not necessary to use the argument e in the function residue. Suppressing it the result is the expected.


original submission:

> Function residue is presenting unexpected results. I have made two examples to show the problem.
>
>


> % Residue test
> %test 1
>
> num1 = [1315.789473684211];
> den1 =  [1   1.100000536842105e+04  1.703789473684211e+03   0];
>
> poles1 = roots(den1)
>
> [r1,p1,k1,e1]=residue(num1,den1)
>
> %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
>
> % test 2
> num2 = [1];
> den2 =  [1   1   1   0];
>
> poles2 = roots(den2)
> [r2,p2,k2,e2]=residue(num2,den2)


>
>
>
> Test 1 presents an wrong result. It was expected p1 equal a poles1 but I got
>
> p1 =
>
>   -7.744602880490632e-02
>   -7.744602880490632e-02
>   -1.099985047636344e+04
>
> poles1 =
>
>   -1.099985047636344e+04
>   -1.548920576098126e-01
>                        0
>
> In test2 the correct result are occurring. p2 equals poles2.
>
> p2 =
>
>   -0.500000000000000 + 0.866025403784438i
>   -0.500000000000000 - 0.866025403784438i
>                    0 +                 0i
>
> poles2 =
>
>   -0.500000000000000 + 0.866025403784438i
>   -0.500000000000000 - 0.866025403784438i
>                    0 +                 0i
>
> I cant find where is the problem. In matlab the code works well.
>
>

Luiz Antonio Maccari Junior <luiz>
Mon 12 Apr 2021 05:15:12 PM UTC, original submission:  

Function residue is presenting unexpected results. I have made two examples to show the problem.



% Residue test
%test 1

num1 = [1315.789473684211];
den1 =  [1   1.100000536842105e+04  1.703789473684211e+03   0];

poles1 = roots(den1)

[r1,p1,k1,e1]=residue(num1,den1)

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% test 2
num2 = [1];
den2 =  [1   1   1   0];

poles2 = roots(den2)
[r2,p2,k2,e2]=residue(num2,den2)




Test 1 presents an wrong result. It was expected p1 equal a poles1 but I got

p1 =

  -7.744602880490632e-02
  -7.744602880490632e-02
  -1.099985047636344e+04

poles1 =

  -1.099985047636344e+04
  -1.548920576098126e-01
                       0

In test2 the correct result are occurring. p2 equals poles2.

p2 =

  -0.500000000000000 + 0.866025403784438i
  -0.500000000000000 - 0.866025403784438i
                   0 +                 0i

poles2 =

  -0.500000000000000 + 0.866025403784438i
  -0.500000000000000 - 0.866025403784438i
                   0 +                 0i

I cant find where is the problem. In matlab the code works well.


Luiz Antonio Maccari Junior <luiz>

 

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

Attach Files:
   
   
Comment:
   

Attached Files
file #51264:  bug60384_mpoles.patch added by mmuetzel (2KiB - application/octet-stream)
file #51257:  tst_residue.m added by rik5 (165B - text/x-matlab)

 

Depends on the following items: None found

Items that depend on this one: None found

 

Carbon-Copy List
  • -email is unavailable- added by mmuetzel (Posted a comment)
  • -email is unavailable- added by rik5 (Updated the item)
  • -email is unavailable- added by luiz (Submitted the item)
  •  

    There are 0 votes so far. Votes easily highlight which items people would like to see resolved in priority, independently of the priority of the item set by tracker managers.

    Only group members can vote.

     

    Follow 11 latest changes.

    Date Changed by Updated Field Previous Value => Replaced by
    2021-04-14 rik5 StatusReady For Test Fixed
        Open/ClosedOpen Closed
    2021-04-14 mmuetzel StatusPatch Submitted Ready For Test
    2021-04-14 mmuetzel Attached File- Added bug60384_mpoles.patch, #51264
        StatusConfirmed Patch Submitted
    2021-04-14 mmuetzel Operating SystemGNU/Linux Any
        SummaryResidue function produces incorrect results when multiplicity of poles &gt; 1 `mpoles` incorrectly classifies distinct poles as unique ones
    2021-04-13 rik5 StatusNone Confirmed
        Release6.1.0 dev
        SummaryResidue Function incorrect results Residue function produces incorrect results when multiplicity of poles > 1
    2021-04-13 rik5 Attached File- Added tst_residue.m, #51257

    Back to the top

    Powered by Savane 3.13-d3ae.
    Corresponding source code