bugGNU Octave - Bugs: bug #59416, ode15s produce results with...

 
 

bug #59416: ode15s produce results with inconsestance size when compared to other ode solver or matlab.

Submitter:  None
Submitted:  Thu 05 Nov 2020 02:48:04 AM UTC
   
 
Category:  Octave Function Severity:  3 - Normal
Priority:  5 - Normal Item Group:  Matlab Compatibility
Status:  Fixed Assigned to:  None
Originator Name:  Bawazeer Originator Email:  -email is unavailable-
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

Tue 10 Nov 2020 04:27:01 PM UTC, comment #24: 

Seems like this has been addressed.  Someone who is interested in Matlab compatibility should probably also test whether the functions ode23.m, ode45.m, ode23s.m also transpose their event outputs when only a single output struct is requested.  It may be that Matlab has changed all of their solvers to do this recently.  But, that would be the subject of a new bug report.  This one has been fixed so I'm closing it now.

Rik <rik5>
Group administrator
Tue 10 Nov 2020 07:55:10 AM UTC, comment #23: 

Rik transposed the fields in the return structure that correspond to events here:
https://hg.savannah.gnu.org/hgweb/octave/rev/2190720bca3e

Can this be closed as fixed?

Markus Mützel <mmuetzel>
Group administrator
Sat 07 Nov 2020 03:49:16 PM UTC, comment #22: 

@Rik: The script fails. But maybe the output when removing the semicolon after the `ode15s` command is good enough:

>> ode15s_event_test2

sol =

  struct with fields:

     solver: 'ode15s'
    extdata: [1×1 struct]
          x: [1×38 double]
          y: [3×38 double]
         xe: [10.0000 10.0000]
         ye: [3×2 double]
         ie: [1 2]
      stats: [1×1 struct]
      idata: [1×1 struct]

Reference to non-existent field 'te'.

Error in ode15s_event_test2 (line 17)
sol.te


Markus Mützel <mmuetzel>
Group administrator
Sat 07 Nov 2020 03:41:59 PM UTC, comment #21: 

@Markus: Since you seem to have access to the newest release (R2020b), could you try the attached ode15s_event_test2.m file?

This tests the question from comment #20 about what the orientation is when ode15s is called with just one struct output.



(file #50221)

Rik <rik5>
Group administrator
Sat 07 Nov 2020 12:48:59 PM UTC, comment #20: 

Thank you very much. With this I withdraw my comment #16.

The only topic that remains is that in R2018b "sol.xe", "sol.ye" in the structure array is a rowvector (comment #9 point 2). The Matlab documentation says that "sol.x" is a rowvector, but it does not say anything about "sol.xe". But it would make sense, if we return "sol.x" as rowvector, that we should do the same for "sol.xe" and "sol.ye".

The same applies to ode15i.

Hg200 <hg200>
Sat 07 Nov 2020 08:31:27 AM UTC, comment #19: 

In Matlab R2020b:

>> ode15s_event_test
te =
   10.0000
   10.0000
ye =
    0.8415    0.0000    0.1585
    0.8415    0.0000    0.1585
ie =
     1
     2


Markus Mützel <mmuetzel>
Group administrator
Sat 07 Nov 2020 02:33:17 AM UTC, comment #18: 

I am the Originator. I run the code in comment #16 on R2019b.


[t, y, te, ye, ie] = ode15s (@rob, [0, 100], [1; 0; 0], opt);
ie
ye
te
size(y)
size(t)


I found the following results
ie =

     1     2

ye =

      0.84147   1.6243e-05      0.15852
      0.84147   1.6243e-05      0.15852

te =

           10
           10

ans =

    38     3

ans =

    38     1

Anonymous
Fri 06 Nov 2020 10:04:48 PM UTC, comment #17: 

I second request by Hg200 to test with a more modern version of Matlab.  Matlab has been known to have incorrect documentation, and thus it can't be fully trusted.  The only true reference for how Matlab behaves is to actually run code on its interpreter.

Rik <rik5>
Group administrator
Fri 06 Nov 2020 07:02:46 PM UTC, comment #16: 

Ok. Step by step. Attached is a minimum example, taken from the Octave BISTs. If i run this example in R2018b a row vector for the event index "ie" is returned. But the documentation says it must be column vector. Can please anybody check with a newer Matlab version?


te =
   10.0000
   10.0000

ye =
    0.8415    0.0000    0.1585
    0.8415    0.0000    0.1585

ie =
    1     2




(file #50217)

Hg200 <hg200>
Fri 06 Nov 2020 05:47:45 PM UTC, comment #15: 

According to Matlab documentation, te and ie are always column vectors.  Rather than transpose the output from ode15s, it would probably be better to figure out why row vectors are getting generated.  If I look at the code in _ode15_.cc I see that te and ie are declared to be of type ColumnVector which appears correct.  What is changing these to RowVector?

Rik <rik5>
Group administrator
Fri 06 Nov 2020 05:22:38 PM UTC, comment #14: 

I forgot one: Events "ie" must be also transposed. This line is also from the BISTs in ode15s:


[t, y, te, ye, ie] = ode15s (@rob, [0, 100], [1; 0; 0], opt);


Matlab:

ie =

   1
   2

Octave:

ie =

   1   2


Hg200 <hg200>
Fri 06 Nov 2020 04:40:59 PM UTC, comment #13: 

@Markus: thank you for fixing up the BIST tests.  The inconsistent behavior to return [y, t] as column vectors (2-output form of calling function), but sol.x, sol.y as row vectors (1-output form of calling function) is required for Matlab compatibility.  See documentation at https://www.mathworks.com/help/matlab/ref/ode15s.html.

I submitted a patch for bug #59063, but again can't test it.  If it works I would be happy to see someone apply it.


Rik <rik5>
Group administrator
Fri 06 Nov 2020 02:25:32 PM UTC, comment #12: 

Hmm, probably better to separate that. ode15i needs to be fixed with respect to event handling bug #59063 due to its BISTs.

Yes, point 3 is invalid.

Hg200 <hg200>
Fri 06 Nov 2020 01:52:25 PM UTC, comment #11: 

@hg200: Please, try to not mix separate issues into a single changeset. I also believe we should continue to test whether passing functions by name works.

I pushed a change that should be fixing the BISTs slightly differently here:
https://hg.savannah.gnu.org/hgweb/octave/rev/bf62eeabf6d1

I don't see what needs to be changed for ode15i.

Re-opening for the point 2.) in comment #9 which seems to be quite similar to the original issue.

Point 3.) in that comment is probably invalid: 163 and 165 are "the same" with a tolerance of 3.

Markus Mützel <mmuetzel>
Group administrator
Fri 06 Nov 2020 01:11:34 PM UTC, comment #10: 

@Markus: Too late for ode15s ;-) Do you want to continue with ode15i? I think we should put everything together (ode15i, ode15s and bug #59063) into one patch, so that we have one clean commit in the mercurial repository?

Hg200 <hg200>
Fri 06 Nov 2020 01:01:25 PM UTC, comment #9: 

I have converted the BISTs of ode15s.m to Matlab and checked with R2018b. I found a few more problems, but i have also good news ;-)

