Excessive delay and memory usage when calling unsupported Numpy functions on TimeSeries.time

Description

Trying to call unsupported Numpy ufuncs or other functions like np.isnan or np.mean on a Time objects raises an exception, offering an (implicit) reminder to directly use the np.array via t.value.
The same operations on an object created as TimeSeries.time for moderately large arrays (length a couple 1000) results in Python running for several minutes and allocating up to several GB (~ len(t) * 200 kB) memory.

Expected behavior

Unsupported functions should raise an exception without undue delay or memory usage of thousands times the original data size, just as for regular Time.

Actual behavior

For a basic TimeSeries like

t = Time(np.arange(10000), format='mjd', scale='tdb')
ts = TimeSeries(time=t)
np.isnan(ts.time)
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
TypeError: ufunc 'isnan' not supported for the input types, and the inputs could not be safely coerced to any supported types according to the casting rule ''safe''

np.abs(ts.time)
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
TypeError: bad operand type for abs(): 'Time'

the last two lines each take 30-40 seconds, ramping memory usage multiple times up to > 2 GB as shown in Activity Monitor or top.

Steps to Reproduce

Running the attached script as (pip install memory_profiler) a function is trying to run a function several times first on the original Time object t, then the same as element ts.time of a TimeSeries, a deepcopy(ts.time) of the latter, and again the original t.

> mprof run python3.8 timeseries_leak_prof.py isnan                                                      
mprof: Sampling memory every 0.1s
1.19.1
4.1rc2
ufunc 'isnan' not supported for the input types, and the inputs could not be safely coerced to any supported types according to the casting rule ''safe''
Filename: timeseries_leak_prof.py

Line #    Mem usage    Increment   Line Contents
================================================
    11     71.9 MiB     71.9 MiB   @profile
    12                             def np_func(x, func):
    13     71.9 MiB      0.0 MiB       try:
    14     73.2 MiB      1.3 MiB           r = func(x)
    15     73.2 MiB      0.0 MiB       except (TypeError, SystemError) as err:
    16     73.2 MiB      0.0 MiB           print(err)
    17     73.2 MiB      0.0 MiB       return


ufunc 'isnan' not supported for the input types, and the inputs could not be safely coerced to any supported types according to the casting rule ''safe''
Filename: timeseries_leak_prof.py

Line #    Mem usage    Increment   Line Contents
================================================
    11     74.1 MiB     74.1 MiB   @profile
    12                             def np_func(x, func):
    13     74.1 MiB      0.0 MiB       try:
    14    167.2 MiB     93.1 MiB           r = func(x)
    15    167.2 MiB      0.0 MiB       except (TypeError, SystemError) as err:
    16    167.2 MiB      0.0 MiB           print(err)
    17    167.2 MiB      0.0 MiB       return


ufunc 'isnan' not supported for the input types, and the inputs could not be safely coerced to any supported types according to the casting rule ''safe''
Filename: timeseries_leak_prof.py

Line #    Mem usage    Increment   Line Contents
================================================
    11    167.2 MiB    167.2 MiB   @profile
    12                             def np_func(x, func):
    13    167.2 MiB      0.0 MiB       try:
    14    185.8 MiB     18.6 MiB           r = func(x)
    15    185.8 MiB      0.0 MiB       except (TypeError, SystemError) as err:
    16    185.8 MiB      0.0 MiB           print(err)
    17    185.8 MiB      0.0 MiB       return


ufunc 'isnan' not supported for the input types, and the inputs could not be safely coerced to any supported types according to the casting rule ''safe''
Filename: timeseries_leak_prof.py

Line #    Mem usage    Increment   Line Contents
================================================
    11    185.8 MiB    185.8 MiB   @profile
    12                             def np_func(x, func):
    13    185.8 MiB      0.0 MiB       try:
    14    176.4 MiB      0.0 MiB           r = func(x)
    15    176.4 MiB      0.0 MiB       except (TypeError, SystemError) as err:
    16    176.4 MiB      0.0 MiB           print(err)
    17    176.4 MiB      0.0 MiB       return

> mprof run python3.8 timeseries_leak_prof.py mean                                                       
mprof: Sampling memory every 0.1s
1.19.1
4.1rc2
Unsupported operand type(s) for +: 'Time' and 'Time'
Filename: timeseries_leak_prof.py

Line #    Mem usage    Increment   Line Contents
================================================
    11     71.5 MiB     71.5 MiB   @profile
    12                             def np_func(x, func):
    13     71.5 MiB      0.0 MiB       try:
    14     77.8 MiB      6.3 MiB           r = func(x)
    15     77.8 MiB      0.0 MiB       except (TypeError, SystemError) as err:
    16     73.1 MiB      0.0 MiB           print(err)
    17     73.1 MiB      0.0 MiB       return


Unsupported operand type(s) for +: 'Time' and 'Time'
Filename: timeseries_leak_prof.py

