ریاضی , علوم پایه 1174 بازدید

در آموزش‌های قبلی از مجموعه مباحث ریاضیات مجله فرادرس، با معادلات دیفرانسیل آشنا شدیم و درباره برخی از روش‌های تحلیلی حل انواع مختلف آن‌ها بحث کردیم. علاوه بر روش‌های تحلیلی، با استفاده از روش‌های عددی نیز می‌توان جواب معادلات دیفرانسیل را محاسبه کرد. روش اویلر، روشی عددی برای حل معادلات دیفرانسیل است. در این آموزش، روش اویلر (Euler Method) یا رونگه-کوتای مرتبه اول (First Order Runge–Kutta) را معرفی می‌کنیم.

روش اویلر

روش اویلر یا رونگه کوتای مرتبه اول، ریاضیات بسیار ساده‌ای دارد و می‌توان آن را به سادگی درک کرد. فرض کنید معادله دیفرانسیل مرتبه اول زیر و شرایط اولیه معلوم $$ y (t _ 0 ) = y _ 0 $$ را داریم:

$$ \large { { d y ( t ) } \over { d t } } = y ^ \prime \left ( t \right ) = f ( y ( t ) , t ) $$

تقریب مشتق را به صورت زیر می‌نویسیم:

$$ \large k _ 1 = y ^ \prime \left ( t \right ) = f ( y ^ * ( t ) , t ) $$

بسط تیلور $$ y (t) $$ را حول $$ t _ 0$$ و با فرض گام زمانی $$h$$ می‌نویسیم:

$$ \large y ( { t _ 0 } + h ) = y ( { t _ 0 } ) + y ^ \prime ( {t _ 0 } ) h + y ^ {\prime \prime } ( { t _ 0 } ) { { { h ^ 2 } } \over 2 } + \cdots $$

در عبارت بالا، از جملات بعد از جمله خطی (جملاتی با توان دو و بالاتر) صرفنظر می‌کنیم. بنابراین، عبارت بالا را به صورت زیر می‌نویسیم:

$$ \large y ( { t _ 0 } + h ) \approx y ( { t _ 0 } ) + y ^ \prime ( { t _ 0 } ) h $$

در نتیجه، می‌توانیم جواب تقریبی را در گام زمانی بعدی به صورت زیر به دست آوریم:

$$ \large y ^ * ( { t _ 0 } + h ) = y ^ * ( { t _ 0 } ) + { k _ 1 } h $$

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

معادله دیفرانسیل خطی همگن

فرض کنید می‌خواهیم با استفاده از کامپیوتر، جواب معادله دیفرانسیل زیر را به دست آوریم:

$$ \large { { d y ( t ) } \over { d t } } + 2 y ( t ) = 0 $$

یا

$$ \large { { d y ( t ) } \over { d t } } = – 2 y ( t ) $$

که شرط اولیه آن، $$y (0 ) = 3 $$ است. در این حالت، جواب دقیق به صورت $$ y ( t) = 3 e ^ { – 2 t } $$ برای $$ t \ge 0 $$ است. این جواب در شکل زیر رسم شده است. با دانستن جواب دقیق، می‌توانیم جواب تقریبی خود را با آن مقایسه کنیم.

نمودار جواب معادله

راه‌های مختلفی برای به دست آوردن جواب تقریبی وجود دارد. در اینجا از بسط تیلور $$ y ( t) $$ حول $$t=0$$ استفاده می‌کنیم (در حالت کلی، این بسط حول نقطه $$ t = t  _ 0 $$ است). بسط تیلور تابع به صورت زیر است:‌

$$ \large y ( t ) = y ( 0 ) + y ^ \prime ( 0 ) t + y ^ { \prime \prime } ( 0 ) { { { t ^ 2 } } \over 2 } + \cdots $$

اکنون جواب را به گام زمانی کوچک $$h$$ بعد از $$t =0 $$ محدود کرده و بسط تیلور را تا مشتق اول در نظر می‌گیریم. بنابراین، داریم:

$$ \large \eqalign {
& y ( h ) = y ( 0 ) + y^ \prime ( 0 ) h + y ^ {\prime \prime } ( 0 ){ { { h ^ 2 } } \over 2 } + \cdots \cr
& y ( h ) \approx y ( 0 ) + y ^ \prime ( 0 ) h \cr } $$

در ادامه، تقریب جواب را $$ y ^ * (h) $$ و مشتق آن را $$ y ^ \prime ( 0 ) =  k _ 1 $$ می‌نامیم.

$$ \large y ^ * ( h ) = y ( 0 ) + { k _ 1 } h $$

شکل زیر، نمودار جواب را برای $$ h = 0.2 $$ نشان می‌دهد.

نمودار جواب تقریبی

