книги / Основы математического моделирования и численные методы
..pdfПо теореме Остроградского - Гаусса имеем
Здесь S - площадь боковой поверхности объема V (рис. 3.16).
Рис. 3.16. Исследуемая область
Аналогичным образом в уравнении (3.42) преобразуется инте
грал
— f x — ^dV ду[ ду,
Объединяя соотношения (3.43) и (3.44) с аналогичными соот ношениями для приведенного выше интеграла и учитывая, что dV = hdA и dS = hdL, получим уравнение (3.42) в виде
ЭI |
Л |
Г, эЛ |
\ |
|
|
/X + |
{ |
КУ |
d L - |
дх J |
дУ) |
/ |
||
à[N T (у ЭЛ Э[N]тг |
dt |
|
||
- J дх Х дх) +' |
ду I |
|
ду - W T9r dA = 0. (3.45) |
|
Здесь предполагается, что толщина элемента h равна единице. Интеграл по контуру в формуле (3.45) может быть выражен че
рез величину dt/dn, где п - внешняя нормаль к границе:
^э[лг]г |
ы |
| a [ y f |
дх |
дх |
ду ду |
Неизвестная функция t в уравнении (3.46) определяется соот ношением
так что
* |
э а в |
. » |
и |
|
|
(3.47) |
дх |
дх |
дх |
|
ду |
ду |
ду |
После подстановки полученных формул (3.47) в первый инте грал (3.46) имеем
£ |
'd \ N ^ |
Э[ЛГ] |
3[W]r 3[W ]' <и{т]. |
(3.48) |
|
дх |
Эх |
ду |
ду |
|
|
Поскольку в уравнении (3.46) интеграл J[iV]r A.-^-d£ содержит
дп
dt
выражение X— , которое определяет плотность теплового потока с
дп
поверхности рассматриваемой области (рис. 3.17) [5], с помощью этого интеграла учитываем граничные условия второго рода:
. dt
и третьего рода:
. dt ( |
ч |
|
~*о) > |
где q - заданная плотность теплового потока; ts - искомая темпе ратура границы тела; t0 - заданная температура окружающей среды.
|
Рис. 3.17. Поток тепла на границе области |
||
Тогда |
|
|
|
Э/ |
|
|
j[N ]Tat0dL.(3A9) |
\[ N ^ X - d L = - \[ N ^ q d L - a \[ N ^ [ N ] { T } d L + |
|||
дп |
L„ |
La |
La |
|
|||
Уравнение (3.46) с учетом выражений (3.48) и (3.49) запишется
как
Э[ЛГ]Г Э М |
Э[ЛГ]г Э[ЛГ] |
dA{T}+ |а [ # ] Г[У]<й{Г}- |
|
||
дх |
дх |
ду |
ду |
|
|
|
|
||||
l[N ]Tq d L -l[ N ] Tat0d L - \ [ N f qvdA = 0. |
(3.50) |
||||
Lq |
|
LQL |
|
A |
|
Перепишем уравнение (3.50) в виде |
|
||||
|
|À [5]r [5]âL4+ |
m = |
|
||
|
A |
|
La |
|
|
r |
J[AT]rqdL + |
|
\ |
|
|
- |
\[N ]T at0dL + $[N]T qvdA |
(3.51) |
|||
L1 |
La |
A |
Выражение (3.51) в матричной форме принято записывать сле дующим образом:
П>Узки), {^} ={/■}„ +{F }, +{*■}»
Здесь |
= |
A |
La |
= - 1М Г 4 d L ; {F}, = { К |
; {F}w = J[W]7q v d A . |
|
А |
Эффективная процедура решения задач МКЭ строится на том, что сначала определяются матрица коэффициентов и вектор-столбец свободных членов для конкретного текущего конечного элемента. Все узловые неизвестные, которые не относятся к этому элементу, исключаются из рассмотрения.
Переходим к матрице по элементам:
!(([*£’+ [* £ > } = { /£ ’+{/}?+{/£*)•
Здесь |
[*£>= Д х » [ л Г [ * Г ) * : |
[* £ * - |
— / « [ л Г -< //;{/}“ =
Будем считать, что X и qv постоянны в пределах элемента. Интерполяционный полином для одномерного линейного эле
мента имеет вид
/ = а д + а д + е д = [ л г ] { т],
где
|
|
|
|
a,= X jY k - i :trt; |
|
1 |
|
|
X |
y |
||
2А |
|
|
' b, = Y j-Y t ; |
или N ,= ± - 1 |
|
Xj |
Yj |
|||||
|
|
ГЬ II |
•h |
J* |
|
1 |
|
xk Yk |
||||
|
|
|
|
|
1 |
|
|
|
|
|
|
|
|
|
|
|
aj = X kYi - X iYk; |
1 |
1 |
|
x , |
Yi |
|||
|
( aj +bJx +cjy)> \b . =Yk - Y ; |
|
|
|
|
|
||||||
Nj =ü |
или N /= — |
1 |
|
X |
У |
|||||||
|
|
|
|
J |
k |
' |
1 |
2A |
1 |
|
|
|
|
|
|
|
и j* l & |
|
|
|
* k |
Yk |
|||
1 |
|
|
|
ak =X |
|
|
|
1 |
|
X, |
Y, |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
Nk = 2А^к +ЬкХ +СкУ); ■bk = Yf - Y j ; |
|
|
1 |
|
X J |
Yj |
||||||
|
|
|
|
t----- |
II |
1 |
|
|
1 |
|
X |
У |
Тогда для треугольного симплекс-элемента матрица градиентов |
||||||||||||
запишется как |
|
|
|
|
|
|
|
|
|
|
||
|
' B ' |
|
|
|
’ЭNt |
BNj |
алг, 1 |
bt |
h |
bk |
||
[5] = |
Эх |
i |
|
II |
dx |
Эх |
dx |
|||||
|
—i |
|
|
|
|
|
|
|
|
|
||
‘ |
Э |
V |
|
|
BNt |
BNj |
dNk 2 А Ci |
cj |
ck |
|||
|
|
|
||||||||||
|
s>y. |
|
|
|
|
|
By |
|
|
|
|
|
|
|
|
_1_ |
Yj~Yk |
Yk -Yt |
Y,-Y, |
|
|
|
|
|
|
|
|
|
2A X k - X j X - X k X j - X t |
|
|
|
|
|
||||
Определим члены уравнения для текущего элемента: |
|
|
|
|||||||||
|
[ i f = f ( 4 * f И * |
= Ч * Г [5] \ds = Х [яГ [В] А = |
|
|||||||||
|
|
х_ |
bfo+cfi, |
bibj+c.Cj |
ЬА +c,ck |
|
|
|
|
|
||
|
|
bA |
+ CjC, |
bjbj+CjCj |
bjbk + cjCk |
|
|
|
|
|
||
|
|
4А |
|
|
|
|
|
|||||
|
|
bkb, +ckCj |
bkbj + ckCj |
bkbk +ckck |
|
|
|
|
|
|||
|
|
|
|
|
|
|
|
|||||
Следующий интеграл определим с учетом выражения для опре деления интеграла для локальных координат треугольного сим
плекс-элемента (3.15) JLiaL2bL3cdA = - а 'Ь-с - |
2А: |
(я + й + с + 2 )! |
|
'N,1
{ / } Ц{ = \ ч у [ ^ dx = qv [<Nj
k
и |
s |
& |
J
A '
&• n
A[A J
* tu.
Y
l > ч.l -
Матрица коэффициентов, учитывающая граничные условия третьего рода для пограничного элемента на рис. 3.17, определяется следующим образом:
" 0 '
Ja[JV]7' [#]«// = a J M 7 |
[N]dl = a j |
Nj |
||
Lo. |
|
Ax |
La Л . |
|
'o |
0 |
0 |
'0 |
0 |
= a j 0 |
NJNJ |
N jN k |
S»ilP«»---- 0 |
AA |
1-0 0 |
NkNj |
NkNk_ |
Ax 0 |
A A |
1 |
v** |
----I |
|
O |
|
0
A A dl.
A A .
Используя |
соотношение |
(3.13) |
? |
a |
L> |
b |
Ct\b\ |
||
IЦ |
|
dx = ------ |
-.-Z, для |
||||||
|
|
|
|
; |
|
|
|
(а+ft + l)! |
|
стороны элемента с узлами j |
и к |
получим |
|
|
|
|
|||
|
|
|
0 |
0 |
0 |
|
|
|
|
|
[ * f |
f. |
0 |
2 |
1 |
|
|
|
|
|
1 Ja |
|
|
|
|
|
|
|
|
|
|
|
0 |
1 |
2 |
|
|
|
|
Здесь Ljk - |
длина стороны треугольника с узлами j |
и к. |
|||||||
Аналогично для двух других сторон |
|
|
|
|
|
||||
|
|
ai,, |
" 2 |
1 |
0 ' |
|
|
|
|
|
|
1 |
2 |
0 |
|
|
|
|
|
|
[*]<*>=:=£ |
|
|
|
|
||||
|
L Ja |
6 |
|
|
|
|
|
|
|
|
|
|
0 |
0 |
0 |
|
|
|
|
|
|
|
м а(«) |
а h i |
2 |
О |
1 |
|
|
|
|
0 |
0 0 |
|
|
||
|
|
|
|
6 |
1 |
|
О 2 |
|
|
|
|
|
|
|
|
||
|
Два |
оставшихся |
интеграла |
в |
векторе |
свободных членов |
||
| |
= |
J ett0[ N fe^Tdl |
и {/}^e) = |
- |
^r[J7V](e)rdl |
вычисляются как |
||
|
4e) |
|
? |
4e) |
|
|
|
|
|
|
|
{/}a(«) |
2 |
1 |
, |
|
|
|
|
|
|
|
0 |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
0 |
|
|
|
|
|
|
atoLjk |
|
|
||
|
|
|
М Г |
2 |
- < 1 |
, |
|
|
|
|
|
|
|
|
1 |
|
|
|
|
|
|
|
|
1 ' |
|
|
|
|
|
{ / |
t |
2 |
-<oL |
|
|
|
|
|
|
|
1 |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
^ - |
|
|
|
|
|
|
|
|
г |
|
|
|
|
|
{ /£ ’ |
g1» |
|
|
||
|
|
|
2 |
< 1 |
’ |
|
||
|
|
|
|
|
|
0 |
|
|
(e) __
{ /Г =
3.13. Решение нестационарной одномерной задачи
теплопроводности МКЭ
Рассмотрим задачу нагревания изолированного провода (рис. 3.18), которая описывается нестационарным уравнением теп лопроводности в одномерной осесимметричной постановке
Э/* |
1 Э |
,d t |
Pc— = |
г or |
Л Э7J + у » |
от |
где t* - температура (искомая функция); X - коэффициент тепло проводности; р - плотность; с - удельная теплоемкость; qv - мощ ность внутреннего источника тепла.
Рис. 3.18. Схема изолированного провода
Граничные условия:
ал |
= 0 |
; |
|
|
- н а оси — |г=0 |
|
|
||
|
|
ЭЛ |
ос/. I |
\ |
-навнеш ней поверхности — |г=п>= |
^ |
~*о)• |
||
В нулевой момент времени температура постоянна: f*(r,0) = tH.
Внутренний источник тепла qv Ф0 . Применяя метод Галеркина, получим
И |
т( dt 1 Э ( |
, d t ' + qy dV = О, |
(3.52) |
I рсэ Г 7 ¥> . |
дг) |
|
где t - приближенное решение. Поскольку
1 |
д |
dt W |
1 д ( . |
a»ï, э[лг]' |
|
|
|
|
|
|
- О Т |
|
Эг |
\& |
|
получим |
|
|
|
|
|
|
|
М |
г 1 Э |
Л Э/') |
1 Э ( т^т\г ( ^ |
dt |
э [ < |
f |
& |
" а |
rX dr |
|
|
dr |
v |
(3.53) |
|
|
г дг |
|
|
dr |
|||
С учетом формулы (3.53) можно записать
|
or у |
|
dV |
По теореме Остроградского - Гаусса имеем |
|
( . |
d t\\ |
} r o r \ |
dr |
Здесь j [ N f |
- интеграл по замкнутой поверхности |
тела.
Тогда для одномерной осесимметричной задачи интеграл через боковую поверхность определится как
+ 2uLzrb[JVf
a объемные интегралы
|
N |
|
rdadzdr = |
0 0 0 |
J |
f 'à[N\T Э(Л rdr;
дг Эr
2TCAr
j[AT] ^ = j j J[w f qyrdadzdr = 2nLz J[w f qvrdr\
V |
|
0 |
0 |
0 |
0 |
|
J [A |
pfc ^ - d V = |
J |
J J|V]7[p c ^ -rd a d zd r = 2nLz J[7V]rp c ^ -rd r . |
|||
V |
* |
0 0 0 |
^ |
0 |
^ |
|
Неизвестная функция / в уравнении (3.52) определяется соот
ношением t = [N]{T] .
Так что
эг э( № } ) _ э[аг]
|
|
У-1 J |
дг |
Э([ЛГ]{Г}) |
а{г}. |
Эт |
Эт |
Эт |
dt |
|
|
дг Г=П> - ? |
( ' U - * ) = ~ |
( [ * ] w U |
С учетом полученных выражений уравнение (3.52) можно запи сать как