OCTAVE求解常微分方程

如何解决OCTAVE求解常微分方程

我必须用倍频程求解以下微分方程

哪里

A是(N-1)X(N-1)的矩阵，
f定义为：f(x,t)=1000*sqrt(|1-t|)，
μ=1，
u(x,0)=1000x(1-x)(1+(3/2)x^3)，
所有L=1的酒吧都有u(0,t)=u(1,t)=0和t长度，
h是1/N（N个间隔的分区数）。

我已经阅读了如何使用该命令，但是在这种特定情况下我不理解如何使用它。如何输入功能du / dt。非常感谢您的帮助。

解决方法

我不是专家，您遗漏了一些信息（例如x的特定性质）。

但是，作为下面的示例，下面的示例作为如何“运行” lsode函数的示例：

您提到了一个长度为1的“条”和N个分区/箱，因此我假设X是一个矢量，表示箱的“边缘”处的对应位置，在这种情况下，您有N + 1条边。由于为您提供了“外部”边缘的解决方案（即x = 0和x = 1），因此我们只关心“内部边缘”，因此我们只需要处理N-1个元素的向量（因此大小矩阵A）。

  N         = 4;   % number of bins
  mu        = 1;
  h         = 1 / N;
  A_matrix  = rand( N - 1,N - 1 );
  t_rvector = [0 : 20].';

  x_cvector = linspace( 0,1,N + 1 ).';   % N+1 edges,including outer ones
  x_cvector = x_cvector(2 : end - 1);      % discard outer edges,keep inner only.


% Create f(t) as defined in the question,but vectorised w.r.t. x
% (i.e. outputting an N-1 element vector that can be added to a corresponding u vector)

  f_function = @(t) 1000 .* sqrt( abs( 1 - t ) ) .* ones( size( x_cvector ) );


% Create a u0 vector (of N-1 elements,corresponding to inner edges on the bar)
  u0_cvector = 1000 .* x_cvector .* (1 - x_cvector) .* (1 + (3 / 2) .* x_cvector .^ 3);

% Create a 'wrapper' function of the form expected by `lsode`
  ode_function = @(u_cvector,t) -mu ./ (h .^ 2) .* A_matrix * u_cvector + f_function(t);

% Run `lsode`
  [u_matrix,ISTATE,MSG] = lsode( ode_function,u0_cvector,t_rvector )

% Pad u_matrix with 0s on left and right,corresponding to x=0 and x=1 values
  TimePoints = length( t_rvector );
  ZeroColumn = zeros( TimePoints,1 );
  u_matrix   = [ ZeroColumn,u_matrix,ZeroColumn ];

% Plot the curves (when plotting a matrix,each column is treated as a separate plot),therefore if you want to plot along the 'components' of u (corresponding to positions in $
  plot ( u_matrix.' );

这可能与您的问题不同，但是上面的代码按预期运行，因此希望它可以作为您的起点。

值得注意的是，传递到u函数中的lsode向量必须采用“列”形式。（时间点列表可以是列或行，没关系）。