کدهای متلب محاسبات عددی | دانلود رایگان مجموعه کامل

در این مطلب روش‌ها و کدهای متلب محاسبات عددی را معرفی می‌کنیم. محاسبات عددی روش‌ها و ویژگی‌های متفاوتی دارند. در این مبحث ضمن معرفی این روش‌ها، کدها و حل دستی برخی از آن‌ها را معرفی می‌کنیم. اگر شما نیز برای حل مسائل ریاضی، مهندسی یا فیزیکی خود به حل عددی نیازمند هستید خواندن این مطلب را از دست ندهید.

محاسبات عددی در متلب

محاسبات عددی روشی برای حل مسائل پیچیده ریاضی است که فقط با استفاده از عملگرهای ساده ریاضی انجام می‌شود. این روش شامل فرمول‌بندی فیزیکی مدل‌های ریاضی و موقعیت‌هایی است که می‌توان با عملگرهای ریاضی حل کرد. برای محاسبات عددی نیاز به توسعه، تجزیه و تحلیل و استفاده از الگوریتم‌ها است.

فیلم آموزش محاسبات عددی با متلب MATLAB در فرادرس

کلیک کنید

این محاسبات شامل تعداد زیادی محاسبات ریاضی هستند و بنابراین به دستگاه‌های محاسباتی سریع و کارآمد احتیاج دارند. انقلاب میکروالکترونیک و متعاقب آن توسعه رایانه‌های شخصی ارزان و کم هزینه تأثیر به سزایی در استفاده از روش‌های عددی برای حل مشکلات علمی داشته است. روش‌های معادلات عددی برای حل حالت‌های زیر اتفاق می‌افتد:

• معادلات غیر خطی تک متغیره که حل تحلیلی و بسته ندارند و برای آن‌ها باید از حل عددی استفاده کرد.
• روش درون‎‌یابی و برازش برای پیدا کردن معادلاتی که داده‌هایی از آن را داریم.
• حل عددی مشتق‌ها و انتگرال‌هایی که حل تحلیلی و دقیقی ندارند.
• حل عددی معادلات دیفرانسیل معمولی یا ODE
• حل عددی دستگاه معادلات خطی یا غیرخطی. البته حل تحلیلی دستگاه معادلات خطی غالباً امکان‌پذیر است اما در حالتی که تعداد متغیرها و معادلات زیاد باشند حل تحلیلی این معادلات زمان زیادی می‌گیرد و در این حالت بهتر است به حل عددی روی بیاوریم.
• حل عددی برای به دست آوردن ویژه مقادیر ماتریس‌ها.
• الگوریتم حداقل مربعات که مجموعه‌ای از داده‌ها را با یک تابع نمایی، هایپربولیک یا ... فیت می‌کند و به این ترتیب اطلاعات بیشتری از یک سیستم به دست می‌آوریم.
• بهینه سازی عددی نیز برای بهینه کردن یک پارامتر در یک سیستم به کار می‌رود.

نکته‌ای که در محاسبات عددی مهم است میزان خطا در این محاسبات است که باید آن را نیز مد نظر قرار داد. در این مطلب چند روش محاسبه عددی در متلب را معرفی می‌کنیم و چند مثال از این روش‌ها را حل می‌کنیم.

حل عددی دستگاه معادلات خطی

در روش حل عددی دستگاه معادلات خطی می‌توان به روش‌های تکراری و غیرتکراری در نرم افزار متلب اشاره کرد که به گروه‌های زیر تقسیم‌بندی می‌شوند:

فیلم آموزش جبر خطی با متلب در فرادرس

کلیک کنید

روش‌های غیرتکراری

• کدهای متلب محاسبات عددی و روش گاوس ناقص
• کدهای متلب محاسبات عددی و روش گاوس یا پیووتینگ
• کدهای متلب محاسبات عددی و روش تجزیه LU

روش‌های تکراری

• کدهای متلب محاسبات عددی و روش ژاکوبی
• کدهای متلب محاسبات عددی و روش گاوس-سیدال

در دستگاه معادلات خطی برای حل یک دستگاه باید تعداد ضرایب معلوم و مجهول با هم برابر باشند. برای حل این دستگاه‌ها ابتدا آن‌ها را به فرمت ماتریسی در می‌آوریم. در حقیقت دستگاه معادلات زیر را در نظر بگیرید:

$$\begin{cases}a_{11}x_{1}+a_{12}x_{2}+\cdots+a_{1n}x_{n}=b_{1} \\ a_{21}x_{1}+a_{22}x_{2}+\cdots+a_{2n}x_{n}=b_{2}\\
\vdots \quad \quad \quad \quad\vdots \quad \quad \quad \quad \quad \quad\vdots\\
a_{m1}x_{1}+a_{m2}x_{2}+\cdots+a_{mn}x_{n}=b_{m}
\end{cases}$$

فرم ماتریسی این دستگاه را می‌توان به شکل خلاصه شده زیر نوشت:

$$Ax=b$$

که A، x و b به ترتیب ماتریس‌های $$m\times n$$، $$1 \times n$$ و $$1\times m$$ هستند و به صورت زیر نمایش داده می‌شوند:

$$\begin{bmatrix}
a_{11}& a_{12}&\cdots & a_{1n}\\
a_{21}&a_{22}&\cdots & a_{2n}\\
\vdots & \vdots & \cdots & \vdots\\
a_{m1}&a_{m2}& \cdots & a_{mn}
\end{bmatrix}\begin{bmatrix}
x_{1}\\
x_{2}\\
\vdots\\
x_{n}\\
\end{bmatrix}=\begin{bmatrix}
b_{1}\\
b_{2}\\
\vdots\\
b_{m}\\
\end{bmatrix}$$

برای حل این دستگاه معادلات خطی همان طور که گفتیم می‎‌توان از روش‌های مستقیم و غیرمستقیم استفاده کرد. در روش‌های مستقیم پاسخ اعداد ثابتی هستند که این اعداد شامل خطاهای حاصل از گردکردن می‌باشند و در روش‌های غیرمستقیم اصل بر حدس‌های پی در پی برای پاسخ به مسئله طرح ریزی شده است. در بین روش‌های مستقیم یا غیرتکراری حل عددی دستگاه معادلات خطی می‌توان به روش‎‌های زیر اشاره کرد:

1. کدهای متلب محاسبات عددی و روش معکوس ماتریس ضرایب
2. کدهای متلب محاسبات عددی و روش کرامر
3. کدهای متلب محاسبات عددی و روش‌های حذفی که شامل روش ساده (ناقص) گاوس و روش گاوس است.
4. کدهای متلب محاسبات عددی و روش تجزیه

نکته مهم در استفاده از این روش‌ها که در بالا بیان شد غیرصفر بودن دترمینان ماتریس ضرایب یا ماتریس A که در بالا معرفی کردیم است. در اینجا به بررسی روش حذف ساده (ناقص) گاوس می‌پردازیم.

روش ساده (ناقص) گاوس

یکی از محدودیت‌های این روش این است که تنها در صورت غیر صفر بودن قطر اصلی ماتریس ضرایب این روش کاربرد دارد. در این روش باید یک ماتریس افزوده تشکیل دهیم که مجموع ماتریس ضرایب و ماتریس معلومات مسئله است. در ماتریس افزوده یک ستون به ماتریس ضرایب اضافه می‌شود که همان ماتریس معلومات مسئله است. برای درک این روش یک مثال را با این روش حل می‌کنیم و سپس آن را در متلب نیز اجرا خواهیم کرد.

فرض کنید چهار مجهول داریم که فرم ماتریسی آن به صورت زیر است:

$$\begin{bmatrix}
1& 1&1 & 1\\
2&3&1 & 5\\
-1 & 1 & -5 & 3\\
3&1&7 & -2\\
\end{bmatrix}\begin{bmatrix}
x_{1}\\
x_{2}\\
x_{3}\\
x_{4}\\
\end{bmatrix}=\begin{bmatrix}
10\\
31\\
-2\\
18\\
\end{bmatrix}$$

می‌خواهیم این دستگاه چهار معادله و چهار مجهول را با روش ساده گاوس حل کنیم. روش کار به این صورت است که اگر $$a_{kk}$$ مخالف صفر باشد آن گاه این روش کاربرد دارد و ضرایب جدیدی تحت عنوان m برای $$k$$اُمین سطر از ماتریس افزوده به صورت زیر تعریف می‌کنیم:

$$m_{ik}=\frac{a_{ik}^{(k)}}{a_{lk}^{(k)}}\Rightarrow\begin{cases}a_{ij}^{(k+l)}=a_{ij}^{(k)}-m_{ik}a_{kj}^{(k)}\\ b_{i}^{(k+l)}=b_{i}^{(k)}-m_{ik}b_{k}^{(k)}\end{cases}$$

حال به ترتیب سطر اول عنصر $$a_{11}$$ را که برابر با 1 است در نظر می‌گیریم و $$m_{21}$$ را محاسبه می‌کنیم که برابر است با:

$$m_{21}=\frac{2}{1}=2$$

بدین ترتیب عناصر سطر دوم به صورت زیر تغییر می‌کنند:

$$\begin{cases}a_{21}^{(2)}=a_{21}^{(1)}-m_{21}a_{11}=0 \\ a_{22}^{(2)}=a_{22}^{(1)}-m_{21}a_{12}=1 \\
a_{23}^{(2)}=a_{23}^{(1)}-m_{21}a_{13}=-1\\
a_{24}^{(2)}=a_{24}^{(1)}-m_{21}a_{14}=3\\
b_{2}^{(2)}=b_{2}^{(1)}-m_{21}b_{1}=11
\end{cases}$$

همچنین برای سطر سوم $$m_{31}$$ برابر است با:

$$m_{31}=\frac{-1}{1}=-1$$

و عناصر سطر سوم به شکل زیر درمی‌‌آیند:

$$\begin{cases}a_{31}^{(2)}=a_{31}^{(1)}-m_{31}a_{11}=0 \\ a_{32}^{(2)}=a_{32}^{(1)}-m_{31}a_{12}=2 \\
a_{33}^{(2)}=a_{33}^{(1)}-m_{31}a_{13}=-4\\
a_{24}^{(2)}=a_{24}^{(1)}-m_{31}a_{14}=4\\
b_{3}^{(2)}=b_{3}^{(1)}-m_{31}b_{1}=8
\end{cases}$$

و در نهایت در این مرحله برای سطر آخر $$m_{41}$$ برابر است با:

$$m_{41}=\frac{3}{1}=3$$

و عناصر این سطر نیز به صورت زیر در می‌آیند:

$$\begin{cases}a_{41}^{(2)}=a_{41}^{(1)}-m_{41}a_{11}=0 \\ a_{42}^{(2)}=a_{42}^{(1)}-m_{41}a_{12}=-2 \\
a_{43}^{(2)}=a_{43}^{(1)}-m_{41}a_{13}=4\\
a_{44}^{(2)}=a_{44}^{(1)}-m_{41}a_{14}=-5\\
b_{4}^{(2)}=b_{4}^{(1)}-m_{41}b_{1}=-12
\end{cases}$$

بعد از انجام این مراحل ماتریس به شکل زیر در می‌آید و داریم:

$$\begin{bmatrix}
1& 1&1 & 1&10\\
0&1&-1 & 3&11\\
0 & 2 & -4 & 4&8\\
0&-2&4 & -5&-12\\
\end{bmatrix}$$

حالا همین مراحل را برای عنصر $$a_{22}$$ در نظر می‌گیریم و به سمت سطر‌های پایین یعنی سطرهای 3 و 4 حرکت می‌کنیم. در این حالت $$m_{32}$$ برابر است با:

$$m_{32}=\frac{2}{1}=2$$

بدین ترتیب عناصر سطر سوم به صورت زیر تغییر می‌کنند:

$$\begin{cases}a_{32}^{(3)}=a_{32}^{(2)}-m_{32}a_{22}=0 \\ a_{33}^{(3)}=a_{33}^{(2)}-m_{32}a_{23}=-2 \\
a_{34}^{(3)}=a_{34}^{(2)}-m_{32}a_{24}=-2\\
b_{3}^{(3)}=b_{3}^{(2)}-m_{32}b_{2}=-14
\end{cases}$$

همچنین عنصر $$m_{42}$$ برابر است با:

$$m_{42}=\frac{-2}{1}=-2$$

و عناصر سطر چهارم به صورت زیر می‌شوند:

$$\begin{cases}a_{42}^{(3)}=a_{42}^{(2)}-m_{42}a_{22}=0 \\ a_{43}^{(3)}=a_{43}^{(2)}-m_{42}a_{23}=2 \\
a_{44}^{(3)}=a_{44}^{(2)}-m_{42}a_{24}=1\\
b_{4}^{(3)}=b_{4}^{(2)}-m_{42}b_{2}=10
\end{cases}$$

