引數空間中的solve_ivp問題

我在創建 Matlab 中給出的示例 ODE 以使用 scipy 的 solve_ivp 時遇到問題。在 Matlab 中，函式定義為

function fixed_point_linear_center()
    clc; clf;

    stepsize=.5;
    xmin=-5;
    xmax=5;
    ymin=-5;
    ymax=5;

    [x,y] = meshgrid(xmin:stepsize:xmax,ymin:stepsize:ymax);

    A = [0 1;-1 0];

    dx = A(1,1)*x   A(1,2)*y;
    dy = A(2,1)*x   A(2,2)*y;

    % Strange scaling for nicer output, only "cosmetics"
    eunorm = ( dx.^2   dy.^2 ).^(0.35);
    dx = dx./eunorm;
    dy = dy./eunorm;

    quiver(x,y,dx,dy);
    axis([xmin xmax ymin ymax]);
    grid on; xlabel('x'); ylabel('y');

    tspan=[0 100];

    x0stepsize=0.25;
    for x0=xmin:x0stepsize:xmax
        hold on
        ic = [x0 0];
        [~,x] = ode45(@(t,x) f(t,x,A),tspan,ic);
        plot(x(:,1),x(:,2),'r');
        hold on
        ic = [0 x0];
        [~,x] = ode45(@(t,x) f(t,x,A),tspan,ic);
        plot(x(:,1),x(:,2),'r');
    end
    hold off
end

function dx = f(~,x,A)
    dx = A*[x(1); x(2)];
end

計算看起來像這樣的解決方案

，但是如果我像這樣在 python 中重新創建函式

def fixed_point_linear_center():
    stepsize   =  0.5
    x0stepsize =  0.25
    xmin       = -5
    xmax       =  5
    ymin       = -5
    ymax       =  5
    
    x    = np.arange(xmin, xmax stepsize, stepsize)
    xval = np.arange(xmin, xmax x0stepsize, x0stepsize)
    y    = np.arange(ymin, ymax stepsize, stepsize)
    yval = np.arange(ymin, ymax stepsize*0.25, stepsize*0.25) # evaluate 4 times for smoothness
    
    [X, Y] = np.meshgrid(x, y)
    
    
    A = np.array([[0,1],[-1,0]])

    dx = A[0,0]*X   A[0,1]*Y # 21x21
    dy = A[1,0]*X   A[1,1]*Y # 21x21
    
    f = lambda t,x,A : np.dot(A,[[x[0]],[x[1]]])

    # Strange scaling for nicer output, but only "cosmetics"
    eunorm = np.float_power(( dx**2   dy**2 ), 0.35)   #( dx**2   dy**2 )**0.35
    eunorm[10,10] = 0.001 # center is 0 which violates division
    dx = dx/eunorm
    dy = dy/eunorm

    plt.figure(figsize = (15,12))
    plt.quiver(X, Y, dx, dy, angles = 'xy', color='#0086b3', width=0.0015)
    plt.grid() 
    plt.xlabel('x') 
    plt.ylabel('y')

    plt.axis([xmin,xmax,ymin,ymax])
    
    tspan=[0,100]
    
    for x0 in xval:
        
        ic = [x0,0]
        #[~,x] = ode45(@(t,x) f(t,x,A),tspan,ic);
        solution = solve_ivp(f, [xmin, xmax], ic, method='RK45', t_eval=yval, dense_output=True, args=(A,))
        #solution = solve_ivp(f, [xmin, xmax], [x0], method='RK45', t_eval=yval, dense_output=False, args=(0,A))
        #solution = solve_ivp(f, [tmin, tmax], [ic], method='RK45', t_eval=tval, args=(A), dense_output=False)
        plt.plot(solution.y[1], solution.y[0],'r')
        
fixed_point_linear_center()

我收到類似的錯誤

ValueError：形狀（2,2）和（2,1,2）未對齊：2（dim 1）！= 1（dim 1）

或類似的，取決于我已經嘗試重寫f的內容。據我了解，solve_ivp 期望 x0 陣列中有一個值，而我回傳一個 2x1 向量。它也不接受向量作為它的 x0 陣列中的值，例如[[x0,0]]

現在我想知道 scipy.solve_ivp 是否能夠像 ode45 那樣對引數空間進行計算（以及我該怎么做），或者我是否必須以其他方式進行計算？

（我已經檢查過，所有其他矩陣和回傳值都與 matlab 計算相同。）

[編輯 2] 好的，現在可以使用了。x 的繪圖引數solution.y[1]當然必須是！

uj5u.com熱心網友回復：

就像 Matlab 求解器一樣，solve_ivp期望狀態是單個向量。改變

f = lambda t,x1,x2, A : np.dot(A,[[x1],[x2]])

到

f = lambda t, x, A : np.dot(A, x)

此外，為確保正確solve_ivp解釋引數的形狀，請將其A傳遞給solve_ivpwith args=(A,)（注意逗號的使用）。

標籤：Python MATLAB scipy 颂

上一篇：神經網路的成本函式計算

下一篇：Matlab：如何通過滑鼠單擊（使用“ButtonDownFcn”）將多個變數點從繪圖中保存到作業區中？