Main Content

ODE 이벤트 위치

이 항목에서는 솔버 함수(ode45, ode15s 등)를 사용하여 ODE를 푸는 동안 이벤트를 감지하는 방법을 설명합니다.

이벤트 위치란?

일부 연립 ODE를 푸는 데 있어 어려운 부분은 풀이 작업을 중지하기에 적절한 시간을 결정하는 것입니다. 적분 구간의 최종 시간이 숫자가 아니라 특정 이벤트에 의해 정의될 수 있습니다. 나무에서 떨어지는 사과를 예로 들 수 있습니다. 사과가 땅에 떨어지면 ODE 솔버가 중지되어야 하지만, 사전에 이 이벤트가 언제 발생할지를 사용자가 알지 못할 수 있습니다. 이와 유사하게, 일부 문제는 풀이 작업을 종료시키지 않는 이벤트를 포함합니다. 행성의 궤도를 공전하는 위성을 예로 들 수 있습니다. 이 경우, 적분을 조기에 중지하기보다는 위성이 행성 궤도를 완전히 한 주기 공전하는 각각의 시점을 감지하기를 원할 것입니다.

ODE를 푸는 동안 특정 이벤트가 발생하는 시점을 감지하려면 이벤트 함수를 사용하십시오. 이벤트 함수는 사용자가 지정한 표현식을 사용하여 해당 표현식이 0인 경우 이벤트가 발생한 것으로 감지합니다. 또한 이벤트 함수는 이벤트를 감지하면 ODE 솔버에 적분을 중단하도록 신호를 보낼 수도 있습니다.

이벤트 함수 작성하기

odeset 함수의 'Events' 옵션을 사용하여 이벤트 함수를 지정합니다. 이벤트 함수의 일반적인 형식은 다음과 같아야 합니다.

[value,isterminal,direction] = myEventsFcn(t,y)

ode15i의 경우, 이벤트 함수는 yp에 대한 세 번째 입력 인수도 받아야 합니다.

출력 인수 value, isterminal 및 direction은 i번째 요소가 i번째 이벤트에 대응하는 벡터입니다.

  • value(i)는 i번째 이벤트를 설명하는 수학 표현식입니다. value(i)가 0일 때 이벤트가 발생합니다.

  • i번째 이벤트가 발생할 때 적분을 종료하려는 경우 isterminal(i) = 1이고, 그렇지 않은 경우에는 0입니다.

  • 모든 0을 찾으려는 경우에는 direction(i) = 0(디폴트 값)입니다. +1 값은 이벤트 함수가 증가하는 곳의 0만 찾고, -1은 이벤트 함수가 감소하는 곳의 0만 찾습니다. 모든 이벤트에 대해 디폴트 값 0을 사용하려면 direction = []을 지정하십시오.

나무에서 떨어지는 사과의 경우를 다시 한 번 살펴보겠습니다. 떨어지는 물체를 나타내는 ODE는 다음과 같습니다.

y''=−1+y' 2 ,

여기서 초기 조건은 y(0)=1과 y'(0)=0입니다. 이벤트 함수를 사용하여 y(t)=0인 시점, 즉 사과가 땅에 떨어지는 시점을 결정할 수 있습니다. 이 문제의 경우, 사과가 땅에 떨어지는 시점을 감지하는 이벤트 함수는 다음과 같습니다.

function [position,isterminal,direction] = appleEventsFcn(t,y)
  position = y(1); % The value that we want to be zero
  isterminal = 1;  % Halt integration 
  direction = 0;   % The zero can be approached from either direction
end

이벤트 정보

이벤트 함수를 지정한 경우 다음과 같이 3개의 추가 출력 인수를 사용하여 ODE 솔버를 호출합니다.

[t,y,te,ye,ie] = odeXY(odefun,tspan,y0,options)

솔버에서 반환하는 3개의 추가 출력값은 감지된 이벤트에 대응됩니다.

  • te는 이벤트가 발생한 시간으로 구성된 열 벡터입니다.

  • ye에는 te의 각 이벤트 시간에 해당하는 해 값이 포함됩니다.

  • ie에는 이벤트 함수에서 반환되는 벡터에 대한 인덱스가 포함됩니다. 이러한 값은 솔버가 감지한 이벤트를 나타냅니다.

또는, 다음과 같이 단일 출력값을 사용하여 솔버를 호출할 수도 있습니다.

sol = odeXY(odefun,tspan,y0,options)

이 경우, 이벤트 정보는 구조체에 sol.te, sol.ye, sol.ie로 저장됩니다.

제한 사항

이벤트 함수와 함께 ODE 솔버에 사용된 근 구하기 메커니즘에는 다음과 같은 제한 사항이 있습니다.

  • 적분의 첫 번째 스텝에서 종료 이벤트가 발생하면 솔버가 이 이벤트를 종료 이벤트가 아닌 것으로 등록하고 적분을 계속합니다.

  • 첫 번째 스텝에서 둘 이상의 종료 이벤트가 발생할 경우 첫 번째 이벤트만 등록되고 솔버가 적분을 계속합니다.

  • 0은 스텝 간 부호 교차 지점으로 결정됩니다. 따라서 스텝 간에 짝수 개수의 교차 지점이 있는 함수의 0은 누락될 수 있습니다.

솔버 스텝이 이벤트를 지나치는 경우 RelTol과 AbsTol을 줄여 정확도를 높여 보십시오. 또는, MaxStep을 설정하여 스텝 크기에 상한을 적용하십시오. tspan을 조정해도 솔버에서 취하는 스텝은 변경되지 않습니다.

단순한 이벤트 위치: 튀어 오르는 공

스크립트 열기

이 예제에서는 ODE 솔버에 사용할 단순한 이벤트 함수를 작성하는 방법을 보여줍니다. 예제 파일 ballode는 튀어 오르는 공의 운동을 모델링합니다. 이벤트 함수는 공이 튀어 오를 때마다 적분을 중단하고, 새로운 초기 조건을 적용해 적분을 다시 시작합니다. 공이 튀어 오를 때마다 적분이 여러 차례 중지되었다가 재시작됩니다.

튀어 오르는 공에 대한 방정식은 다음과 같습니다.

$$\begin{array}{cl} y'_1 &= y_2\\ y'_2 &= -9.8 . \end{array}$$

공의 높이 $y_1(t)$가 점점 감소하다가 0이 될 때 공이 튀어 오릅니다. 이 동작을 코딩하는 이벤트 함수는 다음과 같습니다.

function [value,isterminal,direction] = bounceEvents(t,y)
value = y(1);     % Detect height = 0
isterminal = 1;   % Stop the integration
direction = -1;   % Negative direction only

예제를 실행하고 공이 튀어 오르는 것을 시뮬레이션하는 이벤트 함수를 사용하는 방법을 보려면 ballode를 입력하십시오.

ballode

복잡한 이벤트 위치: 제한된 3체(Three Body) 문제

스크립트 열기

이 예제에서는 이벤트 함수의 방향 성분을 사용하는 방법을 보여줍니다. 예제 파일 orbitode는 한 물체가 그보다 훨씬 큰 두 물체의 궤도를 선회하는 제한된 3체 문제를 시뮬레이션합니다. 이벤트 함수는 궤도를 선회하는 물체의 가장 가까운 지점과 가장 먼 지점을 결정합니다. 궤도상의 가장 가까운 지점과 가장 먼 지점에서 이벤트 함수의 값이 동일하지만, 영점교차 지점의 방향에 따라 두 지점을 구분할 수 있습니다.

제한된 3체 문제의 방정식은 다음과 같습니다.

$$\begin{array}{cl} y'_1 &= y_3\\ y'_2 &= y_4 \\ y'_3 &= 2y_4 + y_1 -
\frac{\mu^* (y_1 + \mu)}{r_1^3} - \frac{\mu(y_1 - \mu^*}{r_2^3}\\ y'_4 &=
-2y_3 + y_2 - \frac{\mu^* y_2}{r_1^3} - \frac{\mu y_2}{r_2^3},
\end{array}$$

여기서

$$\begin{array}{cl} \mu &= 1/82.45\\ \mu^* &= 1-\mu\\ r_1 &=
\sqrt{(y_1+\mu)^2+y_2^2}\\ r_2 &= \sqrt{(y_1-\mu^*)^2 +
y_2^2}.\end{array}$$

처음 2개의 해 성분은 질량이 매우 작은 물체의 좌표이므로, 한 성분을 기준으로 다른 성분을 플로팅하면 물체의 궤도가 나타납니다.

orbitode.m에 중첩되어 있는 이벤트 함수는 두 개의 이벤트를 검색합니다. 한 이벤트는 출발점에서 가장 멀리 떨어진 지점의 위치를 찾고, 다른 이벤트는 우주선이 출발점으로 귀환하는 지점의 위치를 찾습니다. 적분기에서 사용하는 스텝 크기가 이벤트의 위치에 의해 정해지지 않는데도 불구하고 이벤트의 위치가 정확히 찾아집니다. 이 예제에서는 영점교차 지점의 방향을 지정하는 기능이 중요한 역할을 합니다. 출발점으로 귀환하는 지점과 출발점에서 가장 멀리 떨어진 지점 모두 동일한 이벤트 값을 가지며, 교차 지점의 방향에 따라 이러한 지점이 구분됩니다. 이 동작을 코딩하는 이벤트 함수는 다음과 같습니다.

function [value,isterminal,direction] = orbitEvents(t,y)
% dDSQdt is the derivative of the equation for current distance. Local
% minimum/maximum occurs when this value is zero.
dDSQdt = 2 * ((y(1:2)-y0(1:2))' * y(3:4)); 
value = [dDSQdt; dDSQdt];  
isterminal = [1;  0];         % stop at local minimum
direction  = [1; -1];         % [local minimum, local maximum]
end

예제를 실행하려면 orbitode를 입력하십시오.

orbitode

This is an example of event location where the ability to
specify the direction of the zero crossing is critical.  Both
the point of return to the initial point and the point of
maximum distance have the same event function value, and the
direction of the crossing is used to distinguish them.

Calling ODE45 with event functions active...

Note that the step sizes used by the integrator are NOT
determined by the location of the events, and the events are
still located accurately.

참고 항목

|

관련 항목