Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Is it possible to access the partial list content while being Sow'ed or must one wait for it to be Reap'ed?

I've been learning Sow/Reap. They are cool constructs. But I need help to see if I can use them to do what I will explain below.

What I'd like to do is: Plot the solution of NDSolve as it runs. I was thinking I can use Sow[] to collect the solution (x,y[x]) as NDSolve runs using EvaluationMonitor. But I do not want to wait to the end, Reap it and then plot the solution, but wanted to do it as it is running.

I'll show the basic setup example

max = 30;
sol1 = y /. 
   First@NDSolve[{y'[x] == y[x] Cos[x + y[x]], y[0] == 1}, 
     y, {x, 0, max}];
Plot[sol1[x], {x, 0, max}, PlotRange -> All, AxesLabel -> {"x", "y[x]"}]

enter image description here

Using Reap/Sow, one can collect the data points, and plot the solution at the end like this

sol = Reap[
    First@NDSolve[{y'[x] == y[x] Cos[x + y[x]], y[0] == 1}, 
      y, {x, 0, max}, EvaluationMonitor :> Sow[{x, y[x]}]]][[2, 1]];

ListPlot[sol, AxesLabel -> {"x", "y[x]"}]

enter image description here

Ok, so far so good. But what I want is to access the partially being build list, as it accumulates by Sow and plot the solution. The only setup I know how do this is to have Dynamic ListPlot that refreshes when its data changes. But I do not know how to use Sow to move the partially build solution to this data so that ListPlot update.

I'll show how I do it without Sow, but you see, I am using AppenedTo[] in the following:

ClearAll[x, y, lst];
max = 30;
lst = {{0, 0}};
Dynamic[ListPlot[lst, Joined -> False, PlotRange -> {{0, max}, All}, 
  AxesLabel -> {"x", "y[x]"}]]

NDSolve[{y'[x] == y[x] Cos[x + y[x]], y[0] == 1}, y, {x, 0, max}, 
 EvaluationMonitor :> {AppendTo[lst, {x, y[x]}]; Pause[0.01]}]

enter image description here

I was thinking of a way to access the partially build list by Sow and just use that to refresh the plot, on the assumption that might be more efficient than AppendTo[]

I can't just do this:

ClearAll[x, y, lst];
max = 30;
lst = {{0, 0}};
Dynamic[ListPlot[lst, Joined -> False, PlotRange -> All]]

NDSolve[{y'[x] == y[x] Cos[x + y[x]], y[0] == 1}, y, {x, 0, max}, 
 EvaluationMonitor :> {lst = Reap[Sow[{x, y[x]}] ][[2, 1]]; Pause[0.01]}]

Since it now Sow one point, and Reap it, so I am just plotting one point at a time. The same as if I just did:

NDSolve[{y'[x] == y[x] Cos[x + y[x]], y[0] == 1}, y, {x, 0, max}, 
 EvaluationMonitor :> {lst = Sow[{x, y[x]}]; Pause[0.01]}]

my question is, how to use Sow/Reap in the above, to avoid me having manage the lst by the use of AppendTo in this case. (or by pre-allocation using Table, but then I would not know the size to allocate) Since I assume that may be Sow/Reap would be more efficient?

ps. What would be nice, if Reap had an option to tell it to Reap what has been accumulated by Sow, but do not remove it from what has been Sow'ed so far. Like a passive Reap sort of. Well, just a thought.

thanks

Update: 8:30 am

Thanks for the answers and comments. I just wanted to say, that the main goal of asking this was just to see if there is a way to access part of the data while being Sowed. I need to look more at Bag, I have not used that before.

Btw, The example shown above, was just to give a context to where such a need might appear. If I wanted to simulate the solution in this specific case, I do not even have to do it as I did, I could obtain the solution data first, then, afterwords, animate it.

Hence no need to even worry about allocation of a buffer myself, or use AppenedTo. But there could many other cases where it will be easier to access the data as it is being accumulated by Sow. This example is just what I had at the moment.

To do this specific example more directly, one can simply used Animate[], afterwords, like this:

Remove["Global`*"];
max = 30;
sol = Reap[
    First@NDSolve[{y'[x] == y[x] Cos[x + y[x]], y[0] == 1}, 
      y, {x, 0, max}, EvaluationMonitor :> Sow[{x, y[x]}]]][[2, 1]];

Animate[ListPlot[sol[[1 ;; idx]], Joined -> False, 
  PlotRange -> {{0, max}, All}, AxesLabel -> {"x", "y[x]"}], {idx, 1, 
  Length[sol], 1}]

Or, even make a home grown animate, like this

Remove["Global`*"];
max = 30;
sol = Reap[
    First@NDSolve[{y'[x] == y[x] Cos[x + y[x]], y[0] == 1}, 
      y, {x, 0, max}, EvaluationMonitor :> Sow[{x, y[x]}]]][[2, 1]];
idx = 1;
Dynamic[idx];
Dynamic[ListPlot[sol[[1 ;; idx]], Joined -> False, 
  PlotRange -> {{0, max}, All}, AxesLabel -> {"x", "y[x]"}]]

Do[++idx; Pause[0.01], {i, 1, Length[sol] - 1}]

Small follow up question: Can one depend on using Internal``Bag now? Since it is in Internal context, will there be a chance it might be removed/changed/etc... in the future, breaking some code? I seems to remember reading somewhere that this is not likely, but I do not feel comfortable using something in Internal Context. If it is Ok for us to use it, why is it in Internal context then?

(so many things to lean in Mathematica, so little time)

Thanks,

like image 532
Nasser Avatar asked Dec 28 '11 05:12

Nasser


1 Answers

Experimentation shows that both Internal`Bag and linked lists are slower than using AppendTo. After considering this I recalled what Sasha told me, which is that list (array) creation is what takes time.

Therefore, neither method above, nor a Sow/Reap in which the result is collected as a list at each step is going to be more efficient (in fact, less) than AppendTo.

I believe that only array pre-allocation can be faster among the native Mathematica constructs.

Old answer below for reference:


I believe this is the place for Internal`Bag, Internal`StuffBag, and Internal`BagPart.

I had to resort to a clumsy double variable method because the Bag does not seem to update inside Dynamic the way I expected.

ClearAll[x, y, lst];
max = 30;
bag = Internal`Bag[];
lst = {{}};

Dynamic@ListPlot[lst, Joined -> False, PlotRange -> All]

NDSolve[{y'[x] == y[x] Cos[x + y[x]], y[0] == 1}, y, {x, 0, max}, 
 EvaluationMonitor :> {Internal`StuffBag[bag, {x, y[x]}];
                       lst = Internal`BagPart[bag, All];
                       Pause[0.01]}
]
like image 58
Mr.Wizard Avatar answered Dec 29 '22 10:12

Mr.Wizard