迭代 2d 网格插值与洞(缺失值)

2024-06-30 07:31:27 发布

您现在位置:Python中文网/ 问答频道 /正文

我有以下问题:

我有两个网格,s0Bs1B,都是20x20,还有一个函数{},在这些网格点上求值。在

但是,对于某些角点,s_0(s0, s1)无法求值,并返回np.nan。在

我想有一个插值函数,我评估网格点,然后我可以使用在一个根解决问题(不动点在s_0和模拟s_1)。下面是目前的代码草图:

# mesh 1d grids:
ss0 = np.array([ 0.38445455,  0.40388198,  0.42110923,  0.43626859,  0.44949237,
    0.46091287,  0.47066239,  0.47887323,  0.48567771,  0.49120811,
    0.49559675,  0.49897593,  0.50147794,  0.50323509,  0.50437969,
    0.50504404,  0.50536044,  0.50546118,  0.50547859,  0.50554495])
ss1 = np.array([ 0.43490909,  0.50764658,  0.5719651 ,  0.62838234,  0.67741596,
    0.71958365,  0.75540309,  0.78539195,  0.81006791,  0.82994865,
    0.84555185,  0.85739519,  0.86599635,  0.87187299,  0.87554281,
    0.87752348,  0.87833267,  0.87848808,  0.87850736,  0.87890821])

s0B, s1B = np.meshgrid(ss0, ss1, indexing='ij')
# generate interpolated functions
idx = ~np.isnan(s0Max)
if idx.sum() == 0:
    # check if there are any points at all
    raise notImplementedError
s0MaxInterp = interpolate.interp2d(
    s0B[idx], s1B[idx], s0Max[idx],
    fill_value=np.nan, kind='linear')
idx = ~np.isnan(s1Max)
s1MaxInterp = interpolate.interp2d(
    s0B[idx], s1B[idx], s1Max[idx],
    fill_value=np.nan, kind='linear')

def errorFixedPoint(x, s0MaxInterp, s1MaxInterp, grids):
    s0, s1 = x
    # I skip some checks whether s0, s1 are inside the interpolation range
    return np.array([s0 - s0MaxInterp(s0, s1)[0], s1 - s1MaxInterp(s0, s1)[0]])

result = optimize.root(errorFixedPoint, np.array([0.5, 0.9]),
                           args=(s0MaxInterp, s1MaxInterp, (ss0, ss1)), method='lm')

但是,当我尝试这个时,我得到了一个RunTimeWarning

^{pr2}$

插值有时看起来很离谱:

s0Max, grid and interpolated

然而,这可能是因为函数在两个维度上都是非常非线性的:

s0Max in 2d

我在其他地方读到建议使用scipy.interpolate.griddata,而不是{}。但是,我必须迭代地插值(对于根查找过程),而网格数据不支持这一点。在

我如何解决这个问题?在

为了重现性,我将s0Maxs1Max的值粘贴在这里:

