Попытка оптимизировать параметры в модели динамического трения Лугре - PullRequest
Новая
       46

Попытка оптимизировать параметры в модели динамического трения Лугре

1 голос
/ 06 мая 2019

У меня есть данные, собранные в CSV каждого выхода модели трения.модель представляет контакт между поверхностями в виде одномерных щетинок, которые реагируют на изгиб, как пружины при этом отклонении.сила трения моделируется следующим образом: 

FL(V,Z) = sig0*Z +sig1*DZ/Dt +sig2*V

где V - скорость поверхности Z - отклонение щетинок, а DZ / Dt - скорость отклонения, равная: 

DZ/Dt = V + abs(V)*Z/(Fc + (Fs-Fc)*exp(-(V^2/Vs^2))
      = V + abs(V)*Z/G(V)
      = V + H(V)*Z

Где Fc - трение объекта в движении (постоянное), Fs равно Силе, необходимой для приведения объекта в движение (постоянная> Fc), а Vs - общая скорость, необходимая для перехода между доменами.(константа, которую я получил экспериментально).скорость и положение блока указаны в CSV, а также сила трения во времени.Я также создал легко интегрируемое приближение Скорости как функции времени (тригонометрическое).

В связи с проблемой: код не соответствует тому, как я пытаюсь передать списки в функции (я думаю).

Функция передает параметры SEEMS для работы(взят из другого файла, который просто отображает данные), однако я попытался численно интегрировать DZ / Dt и подогнать параметры sig к импортированным данным трения.

Что я импортировал 

import numpy as np
from scipy.integrate import odeint
import matplotlib.pyplot as plt
from scipy import optimize
import pylab as pp
from math import sin, pi, exp, fabs, pow

Параметры 

Fc=2.7  #N
Fs=8.2  #N
Vs=.34  #mm/s

Initial_conditions 

ITime=Time[0]
Iz=[0,0,0]

Построение модели трения 

def velocity(time):
    V=-13/2*1/60*pi*sin(1/60*pi*time+pi)       
    return V

def g(v,vs,fc,fs,sig0):
    G=(1/sig0)*(fc+(fs-fc)*exp(-pow(v,2)/pow(vs,2)))
    return G

def h(v,vg):
  H=fabs(v)/vg
  return H

def findz(z, time, sig):
  Vx=velocity(time)
  VG=g(Vx,Vs,Fc,Fs,sig)
  HVx=h(Vx,VG)
  dzdt=Vx+HVx*z
  return dzdt

def friction(time,sig,iz):
  dz=lambda z,time: findz(z,time,sig)
  z=odeint(dz,iz,time)
  return sig[0]*z+sig[1]*findz(z,time,sig[0])+sig[2]*velocity(Time)

Должно возвращать разницу между построенной функцией и данными и выходомсписок, содержащий оптимизированные параметры 

def residual(sig):
  return Ff-friction(Time,sig,Iz)

SigG=[4,20,1]
SigVal=optimize.leastsq(residual,SigG)

print "parameter values are ",SigVal

Возвращает 

line 56, in velocity
    V=-13/2*1/60*pi*sin(1/60*pi*time+pi)

TypeError: can't multiply sequence by non-int of type 'float'

Это связано с тем, что я передаю списки?

1 Ответ

0 голосов
/ 08 мая 2019

Как я уже упоминал в своем комментарии, Velocity() является причиной ошибки, которая, скорее всего, связана с тем, что она использует значение времени, тогда как вы передаете весь список / массив (с несколькими значениями) в Velocity() когда вы звоните в friction().

Используя некоторые выбранные значения и после сокращения вашего кода и передачи ITime вместо Time код работает правильно, но вам остается только судить, является ли это аналитически тем, чего вы хотели достичь. Ниже мой код: 

import numpy as np
from scipy import optimize
from scipy.integrate import odeint
from math import sin,  pi,  exp,  fabs

# Parameters
Fc = 2.7  #N
Fs = 8.2  #N
Vs = .34  #mm/s

# define test values for Ff and Time
Ff    = np.array([100, 50, 50])
Time  = np.array([10, 20, 30])

# Initial_conditions
ITime = Time[0]
Iz    = np.array([0, 0, 0])

# Building the friction model
V    = lambda                   t: (-13 / 2) * ( 1 / (60 * pi * sin(1 / 60 * pi * t + pi)))      
G    = lambda v, vs, fc, fs, sig0: (1 / sig0) * (fc + (fs - fc) * exp(-v**2 / vs**2))
H    = lambda               v, vg: fabs(v) / vg
dzdt = lambda         z,  t,  sig: V(t) + H(V(t), G(V(t), Vs, Fc, Fs, sig)) * z


def friction(t, sig, iz):
  dz = lambda z, t: dzdt(z, t, sig)
  z  = odeint(dz, iz, t)
  return sig[0]*z + sig[1]*dzdt(z, t, sig[0]) + sig[2]*V(t)

# Should return the difference between the Constructed function and the data
# and yield a list containing the optimized parameters

def residual(sig):

  return Ff-friction(ITime, sig, Iz)[0]

SigG   = np.array([4, 20, 1])
SigVal = optimize.leastsq(residual, SigG, full_output = False)

print("parameter values are ", SigVal )

Выход: 

parameter values are  (array([    4.        ,  3251.47271228, -2284.82881887]), 1)
...