معادلات دیفرانسیل معمولی

[don.Diego]

معادلات دیفرانسیل معمولی به آندسته از معادلات دیفرانسیل اطلاق می شود که تمامی متغیرهای وابسته ( بعنوان مثال x، y و z ) تنها تابعی از یک متغیر مستقل (برای مثال t) می باشند. مثال آشنا در این مورد می­تواند معادله نوسانگر هارمونیک باشد که در آن مکان نوسانگر تنها تابعی از زمان است. در مقابل معادله دیفرانسیل معمولی معادلات دیفرانسیل جزیی قرار می گیرند که نمونه بسیار آشنای آن معادله شرودینگر است که در آن تابع موج (متغیر وابسته) تابعی از مختصات مکانی و زمان (یعنی متغیرهای مستقل ) می باشد.
ابتدا به حل مسایل مقدار اولیه (IVP) با معرفی سه روش اویلر (Euler’s method)، رانگ کوتا (Runge–Kutta method) و تابع آماده در متلب بانام ode45 می پردازیم.
1- روش اویلر
این ساده ترین روش است. با استفاده از تقریب مشتق عددی پیشرو الگوریتم تکرار بسادگی بدست می آید :

که در آن h طول گام است. کافست با توجه به معادله دیفرانسیل تابع f(x,t) مشخص گردد و در رابطه تکرار بالا جایگذاری شود. مثال ساده زیر را در نظر بگیرید.

dy/dt+5y=6 , y(t=0)=2

با توجه به معادله بالا خواهیم داشت : f(x,t)=dy/dt=6-5y
بنابراین رابطه تکرار بصورت زیر خواهد بود،

نهایتا برنامه کوتاه زیر را می­توان نوشت :

clc,clear all
 h=.01;
 t0=0;
 tf=2;
 t=t0:h:tf;
 y(1)=2;
 N=length(t);
 for i=1:N-1
     y(i+1)=y(i)+h*(6-5*y(i));
 end
 plot(t,y)

باید توجه داشته باشیم که هر سه روش بکار گرفته در حل مسایل
مقدار اولیه مبتنی بر حل معادله دیفرانسیل مرتبه یک می باشند در حالیکه بسیاری از معادلات دیفرانسیل مهم در فیزیک مرتبه دو هستند. در چنین مواردی با استفاده از تغییر متغیر هر معادله دیفرانسیل مرتبه دو به دو معادله دیفرانسیل مرتبه یک شکافته می شود و هر کدام از این دو معادله بصورت جداگانه با توجه به شرایط اولیه با استفاده از روش اویلر حل می­شود. بعنان مثال :

نوسانگر هارمونک ساده
با این توضیحات می­توانیم به حل یک مسئله استاندارد یعنی معادله حرکت نوسانگر هارمونیک بپردازیم. جسمی با جرم معین متصل به فنر تحت تاثیر نیروی بازگردانده در یک بعد حرکت نوسانی انجام می دهد با صرفنظر از نیروهای اتلافی (این مورد در روش رانگ کوتا بررسی می شود) معادله حرکت بشکل زیر است :

با تغییر متغیر v=dx/dt جفت روابط زیر خواهیم داشت،

و بنابراین روابط تکرار مطابق مثال قبلی بدست می آیند،

اگر فرض کنیم که نوسانگر در لحظه t=0 بیشترین فاصله از نقطه تعادل را دارد می توانیم شرایط اولیه زیر را در نظر بگیریم،

برنامه زیر (Oscilator1) معادله نوسانگر هارمونیک ساده را به روش اویلر حل می کند. نمودارها x و v برحسب زمان رسم می­شوند و در پایان حرکت واقعی نوسانگر نیز بصورت انیمیشن نمایش داده می­شود.

function Oscilator1
 % simple harmonic Oscilator- Euler's method
 % by saeed babanezhad - solid state physics
 close all
 k=1;
 m=1;
 w=sqrt(k/m);
 h=.1;        % time step
 t0=0;         % initial time
 tf=20;        % final time
 t=t0:h:tf;    % time interval
 N=length(t);
 x=zeros(1,N); % specify size of x matrix
 v=zeros(1,N); % specify size of v matrix
 x(1)=1;       % initial condition for position
 v(1)=0;       % initial condition for speed
 for i=1:N-1
     x(i+1)=x(i)+h*v(i);
     v(i+1)=v(i)-h*w^2*x(i);
 end
 figure(1)
 subplot(3,1,1)
 plot(t,x)
 ylabel('x')
 subplot(3,1,2)
 plot(t,v)
 ylabel('y')
 % animation
 subplot(3,1,3)
 % difine spring
 spring=line('lineStyle','none','Marker','d',...
     'MarkerSize',15,'Color','b','Xdata',[],'ydata',[]);
 % difine body
 body=line('Marker','s','MarkerSize',30,...
     'MarkerFaceColor','r','xdata',[],'ydata',[]);
 xlim([min(x) max(x)])
 for i=1:N
     set(body,'XData',x(i),'YData',0)
     set(spring,'XData',linspace(min(x),x(i),20)...
         ,'YData',zeros(1,20))
     drawnow
 end

