Математичне моделювання процесів теплопровідності в елементах конструкції

Размещено на http://www.allbest.ru/

МІНІСТЕРСТВО ОСВІТИ І НАУКИ УКРАЇНИ

ДВНЗ «Український державний хіміко-технологічний університет»

Кафедра енергетики

КУРСОВА РОБОТА

з дисципліни: «Математичні методи та моделі енергетичного обладнання у розрахунках на ЕОМ»

на тему

«Математичне моделювання процесів теплопровідності в елементах конструкції»

Виконала: ст. гр. 3-ЕТТ-11

Іщук В.І.

Перевірила

доцент Швачич С. В.

Дніпропетровськ 2015



Завдання на курсову роботу

Скласти ММ та реалізувати на ЕОМ находження температурного поля в середині пластини. На зовнішніх поверхнях підтримується постійний тепловий потік, у початковий момент часу T0=const. Задачу вирішити численним методом, заздалегідь задав початкові данні. Для ілюстрації результатів побудувати необхідні графічні залежності.

Початкові дані для пластин приведені у таблиці 1.

Табл.1

№	Матеріал	?, кг/м3	С, кДж/(кг*К)	?, Вт/(м*К)	? м.
1	Сталь 40	7850	0,466	41	0,5

Тепловий потік с однієї сторони q1=3000 Вт/м?, із іншої q2=5000 Вт/м?

T

q1 q2

x=0 x=? X

?

Рис.1 Схема пластини



Реферат

У курсовій роботі розроблена математична модель визначення розподілу температурного поля всередині пластини в залежності від часу.

На поверхнях пластини підтримується постійний тепловий потік.

Розроблена математична модель була реалізована на ЕОМ при використанні чисельного методу кінцевих різниць.



Вступ

Математичне моделювання є одним із основних сучасних методів дослідження. Загалом під моделюванням розуміється процес дослідження реальної системи, який включає побудову моделі, її дослідження та пере- несення одержаних результатів на досліджувану систему.[2] Математичне моделювання з кожним роком знаходить все більш широке застосування в інженерній практиці: при проведенні промислових експериментів, проектуванні і конструюванні технічних систем, управлінні виробництвом і його плануванні. Цим пояснюється особлива актуальність вивчення студентами методів математичного моделювання та набуття навичок їх застосування в теплоенергетиці.[4]



1. Основні положення методу математичного моделювання щодо процесів теплопровідності

теплопровідність пластина час геометричний

Роль моделювання як методу наукового пізнання та методу вирішення технічних завдань завжди оцінювалася досить високо. З розвитком техніки знайшло широке застосування фізичне моделювання споруд, машин і механізмів.

Досягнення математики призвели до поширення математичного моделювання різних об'єктів і процесів.[3]

Слово «модель» походить від латинського modus (копія, образ, обрис). Моделювання -- це заміщення деякого об'єкта А іншим об'єктом Б. Заміщуваний об'єкт А називається оригіналом чи об'єктом моделювання, а що заміщає Б -- моделлю. Іншими словами, модель -- це об'єкт-замінник об'єкта-оригіналу, що забезпечує вивчення деяких властивостей оригіналу.

Метою моделювання є одержання, обробка, представлення та використання інформації про об'єкти, які взаємодіють між собою і зовнішнім середовищем, а модель тут виступає як засіб пізнання властивостей і закономірності поведінки об'єкта.

Моделювання широко використовується в різних сферах людської діяльності, особливо в сферах проектування та управління, де особливими є процеси прийняття ефективних рішень на основі одержуваної інформації.

Зупинимося на одному з найбільш універсальних видів моделювання математичному, що ставить Математичне моделювання -- це засіб вивчення реального об'єкта, процесу або системи шляхом їх заміни математичної моделлю, більш зручною для експериментального дослідження за допомогою ЕОМ.

Математична модель є наближеним поданням реальних об'єктів, процесів або систем, вираженим в математичних термінах, та які зберігають істотні риси оригіналу. Математичні моделі в кількісній формі, за допомогою логіко-математичних конструкцій, описують основні властивості об'єкта, процесу або системи, його параметри, внутрішні та зовнішні зв'язки.[2]

