Кодовое уравнение

Создайте функцию для кодирования уравнений. Эта функция должна иметь подпись dydx = mat4ode(x,y,lambda), где:

Можно записать уравнение Матье как систему первого порядка, используя замены ${\mathit{}}_{\mathrm{y1}}\mathit{=}$ y ${\mathit{и}}_{\mathrm{}}{\mathit{y2}}^{}$= y ′,

${{\mathit{}}_{\mathrm{}}}^{y1\; }={\mathit{}}_{\mathrm{y2}}$,

${{\mathit{}}_{\mathrm{}}}^{}y2\prime =-(\mathrm{}\mathit{}\mathrm{\lambda -2q}\cos \mathit{}({\mathit{2x}}_{\mathrm{)}}$) y1.

Затем выполняется соответствующая функция.

function dydx = mat4ode(x,y,lambda) % equation being solved
dydx = [y(2)
      -(lambda - 2*q*cos(2*x))*y(1)];
end

Примечание.В конце примера все функции включаются как локальные.

Граничные условия кода

Теперь запишите функцию, которая возвращает остаточное значение граничных условий в граничных точках. Эта функция должна иметь подпись res = mat4bc(ya,yb,lambda), где:

Эта задача имеет три граничных условия в интервале $[0,\delta$]. Для вычисления остаточных значений необходимо поместить граничные условия в форму $\mathit{}g\mathit{}(\mathit{x},\mathrm{y}$) = 0. В этой форме граничными условиями являются

${\mathit{y}}^{}\prime \mathrm{}(0\mathrm{)}$= 0,

${\mathit{y}}^{}\prime (\delta \mathrm{)}$= 0,

$\mathit{y}\mathrm{(}0)\mathrm{}-1$= 0.

Затем выполняется соответствующая функция.

function res = mat4bc(ya,yb,lambda) % boundary conditions
res = [ya(2)
       yb(2)
       ya(1)-1];
end

Создать начальное предположение

Наконец, создайте первоначальное предположение о решении. Необходимо предоставить начальное предположение для обоих компонентов решения ${\mathit{}}_{\mathrm{y1}}\mathit{=}\mathit{y}$(x) ${\mathit{}}_{\mathrm{и}}{\mathit{}}^{y2}\mathit{=}$y ′ (x), а также неизвестный $$параметр λ. Можно вычислить только собственные значения и собственные функции, близкие к начальным предположениям. Чтобы увеличить вероятность того, что вычисленная собственная функция соответствует четвертому собственному значению, следует выбрать начальное предположение, которое имеет правильное качественное поведение.

Для этой задачи функция косинуса делает хорошую начальную догадку, поскольку она удовлетворяет трем граничным условиям. Закодируйте начальное предположение для $\mathit{y}$ с помощью функции, которая возвращает предположение для ${\mathit{}}_{\mathrm{y1}}$ и ${\mathit{}}_{\mathrm{y2}}$.

function yinit = mat4init(x) % initial guess function
yinit = [cos(4*x)
        -4*sin(4*x)];
end

Звонить bvpinit используя сетку из 10 точек в интервале $[0,\delta$], функцию начального приближения и предположение 15 для значения $$λ.

lambda = 15;
solinit = bvpinit(linspace(0,pi,10),@mat4init,lambda);

Уравнение решения

Звонить bvp4c с функцией ОДУ, функцией граничного условия и начальным предположением.

sol = bvp4c(@mat4ode, @mat4bc, solinit);

Значение параметра

Печать значения неизвестного параметра $\lambda$, найденного bvp4c. Это значение является четвертым собственным значением ($\mathit{q}\mathrm{=}$5) уравнения Матье.

fprintf('Fourth eigenvalue is approximately %7.3f.\n',...
            sol.parameters)

Fourth eigenvalue is approximately  17.097.

Решение для построения графика

Использовать deval для оценки решения, вычисленного bvp4c в 100 точках в интервале $[0,\delta$].

xint = linspace(0,pi);
Sxint = deval(sol,xint);

Постройте график обоих компонентов решения. На графике показана собственная функция (и её производная), связанная с четвёртым собственным значением ${}_{\mathrm{\lambda \; 4}}\mathrm{=}\mathrm{}$17,097.

plot(xint,Sxint)
axis([0 pi -4 4])
title('Eigenfunction of Mathieu''s Equation.') 
xlabel('x')
ylabel('y')
legend('y','y''')

Локальные функции

Здесь перечислены локальные вспомогательные функции, которые решатель BVP bvp4c вызывает для вычисления решения. Кроме того, эти функции можно сохранить в виде собственных файлов в каталоге по пути MATLAB.

function dydx = mat4ode(x,y,lambda) % equation being solved
q = 5;
dydx = [y(2)
      -(lambda - 2*q*cos(2*x))*y(1)];
end
%-------------------------------------------
function res = mat4bc(ya,yb,lambda) % boundary conditions
res = [ya(2)
       yb(2)
       ya(1)-1];
end
%-------------------------------------------
function yinit = mat4init(x) % initial guess function
yinit = [cos(4*x)
        -4*sin(4*x)];
end
%-------------------------------------------