首页 > 解决方案 > 将一个 numpy 数组传递给一个 numpy 函数数组

问题描述

我正在编写一个旨在解决微分方程组的类,该类是使用与另一个名为 rhs (右侧)的类的组合构成的,该类应该代表给定微分问题的函数。在 ode 系统的情况下,这个 rhs.function 变成了一个函数数组!(出于这个原因,我也使用向量作为解决方案)我遇到的问题是如何将解决方案的数组传递给函数,我将在这里显示代码,Rhs 类是:

 class Rhs:
   '''
      class to define a differential problems 
      contains the dy/dt function and the initial value 
      in order to close the ode problem
   '''
   solution = None

   def __init__(self, fnum : np.ndarray , t0: np.float, tf: np.float, y0 : np.array, n: int , fanal = None ):
      '''   
         Input : 
            - fnum : Function f(t,y(t)) = dy/dt (numeric)
            - 
            - Initial Time t0 
            - Final Time tf
            - Initial value u0_i(t0) = y0_i(t0)  
      '''
      self.func = fnum
      Rhs.solution  = fanal
      self.t0   = t0
      self.tf   = tf
      self.n    = n
      self.u0   = y0


   def createArray(self):
      '''   
         Create the Array time and f(time) 
         - the time array can be create now 
         - the solution array can be just initialize with the IV (array)
      ''' 
      self.t = np.linspace(self.t0, self.tf, self.n )
      self.u = np.array([self.u0  for i in range(self.n) ])


      return self.t,self.u   

   def f(self,ti,ui):
       return  np.array([function(ti,ui) for function in self.func])  

这个类在主文件中被实例化,(problem = Rhs(... args ...))然后传递给初始化类 Euler 在这里报告:

 class Explicit:
   '''
      Solve the ODE using first order Foward difference
      Explicit first order Method (Euler)  
   '''
   solved = False  # True if the method 'solve' was called without error 
   def __init__(self, dydt: rhs.Rhs, save : bool=True, _print: bool=False,  filename : str=None):
      ''' 
         Initialize the Foward Euler solver 
            - dydt (to the super class) : RHS problem dy/dt = f(t,y(t)) t(0)=y0 
            - save : if True returns the 2 vector time and u = du/dt 
      '''
      self._print = _print
      self.save   = save 
      self.file   = filename
      self.dydt   = dydt
      self.dt     = (dydt.tf-dydt.t0)/dydt.n



   def solve(self):
      self.time, self.u = self.dydt.createArray() 
      for i in range(len(self.time)-1):
          #for j in range(len(self.u[0,:])):
              self.u[i+1] = self.u[i] + self.dt*self.dydt.f(self.time[i],self.u[i])


      Explicit.solved = True 

      if self._print:
          with open(self.file,'w') as f:
              for i in range(len(self.time)):      
                 f.write('%.4f %4f %4f %4f\n' %(self.time[i] ,self.u[i,0], self.u[i,1], self.u[i,2]))

      if self.save:
         return self.time,self.u

最后是简单的主要功能:

func1 = lambda t,u : 10   * (u[1] - u[0])
func2 = lambda t,u : 28   * u[0] - u[1] - u[0] * u[2]  
func3 = lambda t,u : -8/3 * u[2] + u[0]*u[1]
def main():

   func = np.array([func1,func2,func3])

   y0 = np.array([1.,0.,0.]) 

   problem3   = rhs.Rhs(func,0.0,100.0,y0,1000) 
   t,u       = problem3.createArray()

   fwdeuler_p1 = euler.Explicit(problem3 , True,_print=True, filename ='lorentz.dat')
   fet,feu = fwdeuler_p1.solve()

我希望 lambda 函数内部使用的语法很清楚:每个都代表一个函数,其中 u[N] 是未知函数。

我在 C++ 中开发了一个类层次结构并且工作正常(如果你有兴趣可以在这里找到它)。

我用 python 程序得到的错误是这样的:

 drive.py:31: RuntimeWarning: overflow encountered in double_scalars
  func2 = lambda t,u : 28   * u[0] - u[1] - u[0] * u[2]
drive.py:32: RuntimeWarning: overflow encountered in double_scalars
  func3 = lambda t,u : -8/3 * u[2] + u[0]*u[1]
drive.py:31: RuntimeWarning: invalid value encountered in double_scalars
  func2 = lambda t,u : 28   * u[0] - u[1] - u[0] * u[2]
/home/marco/Programming/Python/Numeric/OdeSystem/euler.py:77: RuntimeWarning: invalid value encountered in add
  self.u[i+1] = self.u[i] + self.dt*self.dydt.f(self.time[i],self.u[i])

编辑真的很抱歉!纠正现在的错误是正确的。

编辑我解决了使时间步更小:) 感谢所有试图帮助我的人

标签: pythonarraysnumpy

解决方案


推荐阅读