Line #    Mem usage    Increment   Line Contents
================================================
    11     74.5 MiB     74.5 MiB   @profile
    12                             def np_func(x, func):
    13     74.5 MiB      0.0 MiB       try:
    14   2429.7 MiB   2355.2 MiB           r = func(x)
    15   2429.7 MiB      0.0 MiB       except (TypeError, SystemError) as err:
    16    176.6 MiB      0.0 MiB           print(err)
    17    176.6 MiB      0.0 MiB       return


Unsupported operand type(s) for +: 'Time' and 'Time'
Filename: timeseries_leak_prof.py

Line #    Mem usage    Increment   Line Contents
================================================
    11    176.6 MiB    176.6 MiB   @profile
    12                             def np_func(x, func):
    13    176.6 MiB      0.0 MiB       try:
    14   2430.2 MiB   2253.7 MiB           r = func(x)
    15   2430.2 MiB      0.0 MiB       except (TypeError, SystemError) as err:
    16    185.6 MiB      0.0 MiB           print(err)
    17    185.6 MiB      0.0 MiB       return


Unsupported operand type(s) for +: 'Time' and 'Time'
Filename: timeseries_leak_prof.py

Line #    Mem usage    Increment   Line Contents
================================================
    11    185.6 MiB    185.6 MiB   @profile
    12                             def np_func(x, func):
    13    185.6 MiB      0.0 MiB       try:
    14    176.9 MiB      0.0 MiB           r = func(x)
    15    176.9 MiB      0.0 MiB       except (TypeError, SystemError) as err:
    16    176.7 MiB      0.0 MiB           print(err)
    17    176.7 MiB      0.0 MiB       return

Note that

memory_profiler seems to add an additional 5-6 x timing overhead.
The mprof output does show the memory usage for np.mean, but not for ufuncs (because it's all allocated inside compiled C extensions?). However top shows similar usage for the latter, as do the attached plots from mprof plot.
Both for ufuncs and other functions, memory usage failed to trigger the debugger interface running e.g. python -m memory_profiler --pdb-mmem=200.
Each of func(ts.time) and func(copy.deepcopy(ts.time)) is responsible for 3 successive peaks in memory usage.
In a Python session memory generally is not immediately released after running functions like np.mean or np.median, but after a np.ufunc; with older Python or Numpy versions also interpreter segfaults have been observed.

Tested to the same effect with 4.0.1 and current master (including with #10656 merged).

System Details

See above for other Astropy versions tested, also with 4.1rc2 on Pythons 3.6.12, 3.7.9, 3.9.0rc1.

macOS-10.12.6-x86_64-i386-64bit
Python 3.8.5 (default, Jul 22 2020, 22:38:08) 
[Clang 9.0.0 (clang-900.0.39.2)]
Numpy 1.19.1
astropy 4.2.dev560+g679a74a46
Scipy 1.5.2
Matplotlib 3.3.0

timeseries_leak_prof.py.gz

@mhvk

This is clearly a question for @mhvk. Under the hood when this runs I can see that numpy is asking the Time object if it has __array_interface__, __array_struct__, and __array__ attributes. I would have thought it would be trying __array_function__, but this is where I am out of my depth since that is apparently not the case.

What does end up happening is that the __iter__ method in utils/shapes.py gets run and numpy converts the one Time object with 10000 elements to 10000 scalar Times, not good. I'm not sure at all what is going on with TimeSeries that makes this expensive (but doable in a second) operation turn into a showstopper.

@maxnoe

These are ufuncs, so it would try to access __array_ufunc__ -- which Time currently does not have (though @maxnoe is working on it over at #10272). As is, the ufunc creates an object array and then tries to iterate over it.

What I completely fail to understand is why the input time t behaves differently from ts.time. The latter should just be a copy of the former (both are Time instances).

But I don't understand TimeSeries well enough to know whether anything is changed in the instance. cc @astrofrog.

p.s. The difference in timing is also clear in the following in just doing np.array(t) vs np.array(ts.time)

EDIT: and np.array(QTable([t], names=['time'])['time']) is not as slow, so it is TimeSeries specific

EDIT: and np.array(QTable([t], names=['time'])['time']) is not as slow, so it is TimeSeries specific

Yes, I had tried to track it down as far as BaseTimeSeries...

The np.array(t) variation is interesting, as it shows the memory actually seems to be kept in use by the object array. E.g. the memory signature for

tt = copy.deepcopy(ts['time'])
t1 = np.array(tt)
time.sleep(30)
t2 = copy.deepcopy(t1)
time.sleep(30)
del t1
time.sleep(30)
del t2

looks like this

so the excess memory is apparently even copied over.

so it would try to access __array_ufunc__

I did a debug trace and from what I saw it was not looking for __array_ufunc__ but rather those other three I mentioned. It is certainly possible that I am somehow misinterpreting.

Strange, any ufunc should look for __array_ufunc__.

But the strange thing remains the link to TimeSeries - could it be because it is also the index column? It is probably the submodule I know least about.

@mhvk

¡Hola, ustedes! I found some insight that might be helpful.

In response to @mhvk

What I completely fail to understand is why the input time t behaves differently from ts.time. The latter should just be a copy of the former (both are Time instances).

So I wrote a script like so:

from astropy.time import Time; from astropy.timeseries import TimeSeries
import numpy as np
import time as stopwatch

def foo():
	t = Time(np.arange(10000), format='mjd', scale='tdb') # t for time
	
	ts = TimeSeries(time=t) 
	ts_t = ts.time # ts_t for timeseries time
	
	start = stopwatch.time()
	try: 
		[x for x in ts_t]
	except Exception as e:
		pass
	stop = stopwatch.time()
	
	start2 = stopwatch.time()
	try: 
		[x for x in t]
	except Exception as e:
		pass
	stop2 = stopwatch.time()	

	print("ts_t time: {}\nt time: {}".format(str(stop - start), str(stop2-start2)))

>> ts_t time: 7.85418701171875
>> t time: 1.178541660308838

Indeed it appears just iterating over timeseries.time takes longer than time.

I happened to follow the traceback for index.__deepcopy__ it goes something like shapes.__iter__ --> shapes._apply --> time.__apply
and in time._apply, there's this little line
tm.info = self.info which copies the indecies of the time object!

So why is timeseries.time different from time? like @mhvk said, one reason is the index.

print("ts_t.info.indecis:{}\n\n t.info.indecis: {}".format(str(ts_t.info.indices), str(t.info.indices)))
>> ts_t.info.indecis:[ time  rows
------ ----
   0.0    0
   1.0    1
   2.0    2
   3.0    3
   4.0    4
   5.0    5
   6.0    6
   ...  ...
9993.0 9993
9994.0 9994
9995.0 9995
9996.0 9996
9997.0 9997
9998.0 9998
9999.0 9999
Length = 10000 rows]

 t.info.indecis: None

Timeseries adds an index to the time table on init, specifically in the function Timeseries.add_coumn()
Indeed if I remove this line, and make Timeseries.time have no index, I get more similar (though not exact) execution speeds.

So I believe index__deep_copy__ is being called for each iteration of timeseries.time, thus slowing down its execution. The reason Time isn't as slow is because it has no index.

`

@iamsoto

@iamsoto - very nice sleuthing indeed!!

But less sure how to solve this... @taldcroft, any suggestions?

A general question: should the index live on the column rather than on the table? I think it is the only attribute that shares a shape with the host column, which is what makes it problematic - possibly, it needs its own _apply method? EDIT: more basic question: should the indices ever be passed on to the result of _apply?

Perhaps more simply (though at the risk of treating the symptom rather than the cause), since info can (I think) never be relevant for scalars, one could perhaps guard the copy of .info with a shape check, i.e.,

if tm.shape():
    tm.info = self.info

should the index live on the column rather than on the table?

A table index can be on multiple columns, so I think the answer here is no.

@iamsoto

one could perhaps guard the copy of .info with a shape check,

This seems like a sensible idea. So it would appear the memory problems are basically from making N deep copies of an N-long index, which for N=10000 gets bad.

So in time._apply we could do:

        if 'info' in self.__dict__ and tm.shape:
            tm.info = self.info

This of course has precedent from Column which drops info since col[1] just returns the bare numpy scalar or object at index 1.

I just tested this and it passes existing tests and in the functional tests above (the original one and one from @iamsoto), there is no difference in performance with or without index.

Who wants to do this PR? 😄

@iamsoto

@iamsoto - in the code you put in, it appears you are using tabs for indentation instead of 4 spaces. That's a pretty big no-no and will cause all kinds of problems if you try to work with most other code besides your own.

Oh No! sorry.

Confirm that this also resolves the memory blow-up for the numpy ufuncs. Still takes a long time compared to successfully calling the function on ts.time.value – I am wondering if there is no way to hand the exception once it is raised on the first element; but the runtime definitely is no longer a show-stopper.

I think solving the rest will require implementing __array_ufunc__ - could even do __array_ufunc__ = None to just immediately fail. But I suspect that will break code of people relying on things to work via the object-array route (even if very slowly).

I did the proposed fix in #10889.

Navigation Menu

Repository navigation

Excessive delay and memory usage when calling unsupported Numpy functions on TimeSeries.time #10688

Description

Expected behavior

Actual behavior

Steps to Reproduce

System Details

Metadata

Assignees

Labels

Type

Projects

Milestone

Relationships

Development

Excessive delay and memory usage when calling unsupported Numpy functions on TimeSeries.time #10688

Description

Description

Expected behavior

Actual behavior

Steps to Reproduce

System Details

Activity

taldcroft commented on Aug 30, 2020

mhvk commented on Aug 31, 2020

dhomeier commented on Aug 31, 2020

taldcroft commented on Aug 31, 2020

mhvk commented on Sep 1, 2020

i-am-b-soto commented on Oct 9, 2020

mhvk commented on Oct 10, 2020

taldcroft commented on Oct 12, 2020

taldcroft commented on Oct 12, 2020

taldcroft commented on Oct 12, 2020

i-am-b-soto commented on Oct 12, 2020

dhomeier commented on Oct 12, 2020

mhvk commented on Oct 13, 2020

taldcroft commented on Oct 17, 2020

Metadata

Metadata

Assignees

Labels

Type

Projects

Milestone

Relationships

Development

Issue actions