نمودار سمت راست، ناحیه سایه زده شده نمودار سمت چپ را با بزرگنمایی نشان می‌دهد. تقریب تیلور مرتبه اول، یک خط راست است که از مقدار اولیه شروع شده و شیب آن $$-6$$ است (یعنی $$ k _ 1 = y ^ \prime ( 0 ) = -2 y (0) = – 6 $$).

تقریب در $$ t= h = 0.2 $$ برابر با حاصل‌ جمع مقدار اولیه و ضرب شیب در گام زمانی $$h$$ است:

$$ \large y ^ * ( h ) = y ^ * ( 0 .2 ) = y ( 0 ) + { k _ 1} y ^ \prime (0) = 1.8$$

واضح است که هرچه $$h$$ کوچکتر باشد، تقریب بهتر و دقیق‌تر خواهد بود.

برای پیدا کردن تقریب گام بعدی $$ y ^ * ( 2 h ) $$، به سادگی فرایند را تکرار کرده و از مقدار قبلی $$ y ^ * ( h ) $$ برای تخمین مشتق در زمان $$h$$ استفاده می‌کنیم (مقدار دقیق $$ y ( h ) $$ را نمی‌دانیم؛ بنابراین، فقط می‌توانیم مشتق ($$ k _ 1 $$) را محاسبه کنیم). ترتیب انجام مراحل به صورت زیر است:

عبارت دقیق مشتق:     $$ \large y ^ \prime ( t ) = – 2 y ( t ) $$

تقریب مشتق:    $$ \large { k _ 1 } = – 2 y ^ * ( h ) $$

سری تیلور حول $$t = h $$:     $$ \large y ( 2 h ) = y ( h ) + y ^ \prime ( h ) h + y ^ { \prime \prime }( h ) { { { h ^ 2 } } \over 2 } + \cdots $$

بسط تیلور بریده:     $$ \large y ( 2 h ) \approx y ( h ) + y ^ \prime ( 0 ) h $$

جواب تقریبی:    $$ \large y ^ * ( 2 h ) = y ^ * ( h ) + { k _ 1 } h $$

حرکت به سمت جلو از زمان $$  t _0 $$ تا $$ t_0 + h $$ در حالت کلی به صورت زیر خواهد بود:

عبارت دقیق مشتق در $$ t = t_0$$:     $$ \large y ^ \prime ( t_0 ) = – 2 y ( t_0 ) $$

تقریب مشتق با استفاده از تقریب قبلی $$ y (t ) $$:    $$ \large { k _ 1 } = – 2 y ^ * ( t_0 ) $$

سری تیلور حول $$t = t _ 0 $$:     $$ \large y ( t_0 + h ) = y ( t _ 0 ) + y ^ \prime ( t _ 0 ) h + y ^ { \prime \prime }( t _ 0 ) { { { h ^ 2 } } \over 2 } + \cdots $$

بسط تیلور بریده:     $$ \large y ( t _ 0 + h ) \approx y ( t_ 0) + y ^ \prime ( t_ 0 ) h $$

جواب تقریبی:    $$ \large y ^ * ( t _ 0 +h ) = y ^ * ( t _ 0 ) + { k _ 1 } h $$

پیاده‌سازی در متلب

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

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

مقایسه جواب دقیق و تقریبی

خلاصه الگوریتم

برای به دست آوردن جواب تقریبی معادله دیفرانسیلِ

$$ \large { { d y ( t ) } \over { d t } } = f ( y ( t ) , t ) $$

به روش رونگه-کوتای مرتبه اول، از نقطه $$t _0 $$ شوع کرده و با گام زمانی $$ h $$، مراحل زیر را به ترتیب طی می‌کنیم:

تقریب مشتق:    $$ \large { k _ 1 } = f ( { y ^ * } { \rm { ( } } { { \rm{ t }} _ 0 } ) , { t _0 } ) $$

تقریب جواب در گام زمانی بعدی:    $$ \large { y ^ * } ( { t _ 0 } + h ) = { y ^ * } ( { t _ 0 } ) + { k _ 1 } h $$

داشتن مقدار اولیه برای شروع الگوریتم، ضروری است.

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

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

$$ \large y ^ \prime ( t ) + 2 y ( t ) = \cos ( 4 t ) $$

یا

$$ \large y ^ \prime ( t ) = – 2 y ( t ) + \cos ( 4 t ) $$

برای به دست آوردن جواب تقریبی، مراحل زیر را طی می‌کنیم:

عبارت دقیق مشتق در $$ t = t _0 $$:    $$ \large y ^ \prime ( { t _ 0 } ) = – 2 y ( { t _ 0 } ) + \cos ( 4 { t _ 0 } ) $$

تقریب مشتق:    $$ \large { k _ 1 } = – 2 y ^ * ( { t _ 0 } ) + \cos ( 4 { t _ 0 } ) $$