1.)


## Events option add further elements in sol
%!testif HAVE_SUNDIALS
%! opt = odeset ("Events", @feve, "Mass", @massdensefunstate,
%!               "MStateDependence", "none");
%! sol = ode15s (@rob, [0, 100], [1; 0; 0], opt);
%! assert (isfield (sol, "ie"));
%! assert (sol.ie, [0;1]);
%! assert (isfield (sol, "xe"));
%! assert (isfield (sol, "ye"));
%! assert (sol.x(end), 10, 1);


In Matlab "assert (sol.ie, [0;1])" fails because ode15s returns [1;2] instead of [0;1]. But, here we have bug #59063 and Rik already provided a patch for this. Applying Rik's patch --> Problem is solved.

2.)


sol =

  struct with fields:

     solver: 'ode15s'
    extdata: [1x1 struct]
          x: [1x38 double]
          y: [3x38 double]
         xe: [10.0000 10.0000]
         ye: [3x2 double]
         ie: [2x1 double]
      stats: [1x1 struct]
      idata: [1x1 struct]


This is from the example above in Matlab. I think sol.xe and sol.ie must be also transposed -> TODO??

3.)


## Refine
%!testif HAVE_SUNDIALS
%! opt2 = odeset ("Refine", 3, "Mass", @massdensefunstate,
%!                "MStateDependence", "none");
%! opt1 = odeset ("Mass", @massdensefunstate, "MStateDependence", "none");
%! [t, y] = ode15s (@rob, [0, 100], [1; 0; 0], opt1);
%! [t2, y2] = ode15s (@rob, [0, 100], [1; 0; 0], opt2);
%! assert (numel (t2), numel (t) * 3, 3);


This fails in Matlab because


numel (t2)
ans =
   163
numel (t)*3
ans =
   165


But i think this issue is acceptable for Octave?

I have attached a patch for ode15s. The patch will only pass if fix from bug #59063 is applied.

ode15i still open.

