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

function second_oder_ode
 
% SOLVE  d2x/dt2+5 dx/dt - 4 x = sin(10 t)
% initial conditions: x(0) = 0, x'(0)=0
 
t=0:0.001:3;   % time scale
 
initial_x    = 0;
initial_dxdt = 0;
 
[t,x]=ode45( @func, t, [initial_x initial_dxdt] );
 
plot(t,x(:,1));
xlabel('t'); ylabel('x');
 
function dxdt=func(t,x)
    dxdt_1 = x(2);
    dxdt_2 = -5*x(2) + 4*x(1) + sin(10*t);
 
    dxdt=[dxdt_1; dxdt_2];
end
end

نکته ی مهمی که باید به آن توجه کنید این است که وقتی معادله ای با مرتبه بالاتر حل می کنید شرایط اولیه برای مشتق های بالاتر(تا یکی قبل از بالاترین درجه) را باید بدانید.
نمودار زیر نشان دهنده جواب معادله دیفرانسیل بالا است

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

در این شکل m جرم جسم، c سرعت گیر(اصطکاک) و k قدرت فنر است. Wf هم فرکانس نیروی وارده است. با تغییر در این مقادیر میتوانید حالت های مختلف نمودار را ببینید.


@MatlabTips
#For_intermediate, #For_advance
#Differential_equations, #Dynamical_systems, #ode45, #Second_order_ode
نویسنده:(#A_2)

#Video_GUI
ساخت یک شمارنده که با هر بار کلیک رنگ آن بصورت تصادفی تغییر می کند.👆

🔵توابع اینلاین🔵
سطح پیچیدگی:🌕🌑🌑🌑🌑
پیش‌نیاز:(آشنایی مقدماتی با متلب)
برای: مقدماتی

فرض کنید تابع زیر را که بر حسب t است می خواهیم رسم کنیم:

h(t) = cos(8t) + cos(9t);

می توان تابع h را در پنجره فرمان به عنوان یک تابع اینلاین بصورت زیر نمایش داد:

h = inline('cos(8*t) + cos(9*t)')

حالا میخواهیم با نوشتن عبارات مناسب این تابع را رسم کنیم:

x = 0 : pi/40 : 6*pi;
plot(x, h(x)), grid

✔️توجه:
1⃣متغیر tدر تعریف اینلاین تابع h آرگومان ورودی است. این صرفا یک متغیر مشقی است و تنها ورودی یک تابع است. شما می توانید هر طور که مایلید نام متغیرتان را اسم گذاری کنید.

2⃣شما می توانید توابعی با بیش از یک متغیر با فرمت اینلاین تعریف کنید:

f = inline('x.ˆ2 + y.ˆ2', 'x', 'y' );
f(1, 2)
ans =
5

@MatlabTips
#Inline
#For_intermediate
نویسنده: (#A_1)

🔵نقطه تعادل🔵
سطح پیچیدگی:🌕🌕🌕🌕🌑
پیش‌نیاز:(کد نویسی متوسط)
برای:(متوسط و پیشرفته)

درک معادلات دیفرانسیل مخصوصا زمانی که مشتق های بالاتر در آن ظاهر می شود دشوار است. حل کردن برخی از معادلات دیفرانسیل مانند مساله سه جسم ذهن بسیاری را به خود مشغول کرده بود. یک روش که به فهم مساله کمک می کند توسط هنری پوانکاره ابداع شد. در این روش به جای حل معادله بر حسب زمان معادله را در مختصات جدیدی متشکل از متغیرهای مساله به عنوان محورها حل می شود. هر نقطه در این فضای جدید متناظر با یک حالت (state) است زیرا در آن حالت سیستم یک وضعیت مشخص دارد و کاملا می توان به صورت یکتا آن را از حالات دیگر تمیز داد.
برای فهم موضوع یک مثال کلاسیک و ساده را در نظر می گیریم. اگر معادله پاندول را در نظر بگیرید به صورت زیر است:

که در آن ثابت C میزان اصطکاک و k متناظر با عکس طول پاندول است. برای حل این معادله دیفرانسیل کافی است تابع آن را به صورت زیر بنویسیم:

function dx = pendulum(t,x)
global c k
 
dx = zeros(2,1);
dx(1) = x(2);
dx(2) = -c*x(2)-k*sin(x(1));

اگر توجه کنید مقدار dx(1) زاویه پاندول (مکان آن) و dx(2) سرعت آن را مشخص می کند(سرعت مشتق مکان است). حالا با فرض c=1 و k=1 و در بازه زمانی صفر تا 20 ثانیه آن را حل می کنیم:

global c k
c = 0; k = 1;
[t,x] = ode45(@pendulum,[0 20],[1 1]);
plot(x(:,1),x(:,2))

که به صورت زیر تبدیل می شود:

اگر قدری به این شکل دقت کنید متوجه می شوید که بر اساس شرایط اولیه پاندول از نقطه 1 با سرعت صفر(نقطه 1 و 0) شروع کرده است(پاندول با قدری زاویه رها شده است) و اندک اندک از زاویه و سرعتش کاسته شده تا در نقطه صفر و صفر متوقف شده است( صفر اولی مربوط به مختصات اول و زاویه آن است که صفر شده و دومی سرعت آن است که به صفر رسیده است)
حال فرض کنید که اصطکاک نداریم و همچنین سرعت اولیه را برابر با 1.5 در نظر بگیرید(در کد بالا کافی است شرایط اولیه را از [1,1] به [1,1.5] تغییر دهید و c=0) در این صورت معادله به صورت زیر می شود:

این بار مقدار سرعت اولیه را به قدری افزایش می دهیم تا کاملا بچرخد(در کد بالا به به جای [1,1] بگذارید [1,3] و c=0)

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

c = 1; k = 1;
intval = linspace(-1,1,10);
 
figure;hold on;
for ii=1:10,
for jj=1:10,
        [t,x] = ode45(@pendulum,[0 20],[intval(ii) intval(jj)]);
        plot(x(:,1),x(:,2));
end
end

نمودار زیر نشان می دهد که تمامی خطوط به نقطه صفر و صفر ختم می شوند.