проверка близости двух сплайнов в питоне - PullRequest
Новая
       8

проверка близости двух сплайнов в питоне

0 голосов
/ 28 июня 2018

У меня есть необходимость проверить, достаточно ли близки два сплайна с определенным допуском друг к другу. первый сплайн: 

spl_psi=scipy.interpolate.UnivariateSpline(t,psi,k=4,s=0)

- это «ввод данных», используемый для проверки второго, а именно: 

def fr1(x,y,z):
   return z

def fr2(t,r1,r2):
   return ((K_1*spl_delta(t)+(K_1*T_3)*spl_delta1(t)-(T_1*T_2)*r2-r1)/(T_1+T_2))

for i in range(0,len(T2)-1,1):  
    k1=h*fr1(T2[i],r1[i],r2[i])
    l1=h*fr2(T2[i],r1[i],r2[i])
    k2=h*fr1(T2[i]+0.5*h,r1[i]+0.5*k1,r2[i]+0.5*l1)
    l2=h*fr2(T2[i]+0.5*h,r1[i]+0.5*k1,r2[i]+0.5*l1)
    k3=h*fr1(T2[i]+0.5*h,r1[i]+0.5*k2,r2[i]+0.5*l2)
    l3=h*fr2(T2[i]+0.5*h,r1[i]+0.5*k2,r2[i]+0.5*l2)
    k4=h*fr1(T2[i]+h,r1[i]+k3,r2[i]+l3)
    l4=h*fr2(T2[i]+h,r1[i]+k3,r2[i]+l3)
    r1[i+1]=r1[i]+(k1+2*k2+2*k3+k4)/6
    r2[i+1]=r2[i]+(l1+2*l2+2*l3+l4)/6

spl_r1_n2=scipy.interpolate.UnivariateSpline(T2,r1,k=5,s=0) 
spl_psi_n2=scipy.interpolate.UnivariateSpline.antiderivative(spl_r1_n2)

spl_psi_n2 - это мой второй сплайн, который должен быть близок к spl_psi, а T_1, T_2, T_3, K_1 - мои параметры. Этот метод называется Рунге-Кутта, и он используется для решения дифференциального уравнения, такого как мой. Теперь я хочу найти параметры, которые делают лучший сплайн, а затем сохранить их в массиве с именем «res». Это мой код: 

def fr1(x,y,z):
return z

def fr2(x,y,z):
return ((K*spl_delta(x)+(K*(t_3/4))*spl_delta1(x)-((t_1/4)+(t_2/4))*z-y)/(t_1*t_2/(16) ))   

r1_0=spl_psi1(0)
r2_0=spl_psi2(0)
h=T2[2]-T2[1]
r1=np.array( [r1_0]*len(T2))
r2=np.array( [r2_0]*len(T2))
toll=0.5
res=[]

for t_1 in range (2,200,1):
  for t_2 in range (2,50,1):
    for t_3 in range (2,50,1):
        for i in range(0,len(T2)-1,1):
            k1=h*fr1(T2[i],r1[i],r2[i])
            l1=h*fr2(T2[i],r1[i],r2[i])
            k2=h*fr1(T2[i]+0.5*h,r1[i]+0.5*k1,r2[i]+0.5*l1)
            l2=h*fr2(T2[i]+0.5*h,r1[i]+0.5*k1,r2[i]+0.5*l1)
            k3=h*fr1(T2[i]+0.5*h,r1[i]+0.5*k2,r2[i]+0.5*l2)
            l3=h*fr2(T2[i]+0.5*h,r1[i]+0.5*k2,r2[i]+0.5*l2)
            k4=h*fr1(T2[i]+h,r1[i]+k3,r2[i]+l3)
            l4=h*fr2(T2[i]+h,r1[i]+k3,r2[i]+l3)
            r1[i+1]=r1[i]+(k1+2*k2+2*k3+k4)/6
            r2[i+1]=r2[i]+(l1+2*l2+2*l3+l4)/6
        spl_r1_n2=scipy.interpolate.UnivariateSpline(T2,r1,k=5,s=0)
        spl_psi_n2=scipy.interpolate.UnivariateSpline.antiderivative(spl_r1_n2)
        print ([t_1,t_2,t_3])
        for i in range (0,len(T2)-1,1):
            if fabs(spl_psi_n2(T2[i])-spl_psi(T2[i]))<toll:
                print 'ok'
                continue
            elif fabs(spl_psi_n2(T2[i])-spl_psi(T2[i]))>toll:
                print 'nope',i
                break
            elif i==len(T2)-2:
                print 'it's done',t_1,t_2,t_3
                res.append([t_1,t_2,t_3])
                break

Теперь, моя проблема в том, что когда я запускаю это, он дает мне «Все готово» для каждой комбинации t_1, t_2, t_3. Я совершенно уверен, что проблема в куске: 

for i in range (0,len(T2)-1,1):
    if fabs(spl_psi_n2(T2[i])-spl_psi(T2[i]))<toll:
        print 'ok'
        continue
    elif fabs(spl_psi_n2(T2[i])-spl_psi(T2[i]))>toll:
        print 'nope',i
        break
    elif i==len(T2)-2:
        print 'it's done',t_1,t_2,t_3
        res.append([t_1,t_2,t_3])
        break

но не могу обойти меня. Любая помощь?

...