0

我想通过 fmin 优化来解决定点迭代。我不知道为什么我得到一个错误，如：OVERFLOW ENCOUNTERS IN DOUBLE SCALAR: I know the results and the first optimization before printing Status: 0 ([ 2.85475105e-01 3.48080784e-01 -2.54628055e-04]) 是正确的并且第二个值对也是。在 3. 优化循环出现错误..

我只是不知道出了什么问题感谢您的帮助。

import numpy as np
import os
from scipy.optimize import fmin



dt = 0.001                  
t0 = 0                                                                     
tau = 4.0                                                                  
Steps = int(tau / dt)  + 1

def DGLS(t, Y0, Y1, Y2, Y3, Y4):

 return [(-Y0+alpha-Y0*Y4*(Y1/(1+mu2*Y2)+c)-phi*Y0*Y4*(Y2+c)-chi*Y0*Y4*(Y3+c)) ,
         (-Y1+nu*Y0*Y4/(1+mur*Y3)*(Y1/(1+mu2*Y2)+c)) ,
         (-Y2+phi*nu*Y0*Y4/(1+mur*Y3)*(Y2+c)/(1+mu1*Y1/(1+mu2*Y2))) ,
         (-Y3+chi*nu*Y0*Y4*(Y3+c)) ,
         (-Y4*(Y1+Y2+Y3))]


def RK4_end(DGLS, Startwert, Steps, t0, dt):

'''

   Runge Kutta Verfahren 4. Ordnung mit der Stepsize dt 
   Wiedergabe des letzten Wertes 

   Startwert ... im Array Startwert sind sind die Anfangswerte gespeichert
   Steps     ... Steps ist eine Zahl die der Iterationslänge entspricht um bis zum Maximalen Injektionstag in time zu rechnen

''' 

Y = np.zeros(5)

for t in range(t0 , Steps):
    if t == t0 :
        Y[0], Y[1], Y[2], Y[3], Y[4] = Startwert[0], Startwert[1], Startwert[2], Startwert[3], Startwert[4]

    k1 = DGLS( t , Y[0] , Y[1] , Y[2] , Y[3] , Y[4] )
    k2 = DGLS( t + 0.5 * dt , Y[0] + 0.5 * dt * k1[0] , Y[1] + 0.5 * dt * k1[1] , Y[2] + 0.5 * dt * k1[2] , Y[3] + 0.5 * dt * k1[3] , Y[4] + 0.5 * dt * k1[4] )
    k3 = DGLS( t + 0.5 * dt , Y[0] + 0.5 * dt * k2[0] , Y[1] + 0.5 * dt * k2[1] , Y[2] + 0.5 * dt * k2[2] , Y[3] + 0.5 * dt * k2[3] , Y[4] + 0.5 * dt * k2[4] )
    k4 = DGLS( t + dt , Y[0] + dt * k3[0] , Y[1] + dt * k3[1] , Y[2] + dt * k3[2] , Y[3] + dt * k3[3] , Y[4] + dt * k3[4] )


    Y[0] = Y[0] + (dt / 6.)*( k1[0] + 2.0 * k2[0] + 2.0 * k3[0] + k4[0] )
    Y[1] = Y[1] + (dt / 6.)*( k1[1] + 2.0 * k2[1] + 2.0 * k3[1] + k4[1] )
    Y[2] = Y[2] + (dt / 6.)*( k1[2] + 2.0 * k2[2] + 2.0 * k3[2] + k4[2] )
    Y[3] = Y[3] + (dt / 6.)*( k1[3] + 2.0 * k2[3] + 2.0 * k3[3] + k4[3] ) 
    Y[4] = Y[4] + (dt / 6.)*( k1[4] + 2.0 * k2[4] + 2.0 * k3[4] + k4[4] )

return Y
def calc(y):

  Startwert = np.array([10., y[0], y[1], y[2], 1.])

  solve = RK4_end(DGLS, Startwert, Steps, t0, dt)
  #print(solve[1],solve[2],solve[3],solve[4])
  return (solve[1]-y[0])**2 + (solve[2]-y[1])**2 + (solve[3]-y[2])**2


def main():

 accuracy = 1.0

 t1 = -5. * accuracy
 t2 = -5. * accuracy
 tr = -5. * accuracy

 T1v = []
 T2v = []
 Trv = []    

 T1 = []
 T2 = []
 Tr = []

 y = np.zeros(3)

 Status = 0
 while t1 < 0.1 * accuracy:

    while t2 < 0.1 * accuracy:

        while tr < 0.1 * accuracy:

            y[0] = 10**( t1 / accuracy )
            y[1] = 10**( t2 / accuracy )
            y[2] = 10**( tr / accuracy )

            simplex = fmin(calc, y , ftol=1E-15)
            print(simplex)

            if simplex[0]>0 and simplex[1]>0 and simplex[2]>0 and calc[y]<0.01:
                T1v.append(simplex[0])
                T2v.append(simplex[1])
                Trv.append(simplex[2])

            print(Status)
            Status = Status + 1

            tr = tr + 1

        tr = -5 * accuracy
        t2 = t2 + 1

    t2 = -5 * accuracy
    t1 = t1 + 1
4

0 回答 0