nan = np.nan
s0Max = np.array([[        nan,         nan,         nan,         nan,         nan,
                nan,         nan,         nan,         nan,         nan,
         0.51494141,  0.51347656,  0.51210937,  0.51210937,  0.51210937,
         0.51210937,  0.51210937,  0.51210937,  0.51210937,  0.51210937],
       [        nan,         nan,         nan,         nan,         nan,
         0.51337891,  0.5109375 ,  0.51083984,  0.50253906,  0.50175781,
         0.50429687,  0.50957031,  0.50859375,  0.50849609,  0.50839844,
         0.50839844,  0.50839844,  0.50839844,  0.50839844,  0.50839844],
       [        nan,         nan,         nan,  0.51347656,  0.50527344,
         0.50009766,  0.49824219,  0.49746094,  0.49746094,  0.496875  ,
         0.49746094,  0.49941406,  0.50742187,  0.50732422,  0.50732422,
         0.50732422,  0.50732422,  0.50732422,  0.50732422,  0.50732422],
       [        nan,         nan,  0.51347656,  0.5015625 ,  0.49863281,
         0.49746094,  0.49628906,  0.49619141,  0.49570313,  0.49570313,
         0.49619141,  0.49746094,  0.50722656,  0.50634766,  0.50625   ,
         0.50625   ,  0.50625   ,  0.50625   ,  0.50625   ,  0.50625   ],
       [        nan,         nan,  0.50253906,  0.49824219,  0.496875  ,
         0.49619141,  0.49560547,  0.49511719,  0.49511719,  0.49511719,
         0.49550781,  0.49677734,  0.50625   ,  0.50625   ,  0.50625   ,
         0.50625   ,  0.50625   ,  0.50625   ,  0.50625   ,  0.50625   ],
       [        nan,  0.51210937,  0.49941406,  0.49746094,  0.49619141,
         0.49560547,  0.49511719,  0.49453125,  0.49453125,  0.49453125,
         0.49511719,  0.49628906,  0.50625   ,  0.50527344,  0.50527344,
         0.50527344,  0.50527344,  0.50527344,  0.50527344,  0.50527344],
       [        nan,  0.50341797,  0.49804688,  0.49628906,  0.49570313,
         0.49511719,  0.49453125,  0.49414062,  0.49404297,  0.49404297,
         0.49453125,  0.49570313,  0.50527344,  0.50527344,  0.50527344,
         0.50527344,  0.50527344,  0.50527344,  0.50527344,  0.50527344],
       [        nan,  0.50078125,  0.49746094,  0.49619141,  0.49511719,
         0.49453125,  0.49414062,  0.49404297,  0.49394531,  0.49394531,
         0.49453125,  0.49570313,  0.50527344,  0.50527344,  0.50527344,
         0.50527344,  0.50527344,  0.50527344,  0.50527344,  0.50527344],
       [        nan,  0.49941406,  0.496875  ,  0.49570313,  0.49511719,
         0.49453125,  0.49404297,  0.49394531,  0.49355469,  0.49394531,
         0.49414062,  0.49570313,  0.50527344,  0.50527344,  0.50517578,
         0.50517578,  0.50507812,  0.50507812,  0.50507812,  0.50449219],
       [        nan,  0.49873047,  0.49667969,  0.49560547,  0.49462891,
         0.49414062,  0.49394531,  0.49355469,  0.49345703,  0.49355469,
         0.49453125,  0.49628906,  0.50527344,  0.50517578,  0.50429687,
         0.50429687,  0.50429687,  0.50429687,  0.50429687,  0.50429687],
       [ 0.51347656,  0.49863281,  0.49628906,  0.49511719,  0.49453125,
         0.49404297,  0.49394531,  0.49355469,  0.49355469,  0.49404297,
         0.49570313,  0.50527344,  0.50429687,  0.50429687,  0.50429687,
         0.50429687,  0.50429687,  0.50429687,  0.50429687,  0.50429687],
       [ 0.51347656,  0.49814453,  0.49628906,  0.49511719,  0.49453125,
         0.49404297,  0.49394531,  0.49355469,  0.49404297,  0.49570313,
         0.50527344,  0.50429687,  0.50429687,  0.50429687,  0.50351563,
         0.50351563,  0.50341797,  0.50341797,  0.50341797,  0.50341797],
       [ 0.51347656,  0.49814453,  0.49619141,  0.49511719,  0.49453125,
         0.49404297,  0.49394531,  0.49404297,  0.49511719,  0.50527344,
         0.50429687,  0.50429687,  0.50351563,  0.50341797,  0.50341797,
         0.50332031,         nan,         nan,         nan,         nan],
       [ 0.51347656,  0.49814453,  0.49619141,  0.49511719,  0.49453125,
         0.49404297,  0.49404297,  0.49453125,  0.49873047,  0.50507812,
         0.50429687,  0.50351563,  0.50341797,         nan,         nan,
                nan,         nan,         nan,         nan,         nan],
       [ 0.51347656,  0.49814453,  0.49619141,  0.49511719,  0.49453125,
         0.49404297,  0.49404297,  0.49511719,  0.50527344,  0.50429687,
         0.50429687,  0.50341797,         nan,         nan,         nan,
                nan,         nan,         nan,         nan,         nan],
       [ 0.51347656,  0.49804688,  0.49619141,  0.49511719,  0.49453125,
         0.49404297,  0.49453125,  0.49570313,  0.50527344,  0.50429687,
         0.50429687,  0.50341797,         nan,         nan,         nan,
                nan,         nan,         nan,         nan,         nan],
       [ 0.51347656,  0.49804688,  0.49619141,  0.49511719,  0.49453125,
         0.49414062,  0.49453125,  0.49619141,  0.50527344,  0.50429687,
         0.50429687,  0.50341797,         nan,         nan,         nan,
                nan,         nan,         nan,         nan,         nan],
       [ 0.51347656,  0.49804688,  0.49619141,  0.49511719,  0.49453125,
         0.49414062,  0.49453125,  0.49628906,  0.50527344,  0.50429687,
         0.50351563,  0.50341797,         nan,         nan,         nan,
                nan,         nan,         nan,         nan,         nan],
       [ 0.51347656,  0.49804688,  0.49619141,  0.49511719,  0.49453125,
         0.49414062,  0.49453125,  0.49628906,  0.50527344,  0.50429687,
         0.50351563,  0.50341797,         nan,         nan,         nan,
                nan,         nan,         nan,         nan,         nan],
       [ 0.51347656,  0.49804688,  0.49619141,  0.49511719,  0.49453125,
         0.49414062,  0.49453125,  0.49628906,  0.50527344,  0.50429687,
         0.50351563,  0.50332031,         nan,         nan,         nan,
                nan,         nan,         nan,         nan,         nan]])