(file #50215)

Hg200 <hg200>
Fri 06 Nov 2020 12:42:21 PM UTC, comment #8: 

In Matlab R2020b:

>> sol=ode15s(@(t,y)y-t^2+1,[0,2]',1); size(sol.y)
ans =
     1    16
>> [t,y]=ode15s(@(t,y)y-t^2+1,[0,2]',1); size(y)
ans =
    16     1


So the new behavior seems to be Matlab compatible.

@hg200: Don't worry diving too deep into this. I'm having a changeset almost ready...
I'm also not sure if I correctly understand: We need BISTs for functions passed by name and by function handle.

Markus Mützel <mmuetzel>
Group administrator
Fri 06 Nov 2020 10:34:42 AM UTC, comment #7: 

Rik is correct (comment #4) and also Markus is correct (comment #6)

When executed under Matlab R2018b one of the Octave ode15s BISTs returns


ref = [-1.205364552835178, 0.951542399860817];
sol = ode15s (@fpol, [0 -2], [2, 0]);
assert ([sol.x(end), sol.y(end,:)], [-2, ref], 5e-3);
size([sol.x(end), sol.y(end,:)])
ans =
     1    37


So we have to adopt the assert line :


assert ([sol.x(end), sol.y(end,:)], [-2, ref], 5e-3);

must be

assert ([sol.x(end), sol.y(:,end).'], [-2, ref], 5e-3);


And the output of the functions ode15i and ode15s is weird:


[t, x] = ode15i (@f, tspan, x0, x0dot, options);
size(t)
ans =
   557     1
size(x)
ans =
   557     2

But

z = ode15i (@f, tspan, x0, x0dot, options);
  struct with fields:

     solver: 'ode15i'
    extdata: [1x1 struct]
          x: [1x557 double]
          y: [2x557 double]
      stats: [1x1 struct]
      idata: [1x1 struct]


I also see some other issues in the BISTs. References for functions are sometimes "strings" and sometimes "@" (mixed). That should be corrected too. I am going through the BISTs and check for errors.

Hg200 <hg200>
Fri 06 Nov 2020 07:14:39 AM UTC, comment #6: 

I started to have a look at this.

It seems only the built-in self tests (BISTs) where ode15s is called with one return argument - like in `sol = ode15s (...)` - need transposing.
When ode15s is called with two output arguments - like in `[t, y] = ode15s (...)` - the output is still a row vector.

Is this the expected behavior? It looks kind of inconsistent - but might be necessary for Matlab compatibility...

Markus Mützel <mmuetzel>
Group administrator
Thu 05 Nov 2020 09:58:12 PM UTC, comment #5: 
Dmitri A. Sergatskov <dasergatskov>
Thu 05 Nov 2020 09:35:08 PM UTC, comment #4: 

This makes sense.  Previously, the functions returned column vectors for struct outputs and now they return row vectors.  The values will be the same but the orientation is different.

I don't have the necessary libraries to run the tests for ode15s and ode15i, but it should be simple enough to find the failures and then transpose the expected results in the BIST tests.

I can even do it if you post the results from fntests.log that have the failures.

Rik <rik5>
Group administrator
Thu 05 Nov 2020 07:35:29 PM UTC, comment #3: 

(ah my keboard)
make that 12 FAILs for ode15s.m

Philip Nienhuis <philipnienhuis>
Group Member
Thu 05 Nov 2020 07:34:26 PM UTC, comment #2: 

'make check' now yields 15 FAILs for ode15i

Philip Nienhuis <philipnienhuis>
Group Member
Thu 05 Nov 2020 03:54:56 PM UTC, comment #1: 

This was an easy change.  Apparently both ode15s and ode15i were not transposing the solution when the output was a structure.  See this changeset: http://hg.savannah.gnu.org/hgweb/octave/rev/7f103819617d.

Marking as fixed and closing report.

Rik <rik5>
Group administrator
Thu 05 Nov 2020 02:48:04 AM UTC, original submission:  



sol=ode15s(@(t,y)y-t^2+1,[0,2]',1); size(sol.y)
sol=ode45(@(t,y)y-t^2+1,[0,2]',1); size(sol.y)


Octave results
ans =

   26    1

ans =

    1   13

Matlab results

ans =

    1   16

ans =

    1   12

It is clear that ode15s function produced transposed results.

Anonymous

 

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

Attach Files:
   
   
Comment:
   

Attached Files
file #50221:  ode15s_event_test2.m added by rik5 (1KiB - text/x-matlab)
file #50217:  ode15s_event_test.m added by hg200 (1KiB - text/x-objcsrc)
file #50215:  patch_ode15s.diff added by hg200 (6KiB - text/x-patch)

 

Depends on the following items: None found

Items that depend on this one: None found

 

Carbon-Copy List
  • -email is unavailable- added by hg200 (Posted a comment)
  • -email is unavailable- added by mmuetzel (Posted a comment)
  • -email is unavailable- added by dasergatskov (Posted a comment)
  • -email is unavailable- added by rik5 (Posted a comment)
  • -email is unavailable- added by None (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
    2020-11-10 rik5 StatusIn Progress Fixed
        Open/ClosedOpen Closed
    2020-11-07 rik5 Attached File- Added ode15s_event_test2.m, #50221
    2020-11-06 hg200 Attached File- Added ode15s_event_test.m, #50217
    2020-11-06 mmuetzel StatusNeed Info In Progress
        Open/ClosedClosed Open
    2020-11-06 hg200 Attached File- Added patch_ode15s.diff, #50215
    2020-11-06 mmuetzel Operating SystemGNU/Linux Any
    2020-11-06 mmuetzel StatusFixed Need Info
    2020-11-05 rik5 StatusNone Fixed
        Open/ClosedOpen Closed

    Back to the top

    Powered by Savane 3.13-4448.
    Corresponding source code