Диференціальне рівняння теплопровідності є математичною моделлю цілого класу явищ теплопровідності . При інтегруванні диференціального рівняння в приватних похідних отримуємо безліч різних рішень.

Залежність між часовою і просторовою зміною температури в тілі встановлює диференційне рівняння теплопровідності.

Щоб із безмежної кількості явищ виділити одне і дати його повний математичний опис, до диференціального рівняння теплопровідності необхідно добавити умови однозначності, які містять у собі геометричні, фізичні, часові і граничні умови.

Геометричні умови визначають форму і розміри тіла, в якому відбувається досліджуваний процес.

Фізичні умови задаються теплофізичними параметрами тіла ?, cV і розподілом внутрішніх джерел теплоти.

Часові (початкові) умови містять розподіл температури в тілі в початковий момент часу.

Граничні умови визначають особливості проходження процесу на поверхні тіла. Граничні умови можуть задаватися декількома способами.

Граничні умови І роду. У цьому випадку задається розподіл температури на поверхні тіла для кожного моменту часу

Граничні умови ІІ роду. У цьому випадку задається величина теплового потоку для кожної точки поверхні тіла в довільний момент часу

Граничні умови ІІІ роду. У цьому випадку задаються температури середовища tо і умови теплообміну цього середовища з поверхнею тіла.

Процеси теплообміну між середовищем і тілом є винятково складними і залежать від багатьох факторів. Докладніше вони розглянуті у другій часині книги. Для опису інтенсивності теплообміну між поверхнею тіла і середовищем використовується гіпотеза Ньютона-Ріхмана, згідно з якою

qст = ?(tст - tо)

де ? - коефіцієнт пропорціональності, названий коефіцієнтом тепловіддачі, Вт/(м2·К).

У разі, коли коефіцієнт тепловіддачі має великі значення (наприклад, під час кипіння рідини на поверхні тіла), граничні умови ІІІ роду переходять в граничні умови І роду, тому що у цьому випадку температура поверхні тіла стає практично рівною температурі рідини.

Граничні умови ІV роду формуються на основі рівності теплових потоків, які проходять крізь прилеглі одна до іншої поверхні тіл[4]

2. Розроблення математичної моделі та розрахунку

1) Для вирішення даної задачі,використовуємо основне рівняння теплопровідності (без внутрішнього стоку теплоти qv=0).

Так як ми маємо одномірну ,нестаціонарну задачу то рівняння матиме вигляд:

;



2) Початкова умова:

3) Граничні умови:

Оскільки на зовнішніх поверхнях пластини підтримуються постійні теплові потоки(із завдання курсової роботи) q1 та q2 у нас буде два граничних рівняння другого роду:

3. Побудова алгоритму вирішення задач теплопровідності

1. Задаючи число шарів розбиття по координаті N розраховують крок по координаті.

h:=del/(N-1)

2. Враховують умову збіжності.

3. У початковий момент часу ( n = 0) всіх вузлів сітки привласнюють значенняня температури (i=0,1,..., N-1, N) виходячи із заданих початкових умов.

4. Виконують крок по часу і визначають розрахунковий номер моменту часу і поточний час процессу n = n+1;

5. Розраховують температуру на внутрішній і зовнішніх , границях по розрахунковій області за формулами:



6. Проводять аналіз закінчення ітераційного процесу та розрахунок в залежності від заданих умов завдання.

Для рішення одномірної по простору нестаціонарної задачі, скористаємось методом кінцевих різностей. Побудуємо сітку: по координаті візьмемо крок dx = 0,05 м.

Для обчислення кроку за часом скористаємося умовою збіжності:

і виберемо dt = 0,11c

Умова виконується:

Для розв'язку рівнянь заміняємо їх кінцево-різницевим співвідношенням:

де i=1,2..N, j=1,2..M.

Початкова умова:

Граничні умови:

1)

2)

4. Перевірка адекватності отриманих результатів

