Минимизиране на leastsq с деление на нула

Бих искал да приложа минимизиране на функция, която съдържа деление и някаква тригонометрична функция, както можете да видите по-долу с част от кода. Но се появява грешка при разделяне на нула:

RuntimeWarning: среща се деление на нула при деление launch_new_instance()

Lx=2592.
Ly=1936.

YA, XA = np.mgrid[0:Ly, 0:Lx]

XAvect=np.reshape(XA,(Lx*Ly))
YAvect=np.reshape(YA,(Lx*Ly))

#Initialization
x0 = 200.
y0 = 200.
KI = 100000.
T = 20.
A1 = 50.
A2 = 50.


def residual_V2west(vars, XA, YA, donnees):
    KI = vars[0]
    A1 = vars[1]
    A2 = vars[2]
    x0 = vars[3]
    y0 = vars[4]
    modeleV2 = KI/(G)*np.sqrt(np.sqrt((XA-x0)**2 + (YA-y0)**2)/(2*np.pi))*np.sin((np.arctan(((YA-y0)/(XA-x0))))/2)*(Kappa + 1 - 2*np.cos(np.arctan(((YA-y0)/(XA-x0)))/2)**2) + \
    A1*np.sqrt((XA-x0)**2+(YA-y0)**2)*np.cos(np.arctan((YA-y0)/(XA-x0))) + A2 
    return (donnees-modeleV2)

from scipy.optimize import leastsq
vars = [KI, A1, A2, x0, y0]
out_V_west = sco.leastsq(residual_V2west, vars, args=(XAvect, YAvect, Vmvect),epsfcn=0.01)
print out_V_west

numpy divide-by-zero minimization
person user3601754    schedule 27.08.2014    source източник
comment
Опитвам с np.seterr(divide='ignore', invalid='ignore'), но изглежда не работи...   -  person user3601754    schedule 27.08.2014
comment
намерихте ли вече решението? ако да, можете да публикувате като отговор...   -  person Saullo G. P. Castro    schedule 27.08.2014
comment
В крайна сметка не работи във всички случаи...не разбирам...   -  person user3601754    schedule 27.08.2014

Отговори (2)


arrow_upward
1
arrow_downward

Можете да опитате да избегнете делението на нула, като дефинирате знаменателя да бъде малко число вместо нула:

XAmx0 = XA-x0 if abs(XAmx0) < 1e-12: XAmx0 = 1e-12

И след това разделете на XAmx0 вместо на XA-x0

Това регулира израза ви в остатъчната функция и не влияе на резултата ви (ако знаменателят в оптималното решение е от порядъка на 1e-12, вероятно ще имате проблеми с числената точност така или иначе)

person skulda    schedule 29.08.2014

arrow_upward
0
arrow_downward

Постигнах това с помощта на: ma.masked_invalid(...)

person user3601754    schedule 27.08.2014