در مقابل سادگی روش اویلر خطای آن نسبت به سایر روشها بیشتر است ( از مرتبه 2). در بخش بعدی می خواهم به روش رانگ کوتا ی مرتبه چهار بپردازم (خطا از مرتبه 4) که بسیار دقیقتر از روش اویلر است.

[don.Diego]

براي ضرب و تقسيم چند جمله ايها مي توانيد توابع deconv و conv را بكار ببريد. چند جمله ايهاي
و را در نظر بگيريد. حاصلضرب اين دو چند جمله اي به طريق زير بدست مي آيد:
» a=[1 1 1]; b=[1 -1];
» c=conv(a,b)
c =
1 0 0 -1
و تقسيم a/b نيز به صورت زير قابل محاسبه است:
» [q,r]=deconv(a,b)
q =
1 2
r =
0 0 3

[don.Diego]

يك چند جمله اي در MATLAB به صورت يك بردار سطري كه مولفه هاي آن ضرايب چندجمله اي به ترتيب نزولي هستند معرفي مي شود. براي مثال چند جمله ای p(x)= در MATLAB به شكل زير معرفي مي گردد:

» p=[1 0 -2 5];

ريشه هاي يك چند جمله اي را مي توانيد به صورت زير بدست آورد:

» r=roots(p)

r =

-2.0946

1.0473 + 1.1359i

1.0473 - 1.1359i

با دانستن ريشه هاي معادله مي توانيد ضرايب چند جمله اي مربوطه را محاسبه نمائيد:

» p2=poly(r)

p2 =

1.0000 0.0000 -2.0000 5.0000

[don.Diego]

تعريف كردن آرايه ها و عمليات جبري روي آنها
چهار نوع آرايه مي توان تعريف كرد: MATLAB در
١. اعداد اسكالر كه تك عضوي هستند.
٢. بردارها كه شامل يك سطر يا ستون مي باشند (يك بعدي).
٣. ماتريسها كه از اعضاي چيده شده در يك آرايش مربعي تشكيل مي گردند (دو بعدي).
٤. آرايه هاي با ابعاد بيش از دو.
اعضاي يك آرايه مي توانند عدد و يا حرف باشند و تفاوتي بين اعداد صحيح و اعشاري وجود ندارد.
مقدار جايگزين شده را MATLAB ، در صورت جايگزيني يك عدد و يا حرف در يك متغير
خاتمه يابد. semicolon بلافاصله نشان مي دهد مگر آنكه عبارت تعريف متغير با

» a=2.5
a =
2.5000
» a=3.2;
» a
a =
3.2000
» p='hello'
p =
hello
٣

بين حروف كوچك و بزرگ فرق قائل است: MATLAB

» A
??? Undefined function or variable 'A'.

از آنجا كه نشان دادن مقادير به شكل فوق قدري طولاني است معمولا" بهتر است كه در انتهاي
استفاده كرد. در صورتي كه اين عمل را فراموش كنيد و برنامه semicolon دستور معرفي متغير از
را فشار دهيد تا CONTROL C شروع به نشان دادن مقاذير يك آرايه طولاني نمايد كافي است كه
نشان دادن مقادير متوقف گردد. همانطور كه در بالا ديديد هميشه مي توان با نوشتن نام متغير
يك خط فاصله بين دستورها MATLAB مقدار آن را مشاهده نمود. همچنين مشاهده مي كنيد
مي گذارد. براي حذف اين خطوط اضافي مي توانيد از دستور زير استفاده كنيد:

» format compact
اكنون چند بردار تعريف مي كنيم:
» v=[1 2 3]
v =
1 2 3
» w=['abcd' '1234']
w =
abcd1234

براي تعريف بردارهاي عددي حتما" بايد از كروشه استفاده كرد ولي استفاده از آنها براي متغيرهاي
به عنوان جاي خالي استفاده MATLAB حرفي الزامي نيست. حالت خاصي از بردار (كه در توابع
بسياري دارد) عبارتست از بردار تهي كه به صورت [ ] تعريف مي گردد.
نحوه تعريف ماتريسها به صورت زير است:

» m=[1 2 3
4 5 6]
m =
1 2 3
4 5 6
» n=['abcd'
'1234']
n =
abcd
1234

اعضاي يك ماتريس را مي شود بطور جداگانه مشاهده كرد و يا تغيير داد

» m(2,3)
ans =
6
٤
» m(2,3)=7
m =
1 2 3
4 5 7

عمليات ساده جبري روي بردارها و ماتريسها به صورت زير انجام مي شود:
» 2*m
ans =
2 4 6
8 10 14
» m+1
ans =
2 3 4
5 6 8
» n1=[2 5 4
-1 -2 0];
» m+n1
ans =
3 7 7
3 3 7

