前提・実現したいこと
Pythonで微分方程式を解き、グラフに表示させたいです。
dy(i)/dt=-ay(i){b-exp(b-cx(i)**4)}+Σ[d*y(k)exp{(x(i)-x(k)/2)**2/(ex(k)**2)}]
Σの変数はkでx(i)<x(k)<3.0です。
縦軸に変数y、横軸に変数x(0<x<3)
時間は120sec,dt=6.0sec
Σ内のy(k)をどう扱ったらよいか分からず、自分でプログラムを組んでみたのですがエラーが出てしまいました。
Python初心者なので分かりやすく教えていただきたいです。
宜しくお願いします。
発生している問題・エラーメッセージ
--------------------------------------------------------------------------- NameError Traceback (most recent call last) <ipython-input-38-6e757769c90d> in <module> 55 t=np.arange(0,120,dt) 56 ---> 57 Phelps_Model(a,b,c,d,e,t,dt,x,y).calc() <ipython-input-38-6e757769c90d> in calc(self) 29 F_dy=lambda value: -self._a*value[1]*(self._b-math.exp(self._b-self._c*value[0]**4))+sigma(G(k), value[0], 3.0) 30 for j in range(0,len(self._t)): ---> 31 yy=F_dy([self._x,self._y])*self._dt 32 self._y+=yy 33 resultx.append(self._x) <ipython-input-38-6e757769c90d> in <lambda>(value) 27 resultx=[] 28 resulty=[] ---> 29 F_dy=lambda value: -self._a*value[1]*(self._b-math.exp(self._b-self._c*value[0]**4))+sigma(G(k), value[0], 3.0) 30 for j in range(0,len(self._t)): 31 yy=F_dy([self._x,self._y])*self._dt <ipython-input-32-4694edbacd1d> in G(k) 4 5 def G(k): ----> 6 return self._d*value[1]*math.exp((value[0]-k/2)**2)/self._e*k**2 7 k=3 8 NameError: name 'self' is not defined
該当のソースコード
import matplotlib.pyplot as plt import numpy as np import math class Phelps_Model(object): def __init__(self, a, b, c, d, e, t, dt, x, y): self._a=a self._b=b self._c=c self._d=d self._e=e self._t=t self._dt=dt self._x=x self._y=y def G(k): return self._d*value[1]*math.exp((value[0]-k/2)**2)/self._e*k**2 def sigma(func, frm, to): result = 0; for i in range(frm, to+0.1): result += func(i) return result def calc(self): resultx=[] resulty=[] F_dy=lambda value: -self._a*value[1]*(self._b-math.exp(self._b-self._c*value[0]**4))+sigma(G(k), value[0], 3.0) for j in range(0,len(self._t)): yy=F_dy([self._x,self._y])*self._dt self._y+=yy resultx.append(self._x) resulty.append(self._y) self.show(resultx,resulty) def show(self,x,y): fig=plt.figure() ax=fig.add_subplot() ax.plot(x,y,'o') plt.xlim(0, 4.0) plt.ylim(0, 100) ax.set_xlabel('fiber length') ax.set_ylabel('Number of fibers') plt.show() a=0.2 b=1.0 c=1.337 d=0.002 e=0.02 x=0.0 y=0.0 dt=6.0 t=np.arange(0,120,dt) Phelps_Model(a,b,c,d,e,t,dt,x,y).calc()
あなたの回答
tips
プレビュー