بسط تیلور حول $$ t = t _0 $$:    $$ \large y ( { t _ 0 } + h ) = y ( { t _ 0 } ) + y ^ \prime ( { t _ 0 } ) h + y ^ { \prime \prime } ( { t _ 0 } ) { { { h ^ 2 } } \over 2 } + \cdots $$

بسط تیلور بریده:    $$ \large y ( { t _ 0 } + h ) \approx y ( { t _ 0 } ) + y ^ \prime ( { t _ 0 } ) h $$

جواب تقریبی:    $$ \large y ^ * ( { t _ 0 } + h ) = y ^ * ( { t _ 0 } ) + { k _ 1 } h $$

جواب دقیق این مثال $$ y ( t ) = 2.9 { e ^ { – 2 { \kern 1 pt } t } } + 0.1 \cos ( 4 t ) + 0.2 \sin ( 4 t ) $$ است.

پیاده‌سازی در متلب

برای پیاده‌سازی محاسبات بالا، می‌توانیم از نرم‌افزار متلب استفاده کنیم. دقت کنید که در محاسبات، مقدار $$k_1$$ تغییر می‌کند.

شکل زیر، نمودار جواب دقیق و جواب‌های تقریبی را برای مقادیر متفاوت $$h$$ نشان می‌دهد. مانند مثال قبل، می‌بینیم که با کاهش $$h$$، جواب تقریبی به جواب دقیق نزدیک می‌شود.

جواب تقریبی و دقیق

معادله دیفرانسیل مرتبه اول غیرخطی

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

$$ \large { { d y ( t ) } \over { d t } } – y ( t ) ( 1 – 2 t ) = 0 , \quad \quad \quad \quad y ( 0 ) = 1 $$

جواب دقیق این معادله به صورت زیر است:

$$ \large y ( t ) = { e ^ { t – { t ^ 2 } } } $$

برای به دست آوردن جواب تقریبی، عبارت $$ y ^ \prime  (t) $$ را نوشته و از آن برای یافتن $$ k _ 1 $$ استفاده می‌کنیم. محاسبات را به صورت زیر انجام می‌دهیم:

عبارت دقیق مشتق در $$ t = t _ 0 $$:    $$ \large y ^ \prime ( { t _ 0 } ) = y ( { t _ 0 } ) ( 1 – 2 { t _ 0 } ) $$

تقریب مشتق:    $$ \large { k _ 1 } = y ^ * ( { t _ 0 } ) ( 1 – 2 { t _ 0 } ) $$

جواب تقریبی:    $$ \large y ^ * ( { t _ 0 } + h ) = y ^ * ( { t _ 0 } ) + { k _ 1 } h $$

پیاده‌سازی در متلب

مانند مثال‌های قبلی، محاسبات اخیر را می‌توانیم در متلب پیاده‌سازی کنیم. برنامه مربوط به مثال اخیر به صورت زیر است.

نمودار جواب دقیق و جواب‌های تقریبی به ازای مقادیر مختلف $$h $$ در شکل زیر نشان داده شده‌اند.

منحنی جواب تقریبی و دقیق

معادله دیفرانسیل مرتبه بالا

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

$$ \large \displaylines {
{ { { d ^ 3 } y ( t ) } \over { d t } } + 4 { { { d ^ 2 } y ( t ) } \over { d t } } + 6 { { d y ( t ) } \over { d t } } + 4 y ( t ) = \gamma ( t ) \quad \quad \quad  \cr
{ \left . { { {{ d ^ 2 } y ( t ) } \over { d t } } } \right | _ { t = { 0 ^ + } } } = 0 , \quad \quad { \left . { { { d y ( t ) } \over { d t } } } \right | _ { t = { 0 ^ + } } } = – 1 , \quad \quad y ( { 0 ^ + } ) = 0 \quad \cr } $$

که در آن، $$\gamma (t) $$ تابع پله واحد است.

جواب دقیق این معادله برای $$t \ge 0 $$، $$ y ( t ) = { 1 \over 4 } + { { \rm { e } } ^ { – t } } { \mkern 1mu } \left ( { \cos ( t ) – { 5 \over 2 } \sin ( t ) } \right ) – { 5 \over 4 } { { \rm { e } } ^ { – 2 { \kern 1pt } t } } $$ است.

با تعریف متغیرهای جدید $$q _ 1$$، $$ q _ 2$$ و $$q _ 3$$ معادله دیفرانسیل فوق را با سه معادله مرتبه اول توصیف می‌کنیم:

$$ \large \eqalign {
{ q _ 1 } ( t ) & = y ( t ) \cr
{ q _ 2 } ( t ) & = y ^ \prime ( t ) \cr
{ q _ 3 } ( t ) & = y ^ { \prime \prime } ( t ) \cr
\cr
{ q _ 1 } ^ \prime ( t ) & = y ^ \prime ( t ) & = { q _ 2 } (t ) \cr
{ q _ 2 } ^ \prime ( t ) & = y ^ { \prime \prime } ( t ) & = { q _ 3 } ( t ) \cr
{ q _ 3 } ^ \prime ( t ) & = y ^ {\prime \prime\prime } ( t ) & = \gamma ( t ) – 4 y ^ { \prime \prime } ( t ) – 6 y ^ \prime ( t ) – 4 y ( t ) \cr
& & = \gamma ( t ) – 4 { q _ 3 } ( t ) – 6 { q _ 2 } ( t ) – 4 { q _ 1 } ( t ) \cr } $$

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

$$ \large \eqalign {
{ \bf { q } } ^ \prime ( t ) & = \left [ { \matrix {
{ { q _ 1 } ^ \prime ( t ) } \cr
{ { q _ 2 } ^ \prime ( t ) } \cr
{ { q _ 3 } ^ \prime ( t ) } \cr } } \right ] = \left [ { \matrix {
0 & 1 & 0 \cr
0 & 0 & 1 \cr
{ – 4 } & { – 6 } & { – 4 } \cr } } \right ] \left [ { \matrix {
{ { q _ 1 } ( t ) } \cr
{ { q _ 2 } ( t ) } \cr
{ { q _ 3 } ( t ) } \cr } } \right ] + \left [ { \matrix {
0 \cr
0 \cr
1 \cr } } \right ] \gamma ( t ) \cr
{ \bf { q } } ^ \prime ( t ) & = { \bf { A q } } ( t ) + { \bf { B \gamma } } ( t ) \cr
{ \bf { A } } & = \left [ { \matrix {
0 & 1 & 0 \cr
0 & 0 & 1 \cr
{ – 4 } & { – 6 } & { – 4 } \cr } } \right ] , \quad \quad \quad { \bf { B } } = \left [ { \matrix {
0 \cr
0 \cr
1 \cr
} } \right ] \cr } $$

اکنون، مسئله مشابه موارد قبلی است؛ با این تفاوت که متغیرها و مقادیر $$k_1$$، در این حالت بردار هستند:

عبارت دقیق مشتق در $$ t = t_0 $$:    $$ \large { \bf { q } } ^ \prime ( { t _ 0 } ) = { \bf { A q } } ( { t _ 0 } ) + { \bf { B \gamma } } ( { t _ 0 } ) $$

تقریب مشتق:    $$ \large { { \bf { k } } _ 1 } = { \bf { A } } { { \bf { q } } ^ * } ( { t _ 0 } ) + { \bf { B \gamma } } ( { t _ 0 } ) $$

بسط تیلور حول $$ t = t _0$$:    $$ \large { \bf { q } } ( { t _ 0 } + h ) = { \bf { q } } ( { t _ 0 } ) + { \bf { q } } ^ \prime ( { t _ 0 } ) h + { \bf { q } } ^ {\prime \prime } ( { t _ 0 } ) { { { h ^ 2 } } \over 2 } + \cdots $$

بسط تیلور بریده:    $$ \large { \bf { q } } ( { t _ 0} + h ) \approx { \bf { q } } ( { t _ 0 } ) + { \bf { q } } ^ \prime ( { t _ 0 } ) h $$

جواب تقریبی:    $$ \large { \bf { q } } ^ * ( { t _ 0 } + h ) = { \bf { q } } ^ * ( { t _ 0 } ) + { { \bf { k } } _ 1 } h $$

پیاده‌سازی در متلب

مانند مثال‌های قبلی، محاسبات اخیر را می‌توانیم در متلب پیاده‌سازی کنیم. برنامه مربوط به این مثال به صورت زیر است:

نمودار جواب دقیق و جواب‌های تقریبی به ازای مقادیر مختلف $$h $$ در شکل زیر نشان داده شده‌اند.

جواب معادله مرتبه سوم

اگر این مطلب برایتان مفید بوده است، آموزش‌های زیر نیز به شما پیشنهاد می‌شوند:

^^

به عنوان حامی، استارتاپ، محصول و خدمات خود را در انتهای مطالب مرتبط مجله فرادرس معرفی کنید.

telegram
twitter

سید سراج حمیدی

«سید سراج حمیدی» دانش‌آموخته مهندسی برق است. فعالیت‌های کاری و پژوهشی او در زمینه سیستم‌های فتوولتائیک و کاربردهای کنترل در قدرت بوده و، در حال حاضر، آموزش‌های مهندسی برق و ریاضیات مجله فرادرس را می‌نویسد.

بر اساس رای 2 نفر

آیا این مطلب برای شما مفید بود؟

نظر شما چیست؟

نشانی ایمیل شما منتشر نخواهد شد. بخش‌های موردنیاز علامت‌گذاری شده‌اند *