Перевірку на адекватність робимо за допомогою тестового прикладу з підручника.[3]

Умова задачі:Знайти температурне поле в середині пластини.

На поверхнях підтримується постійний тепловий потік(граничні умови 2 роду).

Матеріал пластини сталь 40Х.Початкова температура to=100 Температури стінок Tc1=950 ; Tc2=800.

Теплові потоки q1=3000 Вт/м?, ; q2=6000 Вт/м?, ; Товщина ? =0.33м. Фізичні властивості матеріалу ?=7799, кг/м3 С=494, кДж/(кг*К) ?=48, Вт/(м*К)

Отримані результати:

Memo1

0,112

t = 90

32,38 17,58 11,87 10,33 10,04 10,01 10,07 10,37 11,52 14,48

t = 180

42,57 25,92 16,63 12,33 10,72 10,31 10,50 11,33 13,18 16,51

t = 270

50,23 32,81 21,64 15,34 12,29 11,19 11,30 12,39 14,58 18,05

t = 360

56,66 38,78 26,44 18,73 14,48 12,63 12,44 13,54 15,85 19,38

t = 450

62,33 44,14 31,00 22,26 17,03 14,46 13,86 14,81 17,10 20,64

t = 540

67,47 49,06 35,34 25,81 19,78 16,55 15,52 16,22 18,40 21,92

t = 630

72,22 53,65 39,49 29,33 22,64 18,83 17,35 17,77 19,80 23,27

t = 720

76,68 57,99 43,47 32,81 25,56 21,23 19,33 19,45 21,30 24,72

t = 810

80,91 62,11 47,30 36,23 28,51 23,71 21,42 21,25 22,93 26,29

t = 900

84,94 66,06 51,02 39,59 31,46 26,26 23,61 23,16 24,67 27,97



Порівнюємо отримані результати в одних і тих же вузлах сітки (1 вузол) з курсової роботи та тестового прикладу та визначаємо похибку.

Інтервал виведення результату	Результат розрахунку	Результат тестового прикладу	Відносна похибка
90	30.5	32.38	1.88
180	43.12	42.57	0.55
270	52.37	50.23	2.14
360	59.91	56.66	3.25
450	66.41	62.33	4.08
540	72.23	67.47	4.76
630	77.54	72.22	5.32
720	82.46	76.68	5.78
810	87.07	80.91	6.16
900	91.42	84.94	6.08

Графік залежності

Результати роботи



Результати тестового прикладу.

З графіків ми можемо побачити як розподіляється температура всередині пластини ,а саме перепад між центром та краями



Висновки

У даній курсовій роботі була застосована явна кінцево-різницева схема для знаходження розподілу температур по товщині пластини. При використанні явних кінцево-різницевих схем величина припустимого кроку за часом обмежена і для внутрішніх вузлів залежить від обраного кроку по координаті і температуропровідності матеріалу . З умови збіжності бачимо, як залежить крок по часу від кроку по координаті. Якщо умова збіжності не буде виконуватися ми отримаємо помилкові результати. При проведенні розрахунків, треба насамперед подбати про те, щоб величина кроку задовольняла умовам стійкості системи кінцево-різницевих рівнянь.

Результати обчислень приведені у вигляді таблиці та графіку залежності

Для перевірки моделі на адекватність ми використовували тестовий приклад.. Звідки видно що похибка температур не перевищує 6.16%.



Список літератури

1. Бєляєв Н. М. Основи теплопередачі. К.: Вища школа, 1989.

2. Дульнев Г. Н. Застосування ЕОМ для рішення задач теплообміну. - М.: Вища школа, 1990.

3. Кафаров В. В., Глєбов М. Б. Математичне моделювання основних процесів хімічних виробництв. - М.: Вища школа,1991.

4. Юдаев Б. Н. Теплопередача. - М.: Вища школа, 1981.



Додаток А

Список ідентифікаторів

t[i,j] - температура в вузлах сітки,;

N-число вузлів ;

Del- товщина,м.

Tau-час процесу,с;

Taup-інтервал виведення результату;

T0-початкова температура ;

Tc1 , Тс2 -температури стінок;

Lam-коефіцієнт теплопровідності Вт/мК;

Ro-щільність ,кг/м? ;

TE-теплоємність кДж/кгК;

q1,q2-теплові потоки Вт/м?;

h-крок по координаті;

dt-крок за часом;

t-лічильник часу процеса;

tp-лічильник часу вивода результату;

Bi-критерій Біо;

Fo-критерій Фурье;



Додаток Б

Програмний код

implementation

{$R *.dfm}

const N=10; { число вузлів сітки по х }

Tau=900; { час процесу}

taup=90; { інтервал виведення результату}

del=0.5;

T0=10;

Tc1=950;

Tc2=800;

al1=0;

al2=0;

lam=41;

Ro=7850 ;

TE=466;

q1=3000;

q2=5000;

Procedure Tform1.Button1Click(Sender: TObject);

var

i:integer;

h,dt,x,t,tp,a,Fo, Bi1, Bi2:real;

t1,t2: array [1..N] of real;

ss: string;

begin

h:=del/(N-1); { крок по координаті х}

dt:=10; { крок за часом}

t:=0; { лічильник часу процеса }

tp:=0; { лічильник часу виводів результатів}

a:=lam/(Ro*TE);

Fo:=(a*dt)/(h*h); {формула Фур'є}

Bi1:=al1*h/lam;

Bi2:=al2*h/lam;

ss:= ss+ FloatToStrF(Fo, ffFixed,4,3)+' ';

Memo1.Lines.Add(ss);

{Задаються початкові значення функції в вузлах сітки}

for i:=1 to N do

t1[i]:= T0;

{Розрахунковий цикл}

while t<Tau do

begin

{ Розрахунок значень функції в вузлах сітки на наступному кроці за часом }

for i:= 2 to N-1 do

begin

t2[i]:=t1[i]+Fo*(t1[i-1]-2*t1[i]+t1[i+1]);

end;

t2[1]:=2*Fo*(bi1*Tc1+q1*h/lam+t1[2])+(1-2*Fo*(1+bi1))*t1[1];

t2[N]:=2*Fo*(Bi2*Tc2+q2*h/lam+t1[N-1])+(1-2*Fo*(1+Bi2))*t1[N];

{Лічильник часу процесу та виводу результатів}

t:=t+dt;

tp:=tp+dt;

{ Переприсвоєння значень функції для нового кроку за часом}

for i:=1 to N do

t1[i]:=t2[i];

{Вивід поточних результатів розрахунку}

if tp>=taup then

begin

ss:= 't = ' + FloatToStr(t);

Memo1.Lines.Add(ss);

ss:=' ';

for i:=1 to N do

ss:= ss+ FloatToStrF(t1[i], ffFixed,5,2)+' ';

Memo1.Lines.Add(ss);

tp:=0;

end;

end;

end;

end.



Додаток В

Результати

Memo1

0,036

t = 90

30,50 12,79 10,24 10,01 10,00 10,00 10,00 10,04 10,46 13,42

t = 180

43,12 18,31 11,49 10,20 10,02 10,01 10,03 10,25 11,38 15,52

t = 270

52,37 24,17 13,63 10,73 10,12 10,04 10,12 10,61 12,36 17,06

t = 360

59,91 29,80 16,32 11,64 10,36 10,12 10,28 11,05 13,30 18,32

t = 450

66,41 35,08 19,29 12,89 10,79 10,30 10,51 11,55 14,18 19,40

t = 540

72,23 40,02 22,40 14,40 11,40 10,58 10,81 12,08 15,01 20,37

t = 630

77,54 44,66 25,55 16,11 12,19 10,98 11,18 12,63 15,79 21,26

t = 720

82,46 49,06 28,71 17,98 13,14 11,49 11,62 13,20 16,53 22,09

t = 810

87,07 53,23 31,84 19,95 14,24 12,12 12,12 13,78 17,24 22,86

t = 900

91,42 57,22 34,92 22,00 15,47 12,86 12,69 14,38 17,94 23,61

Размещено на Allbest.ur

...