以及

s1Max = np.array([[        nan,         nan,         nan,         nan,         nan,
                nan,         nan,         nan,         nan,         nan,
         0.99894531,  0.99895996,  0.99895996,  0.99895996,  0.99895996,
         0.99895996,  0.99895996,  0.99895996,  0.99895996,  0.99895996],
       [        nan,         nan,         nan,         nan,         nan,
         0.99895996,  0.99895996,  0.99895996,  0.98726563,  0.98429688,
         0.99896118,  0.99895996,  0.99895996,  0.99895996,  0.99895996,
         0.99895996,  0.99895996,  0.99895996,  0.99895996,  0.99895996],
       [        nan,         nan,         nan,  0.99895996,  0.99896118,
         0.97      ,  0.95265625,  0.93984375,  0.93195313,  0.930625  ,
         0.93929688,  0.96421875,  0.99895996,  0.99895996,  0.99895996,
         0.99895996,  0.99895996,  0.99895996,  0.99895996,  0.99895996],
       [        nan,         nan,  0.99895996,  0.98117188,  0.95429688,
         0.93351563,  0.91695313,  0.90453125,  0.89703125,  0.896875  ,
         0.908125  ,  0.93789063,  0.99895996,  0.99895996,  0.99895996,
         0.99895996,  0.99895996,  0.99895996,  0.99895996,  0.99895996],
       [        nan,         nan,  0.98648438,  0.95296875,  0.9278125 ,
         0.9075    ,  0.89101563,  0.87859375,  0.8715625 ,  0.8725    ,
         0.88617188,  0.92085938,  0.99895996,  0.99895996,  0.99895996,
         0.99895996,  0.99895996,  0.99895996,  0.99895996,  0.99895996],
       [        nan,  0.99895996,  0.9640625 ,  0.93265625,  0.90789063,
         0.88757813,  0.87109375,  0.85875   ,  0.85210938,  0.8540625 ,
         0.86992188,  0.90914063,  0.99895996,  0.99895996,  0.99895996,
         0.99895996,  0.99895996,  0.99895996,  0.99895996,  0.99895996],
       [        nan,  0.99023438,  0.94773438,  0.91703125,  0.89234375,
         0.87203125,  0.85546875,  0.84320313,  0.83695313,  0.83976563,
         0.85757813,  0.9009375 ,  0.99895996,  0.99895996,  0.99895996,
         0.99895996,  0.99895996,  0.99895996,  0.99895996,  0.99895996],
       [        nan,  0.97570313,  0.93515625,  0.90476563,  0.88015625,
         0.8596875 ,  0.843125  ,  0.8309375 ,  0.825     ,  0.82867188,
         0.84820313,  0.8953125 ,  0.99895996,  0.99895996,  0.99895996,
         0.99895996,  0.99895996,  0.99895996,  0.99895996,  0.99895996],
       [        nan,  0.96515625,  0.92546875,  0.89515625,  0.87046875,
         0.85      ,  0.83335938,  0.82140625,  0.81585938,  0.82046875,
         0.84203125,  0.89390625,  0.99895996,  0.99895996,  0.99895996,
         0.99895996,  0.99895996,  0.99895996,  0.99895996,  0.99895996],
       [        nan,  0.95757813,  0.918125  ,  0.88773438,  0.86296875,
         0.8425    ,  0.82609375,  0.8146875 ,  0.81054688,  0.81820313,
         0.84664063,  0.91648438,  0.99895996,  0.99895996,  0.99895996,
         0.99895996,  0.99895996,  0.99895996,  0.99895996,  0.99895996],
       [ 0.99896484,  0.953125  ,  0.913125  ,  0.88242188,  0.8575    ,
         0.83726563,  0.821875  ,  0.81289062,  0.81414063,  0.83429688,
         0.89601563,  0.99895996,  0.99895996,  0.99895996,  0.99895996,
         0.99895996,  0.99895996,  0.99895996,  0.99895996,  0.99895996],
       [ 0.99896484,  0.95078125,  0.90992188,  0.87867188,  0.85390625,
         0.83453125,  0.82164063,  0.81859375,  0.83421875,  0.89546875,
         0.99895996,  0.99895996,  0.99895996,  0.99896484,  0.99125   ,
         0.9909375 ,  0.99078125,  0.99078125,  0.99078125,  0.99070313],
       [ 0.99896484,  0.94945313,  0.9078125 ,  0.87632813,  0.85195313,
         0.83429688,  0.82554688,  0.83273438,  0.87804688,  0.99895996,
         0.99895996,  0.99895996,  0.99117188,  0.99023438,  0.98945313,
         0.98890625,         nan,         nan,         nan,         nan],
       [ 0.99896484,  0.94859375,  0.90640625,  0.875     ,  0.85125   ,
         0.835625  ,  0.83203125,  0.8534375 ,  0.9571875 ,  0.99895996,
         0.99895996,  0.99125   ,  0.98992188,         nan,         nan,
                nan,         nan,         nan,         nan,         nan],
       [ 0.99896484,  0.948125  ,  0.905625  ,  0.87429688,  0.85125   ,
         0.83757813,  0.83914063,  0.87632813,  0.99895996,  0.99895996,
         0.99896484,  0.99046875,         nan,         nan,         nan,
                nan,         nan,         nan,         nan,         nan],
       [ 0.99896484,  0.94789063,  0.90523438,  0.87398438,  0.85148438,
         0.83929688,  0.8446875 ,  0.89546875,  0.99895996,  0.99895996,
         0.99896484,  0.98992188,         nan,         nan,         nan,
                nan,         nan,         nan,         nan,         nan],
       [ 0.99896484,  0.94773438,  0.905     ,  0.87382813,  0.85164063,
         0.84023438,  0.84789063,  0.906875  ,  0.99895996,  0.99895996,
         0.99140625,  0.98960938,         nan,         nan,         nan,
                nan,         nan,         nan,         nan,         nan],
       [ 0.99896484,  0.94773438,  0.90492188,  0.87382813,  0.85171875,
         0.84054688,  0.84890625,  0.9109375 ,  0.99895996,  0.99895996,
         0.99132813,  0.98953125,         nan,         nan,         nan,
                nan,         nan,         nan,         nan,         nan],
       [ 0.99896484,  0.94773438,  0.90492188,  0.87375   ,  0.85171875,
         0.840625  ,  0.84914063,  0.91164063,  0.99895996,  0.99895996,
         0.99132813,  0.98945313,         nan,         nan,         nan,
                nan,         nan,         nan,         nan,         nan],
       [ 0.99896484,  0.94765625,  0.90492188,  0.87375   ,  0.85179688,
         0.84085938,  0.84984375,  0.91445313,  0.99895996,  0.99895996,
         0.99132813,  0.989375  ,         nan,         nan,         nan,
                nan,         nan,         nan,         nan,         nan]])

Tags: 函数网格npnanarray插值idxs1
1条回答
网友
1楼 · 发布于 2024-06-30 07:31:27

根据这个极好的答案here,使用interp2d将不适用于更复杂的数据。建议使用radial basis function interpolation。在

这似乎是可行的,并给出更好的插值结果

import scipy.optimize as optimize
from scipy.interpolate import Rbf

s0MaxInterp = Rbf(s0B[idx], s1B[idx], s0Max[idx], function='cubic')
s1MaxInterp = Rbf(s0B[idx], s1B[idx], s1Max[idx], function='cubic')

def errorFixedPoint(x, s0MaxInterp, s1MaxInterp, grids):
    s0, s1 = x
    # I skip some checks whether s0, s1 are inside the interpolation range
    return np.array([s0 - s0MaxInterp(s0, s1), s1 - s1MaxInterp(s0, s1)])

result = optimize.root(errorFixedPoint, np.array([0.5, 0.9]),
                       args=(s0MaxInterp, s1MaxInterp, (ss0, ss1)), method='lm')

注意,我删除了errorFixedPoint函数中的[0]。在

相关问题 更多 >