이 페이지의 최신 내용은 아직 번역되지 않았습니다. 최신 내용은 영문으로 볼 수 있습니다.

포식자-피식자 방정식 풀기

이 예제에서는 ode23ode45를 모두 사용하여 포식자/피식자 모델을 나타내는 미분 방정식을 푸는 방법을 보여줍니다. ode23과 ode45는 가변 스텝 크기 룽게-쿠타 적분 방법을 사용하여 상미분 방정식의 수치 해를 구하기 위한 함수입니다. ode23은 중간 정도의 정확도를 얻기 위해 단순하게 식의 2차와 3차 항을 사용하고, ode45는 더 높은 정확도를 얻기 위해 4차와 5차 항을 사용합니다.

로트카-볼테라 방정식 또는 포식자-피식자 모델로 알려진 1계 상미분 방정식 쌍을 사용해 보겠습니다.

dxdt=x-αxydydt=-y+βxy.

변수 xy는 각각 피식자 개체군과 포식자 개체군의 크기를 측정합니다. 2차 교차 항은 종 간의 상호 작용을 설명합니다. 포식자가 없으면 피식자 개체군이 증가하고 피식자가 부족하면 포식자 개체군이 감소합니다.

방정식 코딩하기

시스템을 시뮬레이션하기 위해, 상태 값과 시간 값이 제공될 경우 상태 도함수로 구성된 열 벡터를 반환하는 함수를 만듭니다. MATLAB에서 두 변수 xy는 벡터 y의 처음 2개 값으로 나타낼 수 있습니다. 마찬가지로 도함수 역시 벡터 yp의 처음 2개 값입니다. 함수는 ty에 대한 값을 받은 후 yp의 방정식에서 생성한 값을 반환해야 합니다.

yp(1) = (1 - alpha*y(2))*y(1)

yp(2) = (-1 + beta*y(1))*y(2)

이 예제에서 방정식은 lotka.m이라는 파일에 포함되어 있습니다. 이 파일은 파라미터 값 α=0.01β=0.02를 사용합니다.

type lotka
function yp = lotka(t,y)
%LOTKA  Lotka-Volterra predator-prey model.

%   Copyright 1984-2014 The MathWorks, Inc.

yp = diag([1 - .01*y(2), -1 + .02*y(1)])*y;

시스템 시뮬레이션하기

ode23을 사용하여 구간 0<t<15에 대해 lotka에 정의된 미분 방정식을 풉니다. 포식자 개체군과 피식자 개체군이 동일하도록 초기 조건 x(0)=y(0)=20을 사용합니다.

t0 = 0;
tfinal = 15;
y0 = [20; 20];   
[t,y] = ode23(@lotka,[t0 tfinal],y0);

결과 플로팅하기

시간에 대해 결과로 나타나는 개체군을 플로팅합니다.

plot(t,y)
title('Predator/Prey Populations Over Time')
xlabel('t')
ylabel('Population')
legend('Prey','Predators','Location','North')

이제 개체군 간을 비교해 플로팅합니다. 결과로 나타나는 위상 평면 플롯에서는 개체군 간 순환 관계가 훨씬 명확해집니다.

plot(y(:,1),y(:,2))
title('Phase Plane Plot')
xlabel('Prey Population')
ylabel('Predator Population')

여러 솔버의 결과 비교하기

ode23 대신 ode45를 사용하여 다시 시스템을 풉니다. ode45 솔버는 각 스텝마다 시간이 더 오래 걸릴 뿐만 아니라 더 많은 스텝을 진행합니다. 그렇지만 ode45의 출력값은 매끄럽게 표시됩니다. 이는 기본적으로 솔버가 연속적 확장 식을 사용하여 각 스텝 범위 안에서 시간 지점을 4개의 균일한 간격으로 나눠 출력값을 생성하기 때문입니다. ('Refine' 옵션을 사용하면 지점 개수를 조정할 수 있습니다.) 비교를 위해 2개의 해를 모두 플로팅합니다.

[T,Y] = ode45(@lotka,[t0 tfinal],y0);

plot(y(:,1),y(:,2),'-',Y(:,1),Y(:,2),'-');
title('Phase Plane Plot')
legend('ode23','ode45')

이 결과는 서로 다른 수치 계산법을 사용하여 미분 방정식을 풀어도 생성되는 답은 약간의 차이만 있음을 보여줍니다.

참고 항목

|

관련 항목