尋找未知常微分方程

Question

鑑於

d²x/dt² + a·dx/dt + 7.9·x³ = 3.2·sin(xt) 

初始條件

x(0)  = +1.2 
dx/dt(0) = −3.3 
x(2.3) = −0.6 

查找數值的a所有可能的值，每個值精確到至少3個顯著數字。

除了蠻力還有解決這個問題的方法嗎？

Answer 1

據我所知，這是不可能解決這個問題，如上所述。

這是我做的。我實現您的問題在一個合理的一般方式：

%{ 

Find all 'a' for which 

d²x/dt² + a·dx/dt + 7.9·x³ - 3.2·sin(xt) = 0 

with initial conditions 

x(0)  = +1.2 
dx/dt(0) = −3.3 
x(2.3) = −0.6 

%} 

function odetest 

    % See how the function search_a(a) behaves around a = 0: 
    test_as = 0 : 0.1 : 10; 
    da = zeros(size(test_as)); 
    for ii = 1:numel(test_as)   
     da(ii) = search_a(test_as(ii)); end 

    figure(100), clf, hold on 
    plot(test_as, da) 
    axis tight 
    xlabel('a') 
    ylabel('|x(2.3) - 0.6|') 


    % Roughly cherry-pick some positive values, improve the estimate, and 
    % plot the solutions 

    opt = optimset('tolfun',1e-14, 'tolx',1e-12); 

    plot_x(fminsearch(@search_a, 0.0, opt), 1) 
    plot_x(fminsearch(@search_a, 1.4, opt), 2) 
    plot_x(fminsearch(@search_a, 3.2, opt), 3) 

    % Plot single solution 
    function plot_x(a,N) 

     [xt, t] = solve_ode(a); 

     figure(N), clf, hold on 
     plot(t,xt) 
     plot(2.3, -0.6, 'rx', 'markersize', 20) 
     title (['x(t) for a = ' num2str(a)]) 
     xlabel('t') 
     ylabel('x(t)') 
    end 
end 

% Solve the problem for a value a, and return the difference between the 
% actual value and desired value (-0.6) 
function da = search_a(a) 

    a_desired = -0.6; 

    xt = solve_ode(a);  
    da = abs(xt(end) - a_desired); 
end 


% Solve the problem for any given value of a 
function [xt, t] = solve_ode(a) 

    y0  = [1.2 -3.3]; 
    tfinal = 2.3; 

    opt = odeset('AbsTol',1e-12, 'RelTol',1e-6);  
    [t,yt] = ode45(@(y,t) odefun(y,t,a), [0 tfinal], y0, opt);  
    xt  = yt(:,1); % transform back to x(t) 
end 

% Most ODE solvers solve first-order systems. This is not a problem for a 
% second-order system, because if we make the transformation 
% 
% y(t) = [ x (t) 
%   x'(t) ] 
% 
% Then we can solve for 
% 
% y'(t) = [ x' (t) 
%    x''(t) ] <- the second-order homogeneous DE 
% 
function dydt = odefun(t,y,a)  
    dydt = [y(2) 
      -a*y(2) - 7.9*y(1)^3 + 3.2*sin(y(1)*t)];  
end 

，第一部分給了我這個人物：

一些進一步的調查表明，這僅增長較大a。

這個數字引起了最初的估計a = [0, 1.4, 3.2]，我通過fminsearch()改善，繪製的解決方案：

所以，這可能使你在交功課： ）

但是，爲什麼我說這是不可能回答這個問題，因爲這是第一個情節看起來像負a：

的振盪行爲似乎無限期地持續下去，並在零點之間的間距似乎在不可預測的方式來減少。

現在，我的大學日子已經很長時間了，我對ODE理論不再那麼瞭解了。也許有是一個模式，它只是不顯示，因爲數值問題。或者也許振盪在一些值之後停止，永不再返回。或者在a = +1053462664212.25可能會出現另一個零。

我無法證明任何這些東西，我只知道如何蠻力;其餘的由你決定。