و بدین ترتیب ماتریس افزوده به صورت زیر در می‌آید:

$$\begin{bmatrix}
1& 1&1 & 1&10\\
0&1&-1 & 3&11\\
0 & 0 & -2 & -2&-14\\
0&0&2 & 1&10\\
\end{bmatrix}$$

در مرحله بعد $$a_{33}$$ را در نظر می‌گیریم و $$m_{43}$$ را محاسبه می‌کنیم که برابر است با:

$$m_{43}=\frac{2}{-2}=-1$$

حال عناصر سطر چهارم به شکل زیر تغییر می‌کنند و داریم:

$$\begin{cases}a_{43}^{(4)}=a_{43}^{(3)}-m_{43}a_{33}=0 \\ a_{44}^{(4)}=a_{44}^{(3)}-m_{43}a_{34}=-1 \\
b_{4}^{(4)}=b_{4}^{(3)}-m_{43}b_{3}=-4
\end{cases}$$

و ماتریس در نهایت به شکل زیر در می‌آید:

$$\begin{bmatrix}
1& 1&1 & 1&10\\
0&1&-1 & 3&11\\
0 & 0 & -2 & -2&-14\\
0&0&0 & -1&-4\\
\end{bmatrix}$$

اگر این معادلات را به شکل جبری بنویسیم داریم:

$$\begin{aligned}
x_{I}+x_{2}+x_{3}+x_{4} &=10 \\
x_{2}-x_{3}+3 x_{4} &=11 \\
-2 x_{3}-2 x_{4} &=-14 \\
-x_{4} &=-4
\end{aligned}$$

که به راحتی مقادیر $$x_{3}$$، $$x_{2}$$، $$x_{1}$$ و $$x_{4}$$ قابل محاسبه است. همین روش را در متلب می‌توان به صورت زیر نوشت اگر ماتریس ضرایب و مقادیر را مانند مثال بالا به عنوان ورودی بدهیم هر دو حالت جواب‌های یکسانی را نشان می‌دهند.

1function ngaussel(A,b)
2n=length(b);
3x=zeros(n,1);
4fprintf('\n');
5disp('The augmented matrix is')
6augm=[A b]
7for k=1:n-1
8    for i=k+1:n
9        m=A(i,k)/A(k,k);
10        for j=k+1:n
11            A(i,j)=A(i,j)-m*A(k,j);
12        end
13        A(i,k)=m;
14        b(i)=b(i)-m*b(k);
15    end
16end
17x(n)=b(n)/A(n,n);
18for i=n-1:-1:1
19    S=b(i);
20    for j=i+1:n
21        S=S-A(i,j)*x(j);
22    end
23    x(i)=S/A(i,i);
24end
25%print the result
26fprintf('\n');
27disp('The transformed upper triangular augmented matrix C is')
28fprintf('\n');
29for i=1:n
30    for j=1:n
31        if(j<i)A(i,j)=0;end
32    end
33end
34C=[A b]
35fprintf('\n');
36disp('Back substitution gives the vector solution')
37x

برای اجرای دستور گاوس در متلب کد بالا را در یک اسکریپت با نام تابع یعنی ngaussel ذخیره کنید و در صفحه کار متلب ماتریس‌های A و b را تعریف کنید. با فراخوانی تابع این روش برای ماتریس ورودی شما اجرا شده و جواب نمایش داده می‌شود.

درون‌یابی و برازش

درون‌یابی و برازش برای پیدا کردن صورت کلی توابعی است که داده‌هایی از آن را داریم. این داده‌ها می‌توانند از آزمایش یا مشاهده به دست آمده باشند. در حقیقت n زوج داده داریم که به صورت $$(x_{i}, y_{i})$$ هستند و می‌خواهیم $$y(x)$$ را تخمین بزنیم.

این روش در حل عددی بسیار متداول و رایج است و باید دقت کرد که این درون‌یابی با عبور از منحنی یا Curve fitting متفاوت است. در حقیقت در روش عبور از منحنی، منحنی که از بیشترین داده‌ها عبور می‌کند در نظر گرفته می‌شود اما در درون یابی باید تابع به دست آمده شامل تمام نقاط و داده‌ها باشد. انواع روش‌های درون یابی و برازش را می‌توان به صورت زیر دسته‌بندی کرد:

