Filter löschen
Filter löschen

What's wrong with my ode45 function?

1 Ansicht (letzte 30 Tage)
Yuxing Zhang
Yuxing Zhang am 27 Nov. 2018
Kommentiert: Yuxing Zhang am 27 Nov. 2018
I am using the ode45 to solve streamfunction and I compelete my functions and main script but always have wrong in ode 45 function, can anyone tell me where is it wrong, A,B,C are just matrix provided ( they are all correct). Thank you very much!
clear all,clc
% initial conditions
n=64;
tspan= 0:0.5:4;
v= 0.001;
L=20;
x2=linspace(-L/2,L/2,n+1);
xspan=x2(1:n);
dx = xspan(2)-xspan(1);
x=x2(1:n);
y=x;
w0=exp(-x.^2-((y.^2)/20));
wvec=reshape(w0,n,1);
% matrix A,B,C ( you can ignore this)
N=n*n; % total size of matrix
e0=zeros(N,1); % vector of zeros
e1=ones(N,1); % vector of ones
e2=e1; % copy the one vector
e4=e0; % copy the zero vector
for j=1:n
e2(n*j)=0; % overwrite every n^th value with zero
e4(n*j)=1; % overwirte every n^th value with one
end
e3(2:N,1)=e2(1:N-1,1); e3(1,1)=e2(N,1); % shift to correct
e5(2:N,1)=e4(1:N-1,1); e5(1,1)=e4(N,1); % positions
a=spdiags([e1 e1 e5 e2 2*e1 e3 e4 e1 e1], ...
[-(N-n) -n -n+1 -1 0 1 n-1 n (N-n)],N,N);
A=full(a)./((dx)^2);
b=spdiags([e1 -e1 e0 e1 -e1],[-n*(n-1) -n 0 n n*(n-1)],N,N);
B=full(b)./(2*dx);
c=spdiags([e5 -e2 e0 e3 -e4 ],[ -n+1 -1 0 1 n-1],N,N);
C=full(c)./(2*dx);
[t,wtsol]=ode45(@(t,wvec)rhseq(t,wvec,psivec,A,B,C,v),tspan,w0);
function wtsol=rhseq(t,wvec,psivec,A,B,C,v)
psivec=A\wvec;
wtsol=(C*psivec).*(B*wvec)-(B*psivec).*(C*wvec)+(v*A*wvec);
end

Antworten (1)

Stephan
Stephan am 27 Nov. 2018
  1. you want to pass psivec to the function, but it is calculated inside the function.
  2. A is size 4096x4096 and wvec is size 64x1. Due to this A\wvec wont work.

Tags

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by