MATLAB에서 직접 2차 미방을 풀어 진자 운동 구현하기

얼마전에 저는 뭐 아무도 쓸일은 없을것 같았지만 그래도 몇 안되는 저의 취밍이자 흥미있어 하는 것이라 블로깅했던 글이 하나 있는데요. 바로 

[The Robot/Prog.Lang.] - MATLAB에서 4차 Runge Kutta를 이용하여 1차 혹은 2차 미분방정식을 푸는 예제

였습니다. MATLAB의 멀쩡한 Simulink나 ODE 명령이 있음에도 불구하고 과감하게(ㅠㅠ) 글을 올렸죠..ㅎㅎㅎ. 뭐 아무튼 그리고 그 글에 대한 응용 예제로 또 하나 후속글을 올릴려고 합니다. 그 예제로는 유명한 진자를 올릴려고 하죠. 사실 저는 꽤 예전에 단순한 진자(pendulum)를 대상으로 연재도 진했었습니만~~^^

저렇게 생긴 아이가 진자라는 아이죠^^. 저 아이의 동역학을 유도했던것은 [바로가기]부터 시작합니다. 뭐 이 글은 그 과정에는 관심이 없으니.. 살짝 결과만 가져오겠습니다. 그 중에서도 모터는 없이 자유 운동을 한다고 보죠.^^

그러면 수식이 위와 같습니다.^^. 위 수식을 이전에 했던 Runge Kutta로 2차 미분방정식을 풀던 예제[바로가기]에 적용하여 코드를 짜야죠^^

h = 0.01;
t = 0;
x = 30*pi/180;
v = 0;

fm = 0.05;
m = 0.1;
l = 100*0.01;
J = 0.02;
g = 9.8;

func = @(tVal, xVal, vVal) -fm/(m*l^2+J)*vVal-m*g*l/(m*l^2+J)*xVal;

figure;

set(gca, 'Xlim', [-80 80], 'Ylim', [-30,130], 'box', 'off', 'XColor', [0.801,0.801,0.801], 'YColor', [0.801,0.801,0.801], 'Color', [0.801,0.801,0.801])
hold(gca, 'on')

goundLine = line([-100, 100],[100, 100],'Color','k','LineWidth',12);
pendulum = line([0, l*100*sin(x)], [100, 100-l*100*cos(x)], 'Color','b','LineWidth',12);

while(1)
    kx1 = v;
    kv1 = func( t, x, v );
    
    kx2 = v + h*kv1/2;
    kv2 = func( t + h/2, x + h*kx1/2, v + h*kv1/2 );
    
    kx3 = v + h*kv2/2;
    kv3 = func( t + h/2, x + h*kx2/2, v + h*kv2/2 );
    
    kx4 = v + h*kv3;
    kv4 = func( t + h, x + h*kx3, v + h*kv3 );
    
    dx = h*(kx1 + 2*kx2 + 2*kx3 + kx4)/6;
    dv = h*(kv1 + 2*kv2 + 2*kv3 + kv4)/6;
    
    x = x + dx;
    v = v + dv;
    
    t = t + h;
    
    updateX = [0,l*100*sin(x)];
    updateY = [100, 100-l*100*cos(x)];
    
    set(pendulum,'Xdata',updateX,'Ydata',updateY);
    
    drawnow;
end

바로 위의 코드입니다. 너무 간단해서 딱히 설명할 것도 없네요. while문 안에 있는것이 Runge Kutta를 풀고 있는 것이구요. 그로부터 구해진 각도를 업데이트(updatedX/Y)해서 그래프를 그리는데 반영하는 것이지요.

미적감각이 무쟈게 아쉬운 저로서는 그냥 저렇게 선으로만 표현할 뿐이지요^^. 흠흠흠...^^

별것이 없습니다만... 그래도 이전에 다룬 Runge Kutta를 예제를 통해 직접 본다는 것과 MATLAB에서 간편하게 그 결과를 애니메이션을 확인할 때 어떻게 할 것인지에 대해 다루었다고 말하면 좀 억지죠???^^ 뭐 여하튼 전 여기까지 이야기할려구요.. ㅎㅎ^^