Correcting time step width of trajectory plots to equal node_time_int…

[chris-konrad]

Fixing #282

[Peter230655]

I am afraid, this may mess up simulations, where the final time is defined as tf = interval * (num_nodes - 1),
unless I misunderstand something.

[chris-konrad]

I don't think this can happen. The results of arrange and linspace with (num_collocation_nodes-1) are equivalent:

In[30]: np.arange(0, num_collocation_nodes * node_time_interval, node_time_interval)
Out[30]: array([0., 1., 2., 3., 4.])

In[31]: np.linspace(0, (num_collocation_nodes-1) * node_time_interval, num=num_collocation_nodes)
Out[31]: array([0., 1., 2., 3., 4.])

It will not really matter which correction we choose. But the current behaviour is definitely not correct :)

[Peter230655]

I don't think this can happen. The results of arrange and linspace with (num_collocation_nodes-1) are equivalent:

In[30]: np.arange(0, num_collocation_nodes * node_time_interval, node_time_interval)
Out[30]: array([0., 1., 2., 3., 4.])

In[31]: np.linspace(0, (num_collocation_nodes-1) * node_time_interval, num=num_collocation_nodes)
Out[31]: array([0., 1., 2., 3., 4.])

It will not really matter which correction we choose. But the current behaviour is definitely not correct :)

...but if an existing simulation uses tf = interval * (num_nodes - 1), then your change will matter?

[chris-konrad]

I am not sure I understand you correctly. Where would a current simulation use tf = interval * (num_nodes - 1) and how would that not be compatible with the proposed change? The current plot function uses interval * num_nodes without the correction by 1.

[Peter230655]

Virtually all the simulations I know or wrote define the starting time t0, mostly t0 = 0 and the final time tf.
Then interval = (tf - t0) / (num_nodes - 1) [you may look at the simulations in examples-gallery]
I guess the advantage is, that you can play around with num_nodes for convergence without messing up the 'physics' of the problem.

As regards compatibility, I do not know for sure.

[chris-konrad]

I see. I don't think you have to worry about this breaking any existing simulations. Here is why:

What you pass to opty when creating a problem is only num_collocation_nodes=num_nodes and node_time_interval=interval. The actual time interval from t0 to tf is calculated from these two parameters in plot_trajectories(). It currently uses the following line to calculate the time array for the plot.

time = np.linspace(0,
                           self.collocator.num_collocation_nodes *
                           node_time_interval,
                           num=self.collocator.num_collocation_nodes)

According to the definition of np.linspace(), the resulting array begins with time[0]=0 and ends with time[-1]=num_collocation_nodes*node_time_interval. In total it will have a length of num_collocation_nodes. This means that the actual interval of the created time array is not interval, it is slightly larger. For example, for t0=0, tf=5, num_nodes=5, interval is, as per your definition, (tf - t0) / (num_nodes - 1) = 1.25. Calling linspace as above means that it will be called as np.linspace(0, 5 * 1.25 = 6.25, num = 5). The resulting array has the elements: [0. , 1.5625, 3.125 , 4.6875, 6.25 ] with the incorrect interval 1.5625. The current code is not doing what you want :)

It doesn't create a problem for the examples because the error only means that the trajectories in the plot are stretched by one sample over the entire plot. You can actually see that in the second plot here. There, the blue traces in the second axis is slightly stretched compared to the true data in orange.

To fix the problem, we can either change the command above to the exact same definition as you use for your simulations:

time = np.linspace(0,
                           (self.collocator.num_collocation_nodes -1) *
                           node_time_interval,
                           num=self.collocator.num_collocation_nodes)

Or we can replace it with np.arange():

        time = np.arange(0,
                         self.collocator.num_collocation_nodes *
                         node_time_interval,
                         node_time_interval)

Although the input to the function looks slightly different, the results of these two variants are exactly equivalent because the way the two functions interpret an interval are defined differently:

• np.array(start, stop, step): [start,stop) -> from including start to the last value smaller then stop
• np.linspace(start, stop, num: [start,stop] -> from including start to including stop.

If you prefer, we can also correct the linspace command by -1 instead of using arange.

[Peter230655]

I have no preference, it is up to Jason to decide about this.
In the examples I know, num_nodes is in there range 300......3000.
So, in practical terms, I guess it does not matter much?

[chris-konrad]

Actually, @Peter230655, turns out there is a good reason to prefer linspace over arange. The latter is less numerically stable when using non-integer steps. This caused two of the examples to fail after all.

So, we should keep using linspace and add the correction by -1 instead.

        time = np.linspace(0,
                           (self.collocator.num_collocation_nodes-1) *
                           node_time_interval,
                           num=self.collocator.num_collocation_nodes)

[chris-konrad]

But yes indeed. The "error" is purely cosmetic and won't really be visible for large sample numbers.

[Peter230655]

Actually, @Peter230655, turns out there is a good reason to prefer linspace over arange. The latter is less numerically stable when using non-integer steps. This actually caused two of the examples to fail after all.

So, we should keep using linspace and add the correction by -1 instead.

        time = np.linspace(0,
                           (self.collocator.num_collocation_nodes-1) *
                           node_time_interval,
                           num=self.collocator.num_collocation_nodes)

I had no idea about the instabilitiy of arange. I am not the biggest Python expert, so I always use linspace.
Looks to me your suggestion makes sense, based on the definition of tf = interval * (num_nodes - 1) as it is done in all examples I know.

[moorepants]

LGTM!

Thanks for your first contribution :)