• کدهای متلب محاسبات عددی و روش لاگرانژ
• کدهای متلب محاسبات عددی و روش تفاضل‌های تقسیم شونده پیشرونده نیوتن
• کدهای متلب محاسبات عددی و اسپلاین‌های مرتبه دو و سه
• کدهای متلب محاسبات عددی و برازش یا عبور از منحنی

در ادامه این مطلب با حل یک مثال روش لاگرانژ را توضیح می‌دهیم.

درون یابی به روش لاگرانژ

با حل یک مثال روش لاگرانژ را توضیح می‌دهیم. در این مثال می‌خواهیم برای نقاط $${(0,6), (1,0), (2,2)} چند جمله‌ای مناسب را پیدا کنیم. برای به دست آوردن چند جمله‌ای مورد نظر آن را به صورت زیر معرفی می‌کنیم:

$$f(x)=\sum_{i=0}^{n} L_{k}(x) y_{k}$$

ضرایب $$L_{k}(x)$$ به صورت زیر تعریف می‌شوند:

$$L_{k}(x)=\frac{\left(x-x_{0}\right) \cdots\left(x-x_{k-1}\right)\left(x-x_{k+1}\right) \cdots\left(x-x_{n}\right)}{\left(x_{k}-x_{0}\right) \cdots\left(x_{k}-x_{k-1}\right)\left(x_{k}-x_{k+1}\right) \cdots\left(x_{k}-x_{n}\right)}=\prod_{j=0 \atop j=2}^{n} \frac{x-x_{i}}{x_{k}-x_{i}}, k=0,1, \cdots, n$$

به این ترتیب برای داده‌های مثال فوق ضرایب $$L_{k}(x)$$ برابر با روابط زیر به دست می‌آیند:

$$L_{0}(x)=\frac{(x-1)(x-2)}{(0-1)(0-2)}=\frac{x^{2}-3 x+2}{+2}$$

$$L_{1}(x)=\frac{(x-0)(x-2)}{(1-0)(1-2)}=\frac{x^{2}-2 x}{-1}$$

$$L_{2}(x)=\frac{(x-0)(x-1)}{(2-0)(2-1)}=\frac{x^{2}-x}{2}$$

بدین ترتیب داریم:

$$\left\{\begin{array}{l}
f(x)=\sum_{k=0}^{2} L_{k}(x) y_{k}=L_{2}(x) y_{0}+L_{1}(x) y_{1}+L_{2}(x) y_{2} \\
f(x)=\sum_{k=0}^{2} L_{k}(x) y_{k}=6 \times\left(\frac{x^{2}-3 x+2}{2}\right)+0 \times\left(\frac{x^{2}-2 x}{-1}\right)+2 \times\left(\frac{x^{2}-x}{2}\right) \\
f(x)=4 x^{2}-10 x+6
\end{array}\right.$$

روش لاگرانژ در متلب به صورت زیر تعریف می‌شود:

1function lagrange(x,y,a)
2n=length(x);
3p=0;
4for k=1:n
5    b(k)=1;
6    d(k)=1;
7    for j=1:n
8        if j~=k
9            b(k)=b(k)*(x(k)-x(j));
10            d(k)=d(k)*(a-x(j));
11        end
12    end
13    c(k)=y(k)/b(k);
14    p=p+c(k)*d(k);
15end
16fprintf('\n p(a)=%10.6f',p)
17fprintf('\n')

با تعریف بردارهای x ،y و مقدار a تابع به روش لاگرانژ محاسبه می‌شود و مقدار تابع به ازای a در خروجی چاپ می‌شود.

مشتق و انتگرال‌گیری عددی

مشتق و انتگرال‌گیری عددی زمانی مورد استفاده قرار می‌گیرد که تابع در دسترس نباشد و یا با یک تابع پیچیده رو به رو باشیم. روش‌های مشتق گیری عددی شامل موارد زیر هستند:

فیلم آموزش روش های عددی ریشه یابی و حل معادلات و پیاده سازی در متلب MATLAB در فرادرس

کلیک کنید

• بسط تیلور
• فرمول ریچاردسون

روش‌های انتگرال گیری عددی نیز به موارد زیر تقسیم می‌شوند:

• روش نیوتن-کوتز
• روش مجموع ریمان
• روش قاعده ذوزنقه‌ای
• روش سیمپسون
• الگوریتم رومبرگ

در ادامه این مطلب روش بسط تیلور را مورد معرفی و بررسی قرار می‌دهیم.

بسط تیلور

در ابتدا بسط تیلور توابع $$f(x+h)$$ را حول $$h=0$$ بررسی می‌کنیم و داریم:

$$\begin{array}{l}
f(x+h)=f(x)+h f^{\prime}(x)+\frac{1}{2 !} h^{2} f^{\prime \prime}(x)+\frac{1}{3 !} h^{3} f^{\prime \prime}(x)+\cdots \\
\Rightarrow \frac{f(x+h)-f(x)}{h}=f^{\prime}(x)+\left[\frac{1}{2 !} h f^{\prime \prime}(x)\right]+\frac{1}{3 !} h^{2} f^{\prime \prime}(x)+\cdots
\end{array}$$

همین بسط برای تابعی به صورت $$f(x-h)$$ حول $$h=0$$ به صورت زیر است:

$$\begin{array}{l}
f(x-h)=f(x)-h f^{\prime}(x)+\frac{1}{2 !} h^{2} f^{\prime \prime}(x)-\frac{1}{3 !} h^{3} f^{\prime \prime}(x)+\cdots \\
\Rightarrow \frac{f(x)-f(x-h)}{h}=f^{\prime}(x)-\left[\frac{1}{2 !} h f^{\prime \prime}(x)\right]+\frac{1}{3 !} h^{2} f^{\prime \prime}(x)-\cdots
\end{array}$$

اگر دو رابطه بالا را با یکدیگر جمع کنیم جوابی به صورت زیر خواهیم داشت:

$$\frac{f(x+h)-f(x-h)}{2 h}=f^{\prime}(x)+\frac{1}{3 !} h^{2} f^{\prime \prime}(x)-\cdots$$\end{array}$$

بدین ترتیب می‌توان فرمول‌های تقریبی برای مشتق گیری مرتبه اول را به فرم‌های زیر تعریف کرد:

• تقریب اختلاف پیشرونده که به صورت زیر تعریف می‌شود:

$$f^\prime(x)\approx \frac{f(x+h)-f(x)}{h}$$

• همچنین تقریب اختلاف پسرو به شکل زیر است:

$$f^\prime(x)\approx \frac{f(x)-f(x-h)}{h}$$

• و تقریب اختلاف مرکزی نیز شکل زیر را می‌گیرد:

$$f^\prime(x)\approx \frac{f(x+h)-f(x-h)}{2h}$$

با محاسبات می‌توان دید که با توجه به این که $$h$$ همواره مقداری کوچکتر از یک دارد، خطای محاسباتی در روش مشتق گیری مرکزی از دو روش پیشرو و پسرو کمتر خواهد بود. همچنین مشتق دوم تابع $$f(x)$$ با توجه به روش مشتق گیری مرکزی حول نقطه $$h=0$$ به شکل زیر خواهد بود:

$$f^{\prime\prime}(x)\approx \frac{f(x+h)-2f(x)+f(x-h)}{h^{2}}$$

این روش را با یک مثال توضیح داده و آن را در متلب اجرا می‌کنیم. برای تابع $$f(x)=\sin(x)$$ مقدار تخمینی مشتق دوم آن را به ازای مقدار $$x=0.8\quad rad$$ و $$h=0.01$$ محاسبه کنید و خطای محاسبات را به دست آورید.

با استفاده از رابطه مشتق دوم از روش بسط تیلور مقدار تخمینی مشتق برابر با $$-0.6967009$$ خواهد شد. این در حالی است که مقدار واقعی از روابط تحلیلی مقداری برابر با $$-0.6967067$$ دارد. به این ترتیب مقدار خطا برابر با $$5.81e-6$$ خواهد بود.

این مثال را در متلب به صورت زیر می‌نویسیم:

1clear all
2h=0.01;
3t=(0:h:1);
4a=sin(2*pi*1*t);
5d2a=-((2*pi*1)^2)*sin(2*pi*1*t);
6d2f=diff(a,2)/(h*h);
7subplot(211),plot(t,d2a,'b',t(3:end),d2f,'r');
8subplot(212), plot(t(3:end),abs(d2a(3:end)-d2f));
9err=norm(d2f-d2a(3:end),2)/norm(d2a(3:end),2)

با اجرای این برنامه نمودارها به صورت زیر نمایش داده می‌شوند:

همچنین میزان ارور یا خطا نیز برابر با $$0.0622$$ در خروجی چاپ می‌شود.

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

در این بخش حل عددی معادلات دیفرانسیل معمولی را بررسی و معرفی می‌کنیم. به صورت کلی روش‌های حل عددی معادلات دیفرانسیل معمولی به صورت زیر معرفی می‌شوند:

فیلم آموزش معادلات دیفرانسیل با میپل Maple در فرادرس

کلیک کنید

• کدهای متلب محاسبات عددی و روش اویلر (بسط تیلور مرتبه اول)
• کدهای متلب محاسبات عددی و بسط‌های تیلور مرتبه‌های بالاتر از 1
• کدهای متلب محاسبات عددی و روش نقطه میانی
• کدهای متلب محاسبات عددی و روش رانگ-کوتای مرتبه 2 و 4
• کدهای متلب محاسبات عددی و روش آدامز-باشفورث
• کدهای متلب محاسبات عددی و روش آدامز-مولتون

در ادامه روش اویلر را معرفی کرده و با استفاده از یک مثال و بررسی آن به صورت دستی و در متلب این روش را توضیح می‌دهیم.

روش اویلر

ساده‌ترین شیوه تحلیلی از میان همه روش‌های محاسباتی در حل معادلات دیفرانسیل معمولی روش اویلر است. معادله‌ای مانند زیر را در نظر بگیرید که مقدار اولیه آن برابر با $$y(a)=y_{0}$$ است و در بازه تقریبی $$[a,b]$$ تعریف می‌شود و داریم:

$$\frac{dy}{dt}=f(t,y)$$

برای شروع از روش اویلر ابتدا بازه را به $$N$$ زیر بازه مساوی تقسیم کرده و نقاط شبکه‌ای به صورت $$t_{i}=a+ih$$ به ازای $$i=0,1,2,...,N$$ تعریف می‌کنیم که $$h$$ یا گام حرکت برابر با $$h=\frac{b-a}{N}$$ است. به این ترتیب بسط تیلور تابع $$y(t)$$ را حول $$t_{i}$$ می‌نویسیم و داریم:

$$y(t_{i+1})=y(t_{i}+h)=y(t_{i})+hy^{\prime}(t_{i})+\frac{h^{2}}{2}y^{\prime\prime}(\zeta)$$

که $$\zeta$$ مقداری بین $$t_{i}$$ و $$t_{i}+h$$ دارد. با صرف نظر از جملات دوم به بعد داریم:

$$y\left(t_{i}+h\right) \approx y\left(t_{i}\right)+h y^{\prime}\left(t_{i}\right)=y(t_{i})+hf(t_{i},y(t_{i}))$$

و با نمایش $$y(t_{i})$$ به صورت $$y_{i}$$ معادله به صورت زیر در می‌آید:

$$y_{i+1}=y_{i}+hf(t_{i},y_{i})\quad\quad i=0,1,\cdots,N$$

معادله فوق به نام روش اویلر مشهور است. برای درک بیشتر موضوع مثال زیر را حل می‌کنیم و داریم:

$$y^{\prime}=-6y$$

می‌دانیم که جواب این معادله دیفرانسیل برابر با $$y=exp(-6x)$$ است. بر اساس معادله اویلر داریم:

$$y_{i+1}=y_{i}+h(-6y_{i})$$

با نوشتن این تابع در متلب و برای گام‌هایی به طول $$0.1$$ و 10 گام و شرط اولیه $$y(1)=1$$ داریم:

1h=0.1; % step's size
2N=10; % number of steps
3y(1)=1;
4for n=1:N
5y(n+1)= y(n)+h*(-6*y(n));
6x(n+1)=n*h;
7end
8plot(x,y)

با اجرای این مثال نموداری به شکل زیر در خروجی ظاهر می‌شود:

با تغییر طول گام و تعداد گام‌ها شکل نمودار به حالت منحنی یا curve نزدیکتر می‌شود و در این حالت خروجی به شکل زیر نمایش داده می‌شود:

حال اگر حل عددی اویلر را با مشتق مستقیم تابع در همان بازه رسم کنیم و با یکدیگر مقایسه کنیم نتایج به صورت زیر نمایش داده می‌شوند:

1h1=0.01; % step's size
2N1=100; % number of steps
3y(1)=1;
4for n=1:N1
5y(n+1)= y(n)+h1*(-6*y(n));
6x(n+1)=n*h1;
7end
8plot(x,y,'r')
9hold on
10x=0:0.001:1;
11y=exp(-6.*x);
12plot(x,y,'g'

حل عددی معادلات دیفرانسیل غیرخطی

معادلات ديفرانسيل آهنگ تغييرات يک کميت را نمايش می‌دهند. در مطلب حل دستگاه معادلات ديفرانسيل غير خطي در متلب | گام به گام در مورد حل تحليلی اين معادلات در متلب صحبت کرديم در اين بخش از اين مبحث حل عددی معادلات ديفرانسيل غيرخطی در متلب را بررسی و معرفی مي‌کنيم.

به صورت کلی معادلات ديفرانسيل ترکيبی از توابع و مشتقات آن‌ها هستند و بسياری از مسائل فيزيکی و مهندسی را مدل‌سازی می‌کنند. اين معادلات شامل شرايط اوليه يا شرايط مرزی هستند که با توجه با اين شرايط اين مسائل و معادلات مورد بررسی قرار می‌گيرند. انواع روش‌های عددی معادلات ديفرانسيل غير خطی به سه دسته کلی تقسيم می‌شوند:

1. روش‌های ريشه يابی عددی دامنه محدود يا بسته
2. روش‌های ريشه يابی عددی دامنه نامحدود يا باز
3. حل عددی دستگاه معادلات غيرخطی به روش نيوتن

روش‌های ريشه يابی عددی دامنه محدود يا بسته را به دو روش نصف کردن يا Bisection و روش نابجايی يا موقعيت خطا می‌توان انجام داد. همچنين روش‌های ريشه يابی عددی دامنه نامحدود يا باز نيز به سه روش سکانت، نيوتون-رافسون و نقطه ثابت صورت می‌گيرد. در ادامه اين بخش روش نصف کردن يا Bisection را معرفی می‌کنيم و آن را با حل يک مثال و بررسی آن در متلب توضیح خواهیم داد.

روش نصف کردن يا Bisection

روش نصف کردن براي هر تابع پيوسته‌ای مورد استفاده قرار می‌گيرد. برای استفاده از اين روش بايد تابع را در بازه‌ای در نظر بگيريم که ضرب مقادير تابع در آن بازه کوچکتر از صفر و مشتق تابع مخالف صفر باشد. بدين ترتيب بازه تابع را نصف کرده و مقدار تابع را در آن نقطه به دست می‌آوريم و با توجه به شرايط آن را جايگزين مقدار اوليه يا نهايی بازه می‌کنيم يعنی داريم:

$$y=f(x)\quad [a,b]$$
$$c=\frac{(a+b)}{2}$$
$$\begin{cases}if\quad f(a).f(c)<0\rightarrow b=c\\ if \quad f(b).f(c)<0 \rightarrow a=c\end{cases}$$

اين کار را تا جايي ادامه می‌دهيم تا مقدار $$f(c)$$ به صفر نزديک شود. اين روش را برای پيدا کردن ريشه تابع $$x^2-3$$ در بازه $$[1,2]$$ استفاده می‌کنيم. بدين ترتيب داريم:

1function c = bisectionMethod(f,a,b,error)%f=@(x)x^2-3; a=1; b=2; (ensure change of sign between a and b) error=1e-4
2c=(a+b)/2;
3while abs(f(c))>error
4if f(c)<0&&f(a)<0
5a=c;
6else
7b=c;
8end
9c=(a+b)/2;
10end

با تعريف تابع و بازه آن و دقت نزديک بودن جواب به صفر و اجرای این کد ريشه اين تابع در خروجی نمايش داده می‌شود.

جمع بندی

این مطلب به محاسبات عددی در متلب اختصاص داده شده است. در حقیقت روش‌ها و کدهای متلب محاسبات عددی بسیار گسترده هستند و روش‌های متفاوتی برای آن‌ها وجود دارد. در این مطلب این روش‌ها را به صورت فهرست وار معرفی کردیم و از هر کاربرد آن یک نمونه را با مثال و بررسی مثال در متلب توضیح دادیم.