24-08-2020, 12:04 AM
نابراین زمانی که از تفاضل مرکزی برای گسستهسازی معادله نفوذ جابهجایی یک بعدی استفاده میشود، سطر i ام تنها شامل سه ضریب بالا است و باقی ضرایب در سطر i ام برابر با صفر در نظر گرفته میشوند. به صورت کلی میتوان بیان کرد که تعداد عبارات صفر در هر سطر به تعداد نقاطی بستگی دارد که با استفاده از آنها تقریب تفاضل محدود را بیان میکنیم.
رابطه زیر که در بالا نیز مورد مطالعه قرار گرفت را دوباره در نظر بگیرید.
رابطه 29
این رابطه به صورت نیمه گسسته است. زیرا گسسته سازی تنها در مکان انجام شده و در زمان از تقریب تفاضل محدود استفاده نشده است. برای اینکه معادله نفوذ و جابهجایی ارائه شده را به صورت عددی حل کنیم، لازم است که گسسته سازی آن را به صورت کامل انجام دهیم. برای این منظور، با استفاده از روشهایی که در معادلات دیفرانسیل معمولی بیان کردیم سمت چپ معادله فوق را گسسته میکنیم. برای مثال در صورتی که از روش اویلر پیشرو ساده استفاده کنیم، معادله فوق به شکل زیر در میآید.
رابطه 30
استفاده از روش تفاضل محدود مرکزی برای بیان مشتقات مکانی و روش اویلر پیشرو برای بیان مشتقات زمانی یکی از رایجترین، پرکاربردترین و معروفترین روشهای موجود در تفاضل محدود ودینامیک سیالات محاسباتی[url=https://blog.faradars.org/computational-fluid-dynamics/][/url] را نتیجه میدهد. این روش به روش «زمان پیش رو و مکان مرکزی» (Forward Time Central Space) معروف است. در محاسبات عددی و CFD این روش را به صورت اختصاری با نماد FTCS نمایش میدهند.
توجه کنید که این رابطه یک رابطه صریح است و نیاز به بیان پارامتر A به صورت صریح نیست. نکته دیگر این است که برای اصلاح جواب در نقطه i در هر زمان، به شکل زیر عمل میکنیم.
رابطه 31
مثال
رابطه جابهجایی یک بعدی را در نظر بگیرید. در صورتی که روش تفاضل محدود روی این معادله اعمال شود. حل عددی آن را میتوان به شکل زیر محاسبه کرد.
رابطه 32
برای گسسته سازی معادله بالا، از روش تقریب تفاضل مرکزی برای مشتق مکانی و اویلر پیش رو برای مشتق زمانی استفاده میشود. این موضوع در رابطه زیر نشان داده شده است.
رابطه 33
توجه کنید که این حالت همان روش زمان پیش رو و مکان مرکزی که در بالا ارائه شد را نشان میدهد و تنها تفاوت آن این است که ترم نفوذ یا دیفیوژن در آن حذف شده است. در ادامه و برای حل عددی این رابطه، شرط اولیه را برابر با رابطه زیر در نظر میگیریم.
رابطه 34
در این مثال، ناحیه حل را برای تولید شبکه در محدود 0 تا ۱ و شرایط مرزی را پریودیک در نظر میگیریم. کد متلب این مثال در ادامه آورده شده است.
% This Matlab script solves the one-dimensional convection
% equation using a finite difference algorithm. The
% discretization uses central differences in space and forward
% Euler in time.
clear all;
close all;
% Number of points
Nx = 50;
x = linspace(0,1,Nx+1);
dx = 1/Nx;
% velocity
u = 1;
% Set final time
tfinal = 10.0;
% Set timestep
dt = 0.001;
% Set initial condition
Uo = 0.75*exp(-((x-0.5)/0.1).ˆ2)';
t = 0;
U = Uo;
% Loop until t > tfinal
while (t < tfinal),
% Forward Euler step
U(2:end) = U(2:end) - dt*u*centraldiff(U(2:end));
U(1) = U(end); % enforce periodicity
% Increment time
t = t + dt;
% Plot current solution
clf
plot(x,Uo,'b*');
hold on;
plot(x,U,'*','color',[0 0.5 0]);
xlabel('x','fontsize',16); ylabel('U','fontsize',16);
title(sprintf('t = %f\n',t));
axis([0, 1, -0.5, 1.5]);
grid on;
drawnow;
end
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
% This Matlab script solves the one-dimensional convection
% equation using a finite difference algorithm. The
% discretization uses central differences in space and forward
% Euler in time.
clear all;
close all;
% Number of points
Nx = 50;
x = linspace(0,1,Nx+1);
dx = 1/Nx;
% velocity
u = 1;
% Set final time
tfinal = 10.0;
% Set timestep
dt = 0.001;
% Set initial condition
Uo = 0.75*exp(-((x-0.5)/0.1).ˆ2)';
t = 0;
U = Uo;
% Loop until t > tfinal
while (t < tfinal),
% Forward Euler step
U(2:end) = U(2:end) - dt*u*centraldiff(U(2:end));
U(1) = U(end); % enforce periodicity
% Increment time
t = t + dt;
% Plot current solution
clf
plot(x,Uo,'b*');
hold on;
plot(x,U,'*','color',[0 0.5 0]);
xlabel('x','fontsize',16); ylabel('U','fontsize',16);
title(sprintf('t = %f\n',t));
axis([0, 1, -0.5, 1.5]);
grid on;
drawnow;
end
نتایج این حل با استفاده از روش تفاضل محدود در شکلهای زیر برای زمان t=0.5 ، t=0.25 و t=1.0 به تصویر کشیده شده است.
در مباحث بعدی وبلاگ فرادرس نشان داده میشود که الگوریتم FTCS برای هر مقداری از Δt برای معادله «جابهجایی خالص» (Pure Convection) ناپایدار است. این موضوع در این مثال نیز به خوبی نشان داده شده است.
شرایط مرزی
در این بخش، پیاده سازی و اجرای روش تفاضل محدود را روی مرزها مورد بررسی قرار میدهیم. توجه کنید که بررسی تمام شرطهای مرزی در یک مطلب دیگر در وبلاگ فرادرس که مرتبط بادینامیک سیالات محاسباتی و روشهای عددی است به صورت دقیق مورد بررسی قرار میگیرد و در مطلب حاضر تنها شرط مرزی دریکله مطالعه میشود.
شرط مرزی دریکله، شرطی است که در آن، حالت مطلوب در مرز مشخص شده است. برای مثال شرط مرزی دریکله در یک مسئله انتقال حرارت، دما را در مرز ناحیه حل مشخص میکند.
به عنوان مثال دیگر برای شرط مرزی دریکله، یک مسئله یک بعدی نفوذ و جابهجایی را در نظر بگیرید. در این مسئله سرعت ورودی به عنوان شرط مرزی دریکله معرفی میشود. بنابراین برای حل تفاضل محدود، مقدار سرعت را در ورودی برابر با مقداری ثابت در نظر میگیریم و معادلات را در نقاط میانی میدان حل بررسی میکنیم.
رابطه 35
رابطه بالا، سرعت را در نقطه و «نود» (Node) یا گره اول یعنی i=0 نشان میدهد. این نود، مرز ورودی را نشان میدهد. توجه کنید که در حلهای عددی نقاط موجود در شبکه و مش حل را نود یا گره مینامند.
در مسئله یک بعدی نفوذ و جابهجایی که در ابتدای این بخش بیان شد، گسسته سازی در نقطه اول یعنی i=1 با روش تفاضل مرکزی انجام میشود. در نهایت این تقریب تفاضل مرکزی، به شکل زیر در میآید.
رابطه 36
در ادامه و با توجه به آنکه مقدار U0 معلوم و ثابت (شرایط مرزی دریکله) است، مقدار آن را در رابطه بالا جایگذاری میکنیم و در نهایت رابطه تفاضل محدود به شکل زیر در میآید.
رابطه 37
در قسمت قبل و در مسیر اجرای یک حل عددی دیدیم که ما مجهولات مسئله که در اینجا سرعت است را به صورت یک بردار نشان میدهیم. بنابراین در حالتی که شرط مرزی دریکله داریم، یکی از اجزای این بردار که مقدار سرعت در مرز را نشان میدهد، معلوم است. بنابراین در این شرایط، بردار مجهولات را به شکل زیر نمایش میدهیم.
رابطه 38
یکی دیگر از موارد بسیار مهم در مسیر اجرای یک حل عددی تعیین پارامتر b در رابطه ۲۴ است. در واقع در مسائلی که شرط مرزی دریکله در آن به کار رفته، پارامتر b مقدار شرط مرزی را نیز در بر میگیرد. برای مثال رابطه بازنویسی شده به وسیله شرایط مرزی در مسئله نفوذ و جابهجایی (رابطه ۳۷) در نظر بگیرید. در این رابطه مقدار b را میتوان به شکل زیر محاسبه کرد.
رابطه ۳۹
سایر شرطهای مرزی مانند شرط مرزی نیومن و کوشی در سایر مطالب وبلاگ فرادرس به صورت دقیق مورد بررسی قرار میگیرد.
بنابراین همانطور که اشاره شد، روش تفاضل محدود، به منظور محاسبه مشتقهای جزئی، محیط و دنیای پیوسته پیرامون ما را به محیط گسسته تبدیل میکند و با تعریف گره و شبکه حل در حالت یک بعدی دو بعدی و سه بعدی، میتوان حل عددی معادلات دیفرانسیل با مشتقات جزئی موجود در ریاضیات و دینامیک سیالات محاسباتی را با روش تفاضل محدود مورد مطالعه قرار داد.
رابطه زیر که در بالا نیز مورد مطالعه قرار گرفت را دوباره در نظر بگیرید.
رابطه 29
این رابطه به صورت نیمه گسسته است. زیرا گسسته سازی تنها در مکان انجام شده و در زمان از تقریب تفاضل محدود استفاده نشده است. برای اینکه معادله نفوذ و جابهجایی ارائه شده را به صورت عددی حل کنیم، لازم است که گسسته سازی آن را به صورت کامل انجام دهیم. برای این منظور، با استفاده از روشهایی که در معادلات دیفرانسیل معمولی بیان کردیم سمت چپ معادله فوق را گسسته میکنیم. برای مثال در صورتی که از روش اویلر پیشرو ساده استفاده کنیم، معادله فوق به شکل زیر در میآید.
رابطه 30
استفاده از روش تفاضل محدود مرکزی برای بیان مشتقات مکانی و روش اویلر پیشرو برای بیان مشتقات زمانی یکی از رایجترین، پرکاربردترین و معروفترین روشهای موجود در تفاضل محدود ودینامیک سیالات محاسباتی[url=https://blog.faradars.org/computational-fluid-dynamics/][/url] را نتیجه میدهد. این روش به روش «زمان پیش رو و مکان مرکزی» (Forward Time Central Space) معروف است. در محاسبات عددی و CFD این روش را به صورت اختصاری با نماد FTCS نمایش میدهند.
توجه کنید که این رابطه یک رابطه صریح است و نیاز به بیان پارامتر A به صورت صریح نیست. نکته دیگر این است که برای اصلاح جواب در نقطه i در هر زمان، به شکل زیر عمل میکنیم.
رابطه 31
مثال
رابطه جابهجایی یک بعدی را در نظر بگیرید. در صورتی که روش تفاضل محدود روی این معادله اعمال شود. حل عددی آن را میتوان به شکل زیر محاسبه کرد.
رابطه 32
برای گسسته سازی معادله بالا، از روش تقریب تفاضل مرکزی برای مشتق مکانی و اویلر پیش رو برای مشتق زمانی استفاده میشود. این موضوع در رابطه زیر نشان داده شده است.
رابطه 33
توجه کنید که این حالت همان روش زمان پیش رو و مکان مرکزی که در بالا ارائه شد را نشان میدهد و تنها تفاوت آن این است که ترم نفوذ یا دیفیوژن در آن حذف شده است. در ادامه و برای حل عددی این رابطه، شرط اولیه را برابر با رابطه زیر در نظر میگیریم.
رابطه 34
در این مثال، ناحیه حل را برای تولید شبکه در محدود 0 تا ۱ و شرایط مرزی را پریودیک در نظر میگیریم. کد متلب این مثال در ادامه آورده شده است.
% This Matlab script solves the one-dimensional convection
% equation using a finite difference algorithm. The
% discretization uses central differences in space and forward
% Euler in time.
clear all;
close all;
% Number of points
Nx = 50;
x = linspace(0,1,Nx+1);
dx = 1/Nx;
% velocity
u = 1;
% Set final time
tfinal = 10.0;
% Set timestep
dt = 0.001;
% Set initial condition
Uo = 0.75*exp(-((x-0.5)/0.1).ˆ2)';
t = 0;
U = Uo;
% Loop until t > tfinal
while (t < tfinal),
% Forward Euler step
U(2:end) = U(2:end) - dt*u*centraldiff(U(2:end));
U(1) = U(end); % enforce periodicity
% Increment time
t = t + dt;
% Plot current solution
clf
plot(x,Uo,'b*');
hold on;
plot(x,U,'*','color',[0 0.5 0]);
xlabel('x','fontsize',16); ylabel('U','fontsize',16);
title(sprintf('t = %f\n',t));
axis([0, 1, -0.5, 1.5]);
grid on;
drawnow;
end
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
% This Matlab script solves the one-dimensional convection
% equation using a finite difference algorithm. The
% discretization uses central differences in space and forward
% Euler in time.
clear all;
close all;
% Number of points
Nx = 50;
x = linspace(0,1,Nx+1);
dx = 1/Nx;
% velocity
u = 1;
% Set final time
tfinal = 10.0;
% Set timestep
dt = 0.001;
% Set initial condition
Uo = 0.75*exp(-((x-0.5)/0.1).ˆ2)';
t = 0;
U = Uo;
% Loop until t > tfinal
while (t < tfinal),
% Forward Euler step
U(2:end) = U(2:end) - dt*u*centraldiff(U(2:end));
U(1) = U(end); % enforce periodicity
% Increment time
t = t + dt;
% Plot current solution
clf
plot(x,Uo,'b*');
hold on;
plot(x,U,'*','color',[0 0.5 0]);
xlabel('x','fontsize',16); ylabel('U','fontsize',16);
title(sprintf('t = %f\n',t));
axis([0, 1, -0.5, 1.5]);
grid on;
drawnow;
end
نتایج این حل با استفاده از روش تفاضل محدود در شکلهای زیر برای زمان t=0.5 ، t=0.25 و t=1.0 به تصویر کشیده شده است.
در مباحث بعدی وبلاگ فرادرس نشان داده میشود که الگوریتم FTCS برای هر مقداری از Δt برای معادله «جابهجایی خالص» (Pure Convection) ناپایدار است. این موضوع در این مثال نیز به خوبی نشان داده شده است.
شرایط مرزی
در این بخش، پیاده سازی و اجرای روش تفاضل محدود را روی مرزها مورد بررسی قرار میدهیم. توجه کنید که بررسی تمام شرطهای مرزی در یک مطلب دیگر در وبلاگ فرادرس که مرتبط بادینامیک سیالات محاسباتی و روشهای عددی است به صورت دقیق مورد بررسی قرار میگیرد و در مطلب حاضر تنها شرط مرزی دریکله مطالعه میشود.
شرط مرزی دریکله، شرطی است که در آن، حالت مطلوب در مرز مشخص شده است. برای مثال شرط مرزی دریکله در یک مسئله انتقال حرارت، دما را در مرز ناحیه حل مشخص میکند.
به عنوان مثال دیگر برای شرط مرزی دریکله، یک مسئله یک بعدی نفوذ و جابهجایی را در نظر بگیرید. در این مسئله سرعت ورودی به عنوان شرط مرزی دریکله معرفی میشود. بنابراین برای حل تفاضل محدود، مقدار سرعت را در ورودی برابر با مقداری ثابت در نظر میگیریم و معادلات را در نقاط میانی میدان حل بررسی میکنیم.
رابطه 35
رابطه بالا، سرعت را در نقطه و «نود» (Node) یا گره اول یعنی i=0 نشان میدهد. این نود، مرز ورودی را نشان میدهد. توجه کنید که در حلهای عددی نقاط موجود در شبکه و مش حل را نود یا گره مینامند.
در مسئله یک بعدی نفوذ و جابهجایی که در ابتدای این بخش بیان شد، گسسته سازی در نقطه اول یعنی i=1 با روش تفاضل مرکزی انجام میشود. در نهایت این تقریب تفاضل مرکزی، به شکل زیر در میآید.
رابطه 36
در ادامه و با توجه به آنکه مقدار U0 معلوم و ثابت (شرایط مرزی دریکله) است، مقدار آن را در رابطه بالا جایگذاری میکنیم و در نهایت رابطه تفاضل محدود به شکل زیر در میآید.
رابطه 37
در قسمت قبل و در مسیر اجرای یک حل عددی دیدیم که ما مجهولات مسئله که در اینجا سرعت است را به صورت یک بردار نشان میدهیم. بنابراین در حالتی که شرط مرزی دریکله داریم، یکی از اجزای این بردار که مقدار سرعت در مرز را نشان میدهد، معلوم است. بنابراین در این شرایط، بردار مجهولات را به شکل زیر نمایش میدهیم.
رابطه 38
یکی دیگر از موارد بسیار مهم در مسیر اجرای یک حل عددی تعیین پارامتر b در رابطه ۲۴ است. در واقع در مسائلی که شرط مرزی دریکله در آن به کار رفته، پارامتر b مقدار شرط مرزی را نیز در بر میگیرد. برای مثال رابطه بازنویسی شده به وسیله شرایط مرزی در مسئله نفوذ و جابهجایی (رابطه ۳۷) در نظر بگیرید. در این رابطه مقدار b را میتوان به شکل زیر محاسبه کرد.
رابطه ۳۹
سایر شرطهای مرزی مانند شرط مرزی نیومن و کوشی در سایر مطالب وبلاگ فرادرس به صورت دقیق مورد بررسی قرار میگیرد.
بنابراین همانطور که اشاره شد، روش تفاضل محدود، به منظور محاسبه مشتقهای جزئی، محیط و دنیای پیوسته پیرامون ما را به محیط گسسته تبدیل میکند و با تعریف گره و شبکه حل در حالت یک بعدی دو بعدی و سه بعدی، میتوان حل عددی معادلات دیفرانسیل با مشتقات جزئی موجود در ریاضیات و دینامیک سیالات محاسباتی را با روش تفاضل محدود مورد مطالعه قرار داد.