@Allan Wittkopf 

I get it!

This has been a very productive exchange for me,

thanks Allan

@Allan Wittkopf 

Fantastic!

I'm ready to move on and play around with my simulations. Thanks a lot Allan.

I understand that 1e-16 is chosen to be so close to 0 that no event will be triggered between 0 and 1e-16 (though that may not be true in all applications, as you mentioned cases where the event is triggered millions of times!) and being not exactly zero it also sets the direction of integration. On my machine, I found that 1e-323 still did the job, while 1e-324 was too small.

Below some not-terribly urgent thoughts:

Is there another way to reset the integrator? I tried to add lines like forget(dsn) before or after the loop, or to forget(EventDates) applied to the procedure EventDates itself, but that does nothing. When I query the cache of the dsn procedure, it seems to be empty. So how does Maple store/remember previous stages of the integration? I read dsolve/numeric and dsolve/numeric/ivp on that, I don't fully understand how the dsolve procedure works. Here's my guess: Is it the case that every time initial conditions are changed in dsn('initial'=NewInitialConditions) this REDEFINES the procedure altogether, replacing the original initial conditions with the new ones?

@Allan Wittkopf 

Fantastic!

I'm ready to move on and play around with my simulations. Thanks a lot Allan.

I understand that 1e-16 is chosen to be so close to 0 that no event will be triggered between 0 and 1e-16 (though that may not be true in all applications, as you mentioned cases where the event is triggered millions of times!) and being not exactly zero it also sets the direction of integration. On my machine, I found that 1e-323 still did the job, while 1e-324 was too small.

Below some not-terribly urgent thoughts:

Is there another way to reset the integrator? I tried to add lines like forget(dsn) before or after the loop, or to forget(EventDates) applied to the procedure EventDates itself, but that does nothing. When I query the cache of the dsn procedure, it seems to be empty. So how does Maple store/remember previous stages of the integration? I read dsolve/numeric and dsolve/numeric/ivp on that, I don't fully understand how the dsolve procedure works. Here's my guess: Is it the case that every time initial conditions are changed in dsn('initial'=NewInitialConditions) this REDEFINES the procedure altogether, replacing the original initial conditions with the new ones?

@DJ Clayworth 

is that 3D plotting algorithm going to be rewritten in Maple's next release or in the near future?

@Allan Wittkopf 

Thank you very much Allan for your patient explanations. Wanting to make this process more systematic, I have written a loop to replicate the steps you have detailed above.

tf := 100:
interface(warnlevel=0):
A := Array([0]):
while A[-1] < tf do
  dsn(tf):
  A(ArrayNumElems(A)+1) := dsn(eventfired=[1])[1]:
  dsn(eventclear):
  end do:
interface(warnlevel=3):
A := convert(A,list):
A[2..-1];

I have used an array instead of a list to be able to deal with a large number of elements. I have silenced the warnings. I initially created the Array with the 0 elements, which I later remove from the final list.

It's intriguing to see that even with a simple problem such as this roundoff errors lead to the values 17-31 being a little off their exact values (I used default digits), but the other values are quite precisely as expected.

I intend to place the loop above into a procedure, for convenience. But something isn't quite right. If I execute the loop above once, I get the results I reported. But if I execute it twice, instead I get the error message: Error, (in _dat[1]) no events to clear

To avoid the error message, I need to rerun dsolve to redefine dsn every time before I execute the loop. It's as if running the loop altered the dsn procedure. I tried to get around that by executing the loop on a copy of dsn, e.g. sol:=dsn: thinking that this way dsn would remain unchanged, but that didn't work. I did it inside a procedure, but that didn't help either.

Here was my attempt at a procedure, based on the loop above:

EventDates := proc(sol,tf)
  local A;
  description "the first argument is dsolve(), the second argument is the end-date";
  interface(warnlevel=0);
  A := Array([0]);
  while A[-1] < tf do
    sol(tf);
    A(ArrayNumElems(A)+1) := sol(eventfired=[1])[1];
    sol(eventclear);
  end do;
  A := convert(A,list);
  A[2..-1];
end proc:

A := EventDates(dsn,100);

As with the loop above, works once, not twice.

It does work once, so I can live with it, you know. Thanks!

@Allan Wittkopf 

Thank you very much Allan for your patient explanations. Wanting to make this process more systematic, I have written a loop to replicate the steps you have detailed above.

tf := 100:
interface(warnlevel=0):
A := Array([0]):
while A[-1] < tf do
  dsn(tf):
  A(ArrayNumElems(A)+1) := dsn(eventfired=[1])[1]:
  dsn(eventclear):
  end do:
interface(warnlevel=3):
A := convert(A,list):
A[2..-1];

I have used an array instead of a list to be able to deal with a large number of elements. I have silenced the warnings. I initially created the Array with the 0 elements, which I later remove from the final list.

It's intriguing to see that even with a simple problem such as this roundoff errors lead to the values 17-31 being a little off their exact values (I used default digits), but the other values are quite precisely as expected.

I intend to place the loop above into a procedure, for convenience. But something isn't quite right. If I execute the loop above once, I get the results I reported. But if I execute it twice, instead I get the error message: Error, (in _dat[1]) no events to clear

To avoid the error message, I need to rerun dsolve to redefine dsn every time before I execute the loop. It's as if running the loop altered the dsn procedure. I tried to get around that by executing the loop on a copy of dsn, e.g. sol:=dsn: thinking that this way dsn would remain unchanged, but that didn't work. I did it inside a procedure, but that didn't help either.

Here was my attempt at a procedure, based on the loop above:

EventDates := proc(sol,tf)
  local A;
  description "the first argument is dsolve(), the second argument is the end-date";
  interface(warnlevel=0);
  A := Array([0]);
  while A[-1] < tf do
    sol(tf);
    A(ArrayNumElems(A)+1) := sol(eventfired=[1])[1];
    sol(eventclear);
  end do;
  A := convert(A,list);
  A[2..-1];
end proc:

A := EventDates(dsn,100);

As with the loop above, works once, not twice.

It does work once, so I can live with it, you know. Thanks!

@Joe Riel 

This is a very clear explanation Joe, thanks a lot. I had read the original exchange I quoted above and got a sense that the matter was more subtle than my characterization, but your description has cleared it all up.

The use of the round parentheses (instead of square brackets) had been lost on me, actually. I only now notice it.

Thanks Joe!

I wish I could 1) upvote your comment, 2) make it a favorite, 3) recommend your post as inclusion a list of FAQs, but currently I can do none of these 3 things.

@Joe Riel 

This is a very clear explanation Joe, thanks a lot. I had read the original exchange I quoted above and got a sense that the matter was more subtle than my characterization, but your description has cleared it all up.

The use of the round parentheses (instead of square brackets) had been lost on me, actually. I only now notice it.

Thanks Joe!

I wish I could 1) upvote your comment, 2) make it a favorite, 3) recommend your post as inclusion a list of FAQs, but currently I can do none of these 3 things.

can you copy-paste your input into this box, rather than attaching a document?

@acer 

oh interesting, I didn't know that, thanks, no wonder it's hard to find as this list appears under maplesoft rather than mapleprimes...

@Allan Wittkopf 

I've been trying to apply your suggestions and have a couple of follow-up questions if I may.

Above you wrote

> dsn(6.1);
                [t = 6.1, y(t) = 0.900000000000001, s(t) = 6.]

> dsn(eventfired=[1]);
                              [6.00000000000000]

the query dsn(eventfired=[1]); depends upon the previous call dsn(6.1); and returns ONLY THE LAST DATE at which the event fired on the interval [0,6.1]. I would like a list of EVERY DATE at which the event has triggered on that interval. Feasible?

Related to my failed attempts in generating a list, a question came up. The help page states that one can use eventfired without previously requiring a solution value, but one needs to specify the direction of integration. No example is provided. How do you do that? I tried to add 'direction'=1 to the call to dsolve. I also tried to write things like Sol(direction=1); and variations thereof. Obviously that's not how to set the direction of integration. How do you set "direction=1" ?

To sum up, my questions are:

1. Can I use eventfired to obtain EVERY DATE at which the event has triggered on the interval [0,6.1] ?
2. How can I set the direction of integration without requesting a specific solution value?

Thanks!

Edit. The help page has a paragraph which provides a hint that it may be feasible, namely:

The calling sequence eventfired=[n1=v1,n2=v2,n3=v3] can be used to specify when each event has last fired (and as a result the events will not fire until one of the triggers fire pas the relevant specified value). When this form is used, the prior values of eventfired are returned in a list.

In bold, the bit that seems to suggest that values of eventfired prior to the last one may be returned. Unfortunately I'm too thick to implement the suggestion. The n1 stands for the 1 we have in eventfired=[1], I think. But not sure where to get the v1 ...

@Allan Wittkopf 

I've been trying to apply your suggestions and have a couple of follow-up questions if I may.

Above you wrote

> dsn(6.1);
                [t = 6.1, y(t) = 0.900000000000001, s(t) = 6.]

> dsn(eventfired=[1]);
                              [6.00000000000000]

the query dsn(eventfired=[1]); depends upon the previous call dsn(6.1); and returns ONLY THE LAST DATE at which the event fired on the interval [0,6.1]. I would like a list of EVERY DATE at which the event has triggered on that interval. Feasible?

Related to my failed attempts in generating a list, a question came up. The help page states that one can use eventfired without previously requiring a solution value, but one needs to specify the direction of integration. No example is provided. How do you do that? I tried to add 'direction'=1 to the call to dsolve. I also tried to write things like Sol(direction=1); and variations thereof. Obviously that's not how to set the direction of integration. How do you set "direction=1" ?

To sum up, my questions are:

1. Can I use eventfired to obtain EVERY DATE at which the event has triggered on the interval [0,6.1] ?
2. How can I set the direction of integration without requesting a specific solution value?

Thanks!

Edit. The help page has a paragraph which provides a hint that it may be feasible, namely:

The calling sequence eventfired=[n1=v1,n2=v2,n3=v3] can be used to specify when each event has last fired (and as a result the events will not fire until one of the triggers fire pas the relevant specified value). When this form is used, the prior values of eventfired are returned in a list.

In bold, the bit that seems to suggest that values of eventfired prior to the last one may be returned. Unfortunately I'm too thick to implement the suggestion. The n1 stands for the 1 we have in eventfired=[1], I think. But not sure where to get the v1 ...

And here the answer I was asking myself: why create an Array that is later converted to a list

http://www.mapleprimes.com/questions/39483-Memory-Errors-List-Vs-Array#comment72012

the reason: lists are limited to 100 elements.

And here the answer I was asking myself: why create an Array that is later converted to a list

http://www.mapleprimes.com/questions/39483-Memory-Errors-List-Vs-Array#comment72012

the reason: lists are limited to 100 elements.

It looks like this test is just measuring the distance from Rubi's results, where Rubi is taken as the cannon.

Isn't the road to success in these kind of tests simply to define new functions? Would a solution in terms of the Lambert-W function be deemed more simple than an equivalent representation just because it has been named?