寻找未知常微分方程

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可能会出现另一个零。

我无法证明任何这些东西，我只知道如何蛮力;其余的由你决定。