2016-06-01 1 views
1

Я пытаюсь использовать fsolve но он возвращает результат от вызова функции не является собственно массив поплавки здесь кодРезультата от вызова функции не является надлежащим массив поплавков, fsolve, SymPy

Если я сделайте это так, это не сработает.

from sympy import * 
    B=symbols('B') 
    L=0.97 
    a = 0.04 
    b= 0.029 
    m = 4.8 
    I = 3.14*(a/2)*(b/2)**3/4 
    E = 10**10 
    ro = 900 
    A = 3.14*(a/2)*(b/2) 
    c = sqrt(E*I/(ro*A)) 
    AA=Matrix([[0,-B**3,0,B**3],[-B**2,0,B**2,0],[B**4*c**2*m*(cos(B*L)) + E*I*(B**3*sin(B*L)), B**4*c**2*m*(sin(B*L))+ E*I*(- B**3*cos(B*L)),B**4*c**2*m*(cosh(B*L))+ E*I*(B**3*sinh(B*L)),B**4*c**2*m*(sinh(B*L))+ E*I*(B**3*cosh(B*L))],[-B**2*cos(B*L),- B**2*sin(B*L),B**2*cosh(B*L),B**2*sinh(B*L)]]) 

    detA=AA.det() 

    from math import cosh, sinh, cos, sin 
    from scipy.optimize import fsolve 
    def func(B): 
     return detA 

    for i in range(1,15): 
     B0=fsolve(func,i) 
     print("B0=",B0) 
     print("w0=",B0**2*c) 

Но если я распечатаю detA, а затем скопирую/вставьте его в функцию, он будет работать нормально.

from math import cosh, sinh, cos, sin 
    from scipy.optimize import fsolve 
    def func(B): 
     return 5606.66666666667*B**11*sin(0.97*B)*cosh(0.97*B) - 5606.66666666667*B**11*cos(0.97*B)*sinh(0.97*B) + 478.634125*B**10*sin(0.97*B)**2 + 478.634125*B**10*cos(0.97*B)**2 - 957.26825*B**10*cos(0.97*B)*cosh(0.97*B) - 478.634125*B**10*sinh(0.97*B)**2 + 478.634125*B**10*cosh(0.97*B)**2 

    for i in range(1,15): 
     B0=fsolve(func,i) 
     print("B0=",B0) 
     print("w0=",B0**2*c) 

Я предполагаю, что, поскольку я использовал SymPy, определитель не плывут, как я мог бы превратить его плавать?

ответ

1

Использование lambdify преобразовать символическое выражение в функцию, которая может быть оценена численно:

func = sy.lambdify(B, detA) 

import sympy as sy 
import math 
from scipy.optimize import fsolve 

B = sy.symbols('B') 
L = 0.97 
a = 0.04 
b = 0.029 
m = 4.8 
I = 3.14 * (a/2) * (b/2)**3/4 
E = 10**10 
ro = 900 
A = 3.14 * (a/2) * (b/2) 
c = sy.sqrt(E * I/(ro * A)) 
AA = sy.Matrix([[0, -B**3, 0, B**3], [-B**2, 0, B**2, 0], [B**4 * c**2 * m * (sy.cos(B * L)) + E * I * (B**3 * sy.sin(B * L)), B**4 * c**2 * m * (sy.sin(B * L)) + E * I * (- B**3 * sy.cos(B * L)), B**4 * c**2 * m * (sy.cosh(B * L)) + E * I * (B**3 * sy.sinh(B * L)), B**4 * c**2 * m * (sy.sinh(B * L)) + E * I * (B**3 * sy.cosh(B * L))], [-B**2 * sy.cos(B * L), - B**2 * sy.sin(B * L), B**2 * sy.cosh(B * L), B**2 * sy.sinh(B * L)]]) 
detA = AA.det() 

func = sy.lambdify(B, detA) 

def func2(B): 
     return 5606.66666666667*B**11*math.sin(0.97*B)*math.cosh(0.97*B) - 5606.66666666667*B**11*math.cos(0.97*B)*math.sinh(0.97*B) + 478.634125*B**10*math.sin(0.97*B)**2 + 478.634125*B**10*math.cos(0.97*B)**2 - 957.26825*B**10*math.cos(0.97*B)*math.cosh(0.97*B) - 478.634125*B**10*math.sinh(0.97*B)**2 + 478.634125*B**10*math.cosh(0.97*B)**2 

assert abs(func(1) - func2(1)) < 1e-8 

for i in range(1, 15): 
    B0 = fsolve(func, i) 
    print("B0=", B0) 
    print("w0=", B0**2 * c) 
Смежные вопросы