[Numpy-discussion] Interpolation question

Andrea Gavana andrea.gavana at gmail.com
Sun Mar 28 16:47:31 EDT 2010

HI All,

On 28 March 2010 19:22, Robert Kern wrote:
> On Sun, Mar 28, 2010 at 03:26, Anne Archibald <peridot.faceted at gmail.com> wrote:
>> On 27 March 2010 20:24, Andrea Gavana <andrea.gavana at gmail.com> wrote:
>>> Hi All,
>>>    I have an interpolation problem and I am having some difficulties
>>> in tackling it. I hope I can explain myself clearly enough.
>>> Basically, I have a whole bunch of 3D fluid flow simulations (close to
>>> 1000), and they are a result of different combinations of parameters.
>>> I was planning to use the Radial Basis Functions in scipy, but for the
>>> moment let's assume, to simplify things, that I am dealing only with
>>> one parameter (x). In 1000 simulations, this parameter x has 1000
>>> values, obviously. The problem is, the outcome of every single
>>> simulation is a vector of oil production over time (let's say 40
>>> values per simulation, one per year), and I would like to be able to
>>> interpolate my x parameter (1000 values) against all the simulations
>>> (1000x40) and get an approximating function that, given another x
>>> parameter (of size 1x1) will give me back an interpolated production
>>> profile (of size 1x40).
>> If I understand your problem correctly, you have a function taking one
>> value as input (or one 3D vector) and returning a vector of length 40.
>> You want to know whether there are tools in scipy to support this.
>> I'll say first that it's not strictly necessary for there to be: you
>> could always just build 40 different interpolators, one for each
>> component of the output. After all, there's no interaction in the
>> calculations between the output coordinates. This is of course
>> awkward, in that you'd like to just call F(x) and get back a vector of
>> length 40, but that can be remedied by writing a short wrapper
>> function that simply calls all 40 interpolators.
>> A problem that may be more serious is that the python loop over the 40
>> interpolators can be slow, while a C implementation would give
>> vector-valued results rapidly.
> 40 is not a bad number to loop over. The thing I would be worried
> about is the repeated calculation of the (1000, 1000) radial function
> evaluation. I think that a small modification of Rbf could be made to
> handle the vector-valued case. I leave that as an exercise to Andrea,
> of course. :-)

It seems like this whole interpolation stuff is not working as I
thought. In particular, considering scalar-valued interpolation (i.e.,
looking at the final oil recovery only and not the time-based oil
production profile), interpolation with RBFs is giving
counter-intuitive and meaningless answers. The issues I am seeing are
basically these:

# Interpolate the cumulative oil production
rbf = Rbf(x1, x2, x3, x4, x5, x6, final_oil_recovery)

# Try to use existing combination of parameters to get back
# the original result (more or less)
possible_recovery = rbf(x1[10], x2[10], x3[10], x4[10], x5[10], x6[10])

1) If I attempt to use the resulting interpolation function ("rbf"),
inputting a single value for each x1, x2, ..., x6 that is *already
present* in the original x1, x2, ..., x6 vectors, I get meaningless
results (i.e., negative oil productions, 1000% error, and so on) in
all cases with some (rare) exceptions when using the "thin-plate"
interpolation option;
2) Using a combination of parameters x1, x2, ..., x6 that is *not* in
the original set, I get non-physical (i.e., unrealistic) results: it
appears that RBFs think that if I increase the number of production
wells I get less oil recovered (!), sometimes by a factor of 10.

I am pretty sure I am missing something in this whole interpolation
thing (I am no expert at all), but honestly it is the first time I get
such counter-intuitive results in my entire Python numerical life ;-)


"Imagination Is The Only Weapon In The War Against Reality."

==> Never *EVER* use RemovalGroup for your house removal. You'll
regret it forever.
http://thedoomedcity.blogspot.com/2010/03/removal-group-nightmare.html <==

More information about the NumPy-Discussion mailing list