Matematické Fórum

Zdravím,

mám za úkol naprogramovat metodou AM6 následující rovnici





metodu mám naprogramovanou zde:

function [t,y]=method(ode,tspan,y0,Q)
%% metoda AM6
%
% vstupni data:
%
% ode   ... prava strana diferencialni rovnice, funkce promennych (t,y),
%           pole typu (d,1)
% tspan ... interval, v nemz pocitame reseni, a=tspan(1) je pocatek,
%           b=tspan(end) je konec intervalu integrace, vektor delky d
% y0    ... pocatecni podminka, vektor delky d, y0(i) je pocateni
%           podminka pro i-tou slozku
% Q     ... pocet dilku, t(n)=a+n*tau, n=0,1,...,Q, kde tau=(b-a)/Q
%
% vystupni data:
%
% t ... deleni, pole typu (1,Q+1)
% y ... reseni, pole typu (Q+1,d), y(n,i) je i-ta slozka reseni v uzlu t(n)
%         
% vzorove volani:
%
% pomoci programu go.m:
% go(1) % ==> prvni testovaci priklad
% go(2) % ==> druhy testovaci priklad
%%
k=6;
milne=true;
if nargin<6, milne=true; end;

d=length(y0);        % pocet rovnic
a=tspan(1);          % zacatek integrace
b=tspan(end);        % konec integrace
tau=(b-a)/Q;         % delka kroku
t=linspace(a,b,Q+1); % deleni
y=zeros(Q+1,d);      % pseudodeklarace
f=zeros(Q+1,d);      % pseudodeklarace derivace reseni
y(1,:)=y0(:)';       % pocatecni podminka do prvniho radku
p(1,:)=ode(t(1),y(1,:)')';    % vlozeni odpovidajicich hodnot prave strany
atol=1e-14;          % absolutni tolerance
rtol=1e-12;          % relativni tolerance
tol=[atol rtol];     % tolerance pro solver nsoli

%startovací hodnoty dle PDF

tspan=[t(1),t(1)+(k-1)*tau]; % startovaci uzly
opt=odeset('AbsTol',1e-14,'RelTol',1e-11);
sol=ode45(ode,tspan,y0(:),opt);
yy=deval(sol,tspan(1):tau:tspan(2))';
y(1:k,:)=yy; %y(2:k,:)=yy(2:k,:);
for j=1:k,
     f(j,:)=ode(t(j),y(j,:)')';
end;
  atol=1e-12;      % absolutni tolerance pro solver nsoli
  rtol=1e-9;       % relativni tolerance pro solver nsoli
  tol=[atol rtol]; % tolerance pro solver nsoli

%% vypocet y_n, n=k+1,...,Q+1
est=0;  % pocatecni nastaveni pro urceni maximalni lokalni chyby
for n=k:Q     % cyklus pres jednotlive uzly deleni
  tn=t(n);    % uzel t_n                   
  yn=y(n,:)'; % reseni v uzlu t_n
  fn=f(n,:)'; % derivace v uzlu t_n

F=@(z) -z+yn+1/1440*tau*(475*f(tn,z)+1427*f(n,:)'-798*f(n-1,:)'+482*f(n-2,:)'-173*f(n-3,:)'+27*f(n-4,:)'); % C: AM6
z=nsoli(y(n,:)',F,tol)';           % reseni F(z)=0
y(n+1,:)=z; 
                                            %y(n+1,:)=nsoli(yn,F,tol)';                         % vlozeni ynew           
f(n+1,:)=ode(tn+tau,y_new)';             % E:
end

kdy mi ale dělá problém F=@(z), která musí projít přes funkci nsoli, čehož nemůžu dosáhnout. Pořád mi to vypisuje Subscript indices must either be real positive integers or logicals. Neví si s tím někdo rady?

Děkuji.