Investigate returning np arrays as arguments

[MaxvdKolk]

Either we find a way that returning the arguments writes their result into the numpy arrays, similar to x[:] = ..., or pass the arrays into the routine as part of the arguments (oldschool style) and update them therein.

Any alternatives might work too, e.g. storing attributes in the classes, maybe using setter/getters.

[MaxvdKolk]

Something I noticed, the numpy functions seem to provide an additional optional argument out=None that allows to pass an output array for the routine, e.g. np.dot.

Output argument. This must have the exact kind that would be returned if it was not used. In particular, it must have the right type, must be C-contiguous, and its dtype must be the dtype that would be returned for dot(a,b). This is a performance feature. Therefore, if these conditions are not met, an exception is raised, instead of attempting to be flexible.

We could look in their source see how the implement these checks and adopts something similar for the routines where we noticed many new np.empty constructions.

Having the optional argument only makes sense if we want to have both options, otherwise we might just require that the output is part of the required set of arguments.

[artofscience]

@Giannis1993 can you remove the np.empty calls (store data in attributes or some other way) and run the profile again. What was the overhead of these calls?

[aatmdelissen]

@MaxvdKolk @Giannis1993 I guess it should look like this:

def g(self, x, out=None):
    if out is None:
        out = np.empty(...)
    ...
    return out

[aatmdelissen]

Update: This functionality has been added in #57. The next thing to do is to use these in the function calls of the subsolvers (@artofscience).

[artofscience]

@aatmdelissen use them as: f = zeros(n,1) g(x, out=f) ?

[aatmdelissen]

@artofscience yes :) can be used as

f = g(x)  # Create new array
g(x, out=f)  # Reuse the array
f = g(x, out=f)  # Also reuse the array, and get reference back

[artofscience]

Finished