لازم به ذكر است كه اعضاي يك سطر ماتريس را مي توان هم با فاصله و هم با ويرگول از هم جدا
در تعريف يك ماتريس به معناي انتقال به سطر بعدي مي باشد: semicolon

كرد. بكار بردن
» q=[1, 2, 3
4 5 6; 7 8 9]
q =
1 2 3
4 5 6
7 8 9

عملگر دو نقطه ( : ) كاربرد زيادي در رجوع به سطرها، ستونها و يا بخشي از آرايه دارد:

» q(1,:)
ans =
1 2 3
» q(:,2)
ans =
2
5
8
» q(1:2,2:end)
ans =
2 3
5 6
٥

استفاده نماييد: who اگر بخواهيد نام متغيرهاي ايجاد شده را ببينيد مي توانيد از دستور

» who
Your variables are:
a n q w
m p v

را بكار whos براي مشاهده نام متغيرهاي موجود به همراه اطلاعات اضافه تر در مورد آنها دستور
ببريد:

» whos
Name Size Bytes Class
a 1x1 8 double array
m 2x3 48 double array
n 2x4 16 char array
p 1x5 10 char array
q 3x3 72 double array
v 1x3 24 double array
w 1x8 16 char array
Grand total is 31 elements using 122 bytes

براي توليد بردارهاي عددي كه اعضاي آن به فاصله مساوي از هم قرار دارند روش ساده اي در
برداري باشد كه عضو اول آن ٠، عضو آخر آن ٢ و اعضاي t وجود دارد. فرض كنيد كه MATLAB
٠ از يكديگر باشند / آن به فاصله مساوي ٥

» t=0:.5:2
t =
0 0.5000 1.0000 1.5000 2.0000

آرايه هاي چند بعدي (آرايه هايي كه بيش از دو بعد دارند) از امكانات جديد پيش بيني شده در
كه قبلا" m هستند. به عنوان مثال مي توان بعد سوم را به شكل زير به ماتريس MATLAB 5
تعريف شده افزود:

» m(:,:,2)=ones(2,3)
m(:,:,1) =
1 2 3
4 5 7
m(:,:,2) =
1 1 1
1 1 1

افزودن بعدهاي چهارم و بيشتر نيز به طريق مشابه امكان پذير است. اصطلاحا" به بعد سوم صفحه
گفته مي شود ولي نام خاصي براي ابعاد چهارم به بعد وجود ندارد.
٦
استفاده كنيد: length براي بدست آوردن طول يك بردار مي توانيد از دستور

» length(t)
ans =
5

تعداد سطرها و ستونهاي يك ماتريس را نمايش مي دهد: size دستور

» size(n)
ans =
2 4

در مورد آرايه هاي چند بعدي برداري را مي دهد كه مولفه هاي آن طول آرايه در size استفاده از
هر يك از ابعاد آن است.

برخي از توابعي كه در ساختن آرايه ها بكار مي روند عبارتند از:

ones(2)

یک ماتریس 2*2 با مولفه هاي ١ ايجاد مي كند

ones(2,3)

یک ماتریس 3*2 با مولفه هاي ١ ايجاد مي كند

zeros(2)

یک ماتریس 2*2 با مولفه هاي 0 ايجاد مي كند

eye(3)

یک بردار یکه ایجاد می کند

linspace(- برداري با ٧ مولفه با فواصل مساوي بين ١- و ٥ ايجاد مي كند ( 1,5,7
linspace(- ١٠ و ١٠٢ ايجاد مي كند ( 1,2,8 - برداري با ٨ مولفه با فواصل لگاريتمي مساوي بين ١
تعدادي از توابعي كه روي آرايه ها عمل مي كنند عبارتند از:
sum(x) x حاصل جمع مولفه هاي
cumsum(x) از اول تا هر مولفه x حاصل جمع مولفه هاي
prod(x) x حاصلضرب مولفه هاي
cumprod(x) از اول تا هر مولفه x حاصلضرب مولفه هاي
max(x) را پيدا مي كند x بزرگترين مولفه
max(x) را پيدا مي كند x كوچكترين مولفه
sort(x) را مرتب مي كند x مولفه هاي
mean(x) x ميانگين حسابي مولفه هاي
std(x) x انحراف معيار مولفه هاي

[don.Diego]

مشتق چند جمله اي را مي توانيد با بكار بردن تابع polyder محاسبه كنيد.
» c=polyder(a)
c =
2 1
مشتق حاصلضرب دو چند جمله اي a*b را مي توانيد به صورت زير بدست آوريد:
» d=polyder(a,b)
d =
3 0 0
در صورتي كه تعداد آرگومانهاي خروجي تابع polyder برابر ٢ باشد، تابع مشتق تقسيم دو چند جمله اي جمله اي را تعيين مي نمايد:
» [q,d]=polyder(a,b)
q =
1 -2 -2
d =
1 -2 1



منبع : mathworks.ir