Could anyone check my ode45() code please?
조회 수: 3 (최근 30일)
이전 댓글 표시
NB
Moved the following from title to body but have no idea what the real question is...dpb
Because nu represents the mean value so there is someone who expected to got my frequency at 1!!
ORIGINAL posting picks up here again with code...end dpb edits
My system is
function dgdt=stochasticnu(t,g,S0,Sr,nu)
% The system consists of two competing species,Grass density and Woody-vegetation density with cattle-stocking rate
dgdt(1)=1.5*g(1)*(1-(S0+Sr*cos(nu*t))-0.7*g(1)-g(2));
dgdt(2)=0.03 + g(2)*(1-2*g(1)-1.03*g(2));
dgdt=[dgdt(1) dgdt(2)]';
end
My code is:
function RunStocnnu
% Finding the DFT(Discrete Fourier Transform ) using FFT procedure to get
% the spectrum plot of the system of two competing species, Grass
% density and Woody-vegetation density after adding a small perturbation
% to the cattle-stocking rate.
tinit=0; % The initial time
tf=50; % The final time
h=0.1; % Step size
tspan=tinit:h:tf; % The integration time vector
S0=.3;
Sr=.1;
g0=[0.1,0.1]; % The initial conditions
nu=1;
[t,g]=ode45(@stochasticnu,tspan,g0,[],S0,Sr,nu); % ODE solver
% Begin a new figure
figure(01)
plot(g(:,1),g(:,2),'k*-') % plotting the phase plane of the two species
xlabel('Grass density') % The x axis
ylabel('Woody-vegetation density') % The y axis
title('A trajectory of stochastic system') % The title of the figure
n=length(g); % The length of the input data vector
dt=t(end)/(length(g)-1);
Fs=1/dt;
% T=1/Fs;
NFFT =n;
y=fft(g,NFFT)/n;
y=y(1:NFFT/2);
mx=abs(y);
f=(0:NFFT/2-1)*(Fs/NFFT); % The frequency range
figure(02)
plot(f(1:50),y(1:50)); % FFT is symmetric, throw away second half
xlabel('Frequency')
figure(03)
plot(f,mx)
댓글 수: 2
dpb
2014년 9월 1일
Please recast the question more clearly and provide an appropriate question title. I moved the long semi-question into body but can't interpret the issue well enough to try to take a stab at an answer, sorry...
Image Analyst
2014년 9월 1일
I have no idea what the question is. One good point though is I was able to copy and paste the code and run it, and it did successfully produce these plots and didn't give any error messages.
답변 (0개)
참고 항목
카테고리
Help Center 및 File Exchange에서 Ordinary Differential Equations에 대해 자세히 알아보기
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!