حل دستگاه معادلات دیفرانسیل به روش ODE در متلب
در پست «حل عددی معادله دیفرانسیل در متلب به کمک ode» فراگرفتیم که چگونه می توان یک معادله دیفرانسیل را به کمک ode در متلب به صورت عددی حل نمود. در این پست می خواهیم فرا بگیریم که چگونه می توان یک دستگاه معادله دیفرانسیل را به کمک ode حل نمود.
پس قبل از اینکه این پس رو مطالعه بفرمائید حتما پس مربوط به «حل عددی معادله دیفرانسیل در متلب به کمک ode» را بخوانید.
فرض کنید که دستگاه معادله دیفرانسیل مقابل با شرایط مرزی زیر داده شده است:
y\" z-z\'=cos(x)
z\" z\' y\'=-sin(x)
y(0)=0
y\'(0)=1
z(0)=1
z\'(0)=0
همانطور که می دانید ode در متلب تنها می تواند معادلات دیفرانسیل مرتبه اول را حل کند. پس اولین کار قبل از رفتن به متلب این است که معادلات یا دستگاه معادلات را به معادلات مرتبه اول تبدیل نمائیم. پس در این مقاله معادله اول تبدیل به دو معادله و معادله دوم نیز تبدیل به دو معادله خواهد شد. این تبدیلات در زیر آمده است:
y=y1
y\'1=y2
y\'2=cos(x) z\'-z
-----------------------
z=z1
z\'1=z2
z\'2=-sin(x)-y\'-z\'
-----------------------
با تلفیق این دو دسته معادلات، دستگاه معادلات جدید به صورت زیر بدست خواهد آمد:
y\'1=y2
y\'2=cos(x) z2-z1
z\'1=z2
z\'2=-sin(x)-y2-z2
شرایط مرزی نیز به صورت زیر تبدیل خواهد شد:
y1(0)=0
y2(0)=1
z1(0)=1
z2(0)=0
حال باید دستگاه معادله 4 گانه تبدیل شده را در یک function در متلب بنویسیم:
function out=fun(x,s)
y1=s(1)
y2=s(2)
z1=s(3)
z2=s(4)
out=[y2;cos(x) z2-z1;z2;-sin(x)-y2-z2];
توجه کنید که متغیرهای ما تبدیل به 4 متغیر y1 و y2 و z1 و z2 شده اند. اما همانطور که می دانید در ode در هنگام نوشتن function شما تنها یک متغیر مستقل (در اینجا x) و یک متغیر وابسته (در اینجا s) دارید. به همین دلیل s دارای 4 مولفه خواهد بود که در function آمده است. یعنی به طور مثال اگرکه یک دستگاه معادله را شما تغییر دادید به معادلات مرتبه اول و 10 معادله مرتبه اول شد، آنگاه s دارای 10 مولفه خواهد بود.
پس از اینکه نوشتن function تمام شد حالا در m فایل دیگری و یا در command window دستور ode45 را اجرا می کنیم:
>>x_span=[0 2*pi];
>>IC=[0 1 1 0];
>>[X OUT]=ode45(@fun,x_span,IC);
با اجرای این دستور معادلات برای x بین صفر تا 2*pi حل می گردد. توجه نمائید که بردار IC بردار شرایط اولیه خواهد بود که به ترتیب معادلاتی که در function نوشته اید باید این شرایط اولیه را بنویسید. خروجی این دستورات بردار X و ماتریس OUT می باشد. همانطور که میبینید ماتریس OUT دارای 4 ستون می باشد. به ترتیبی که معادلات را در function نوشته اید ستون اول y1 و ستون دوم y2 و ستون سوم z1 و ستون چهارم z2 می باشد. با توجه به تبدیلاتی که انجام داده ایم ما به دنبال جواب y و z هستیم که در اینجا می شود y1 و z1. پس اگر که بخواهیم جواب ها را با همدیگر رسم نمائیم کافیست دستورات زیر را اجرا کنیم:
>>plot(x,OUT(:,1),\'r\',x,OUT(:,3),\'b\')
>>legend(\'y\',\'z\')
همانطور که از شکلها هم می توان فهمید جواب y همان نمودار قرمز رنگ یعنی sin و جواب z همان نمودار آبی رنگ یعنی cos می باشد.
