Курс лекций - Математическое моделирование технических объектов (1075784), страница 11

Просмтор этого файла доступен только зарегистрированным пользователям. Но у нас супер быстрая регистрация: достаточно только электронной почты!

Текст из файла (страница 11)

Такое выделение позволяет записать вклад отдельного конечного элемента в общую сумму (12.1) в более общем матричном виде:Введем обозначения:KV[e] =  [B[e] T] [B[e]]dV(12.8)V[e]KS[e] = h [N[e] T] [N[e]]dS(12.9)[e] T(12.10)S2 [e]FS1[e] = q [N] dSS1[e]FS2[e] = hTOC[N[e] T] dS(12.11)S2 [e]и запишем в матричной форме соотношение, представляющее вклад отдельного конечногоэлемента в общую сумму (12.1):д[e]= ( [KV[e] ] + [KS[e] ] ) {T} + {FS1(e)} + {FS2(e)} = [K[e] ] {T} + {F(e)}(12.12)дTЗдесь матрица теплопроводности конечного элемента[ K(e)] и его вектор нагрузки { F(e)}называются далее матрицами е-го конечного элемента.

Приравнивая данное выражение нулю, запишем окончательную систему уравнений в краткой форме:[ K]  { T } = { F}где:Е(12.13)49[K] =[K[e]](12.14)е=1и:Е[F] = -{F[e]}(12.15)е=1Рассмотрим методику получения матриц конечных элементов на нескольких примерах решения задачи переноса тепла в стержне.Пример 12.1. Одномерный случай переноса теплаТребуется вычислить температуру TХ на правом конце стержня, если температура еголевого конца поддерживается равной величине T1=150оС.

Радиус стержня R=1(см), длинаL=7,5 (см). Коэффициенты теплопроводности материала стержня и конвективного теплообмена по всей поверхности стержня соответственно равны: =75 Вт/(см  ОС), h = 10 Вт/(см2 ОС). Температура окружающей среды равна ТОС=40 ОС.Рис. 12.1Рис.12.2Решение.Разобьем стержень на три одинаковых одномерных линейных элементов длиной 2.5см., как это показано на рисунке (12.1), и запишем интерполяционный полином для первогоэлемента:T = N1T1+N2T2Совместим начало координат с 1-м узлом, тогда ФФ примут вид:N1 = (1–x/L) и N2= (x/L)Запишем выражения для матриц [N(e)] и [B(e)]:[N(e)] = [(1–x/L) (x/L)][B(e)] = [ (-1/L) (1/L)]Согласно (12.8) матрица KV(1) примет вид:LKV[1] =0-1L1LL[-11LL] Adx =АL21-1-11dx0После вычисления интеграла окончательно имеем:KV[1] =АL1-1-11(12.16)50Для определения матрицы KS(1) рассмотрим все поверхности конечного элемента 1, обозначенные на рисунке 12.2 через S1, S2 и S3. Через эти поверхности конечный элемент теряеттепло за счет конвекции (h).

Через поверхность S1 конвективного обмена с окружающей средой нет, так как здесь по всей поверхности поддерживается постоянная температура 150 оС.Через поверхность S3 конвективный обмен у первого элемента также отсутствует. То естьдолжна учитываться только поверхность S3. Диаметр стержня не изменяется по оси ОХ, поэтому дифференциал dS в (12.8) примет вид: dS = (Pdx), где Р – периметр, и:LKS[1] = hP1-0xL1L[(1-xL)1L] dx =hPL62112(12.17)Складывая эти матрицы согласно (12.10), получим матрицу [K(1)] теплопроводностидля первого конечного элемента:K[1] =АL1-1-11hPL6+2112(12.18)Матрица теплопроводности второго элемента идентична (12.18). Матрица [K(3)] отличается от (12.18) дополнительным членом, описывающим конвективный обмен со средой поповерхности S3.

Вычислим этот дополнительный поверхностный интеграл, используя выражение (12.9):(12.19)При вычислении интеграла (12.19) учитывалось, что на всей поверхности S3 Ni=0,Nj=1, поскольку эта поверхность является j-м узлом в 3-ем конечном элементе.Рассмотрим теперь интегралы вектора нагрузки.

Начнем с первого конечного элемента. Составляющие вектора нагрузки описывают действие внешних тепловых источников истоков тепловой энергии. Поскольку в нашем примере вообще нет никаких источников тепла,то составляющая q в выражении (12.10) равна нулю и составляющая вектора нагрузки первого элемента FS1 описывается только выражением (12.11) и зависит от величины поверхностиS2:L{FS2(3)} = hTOC PL(1-x/L)x/Ldx =hTOCA211(12.20)0Вектор нагрузки для второго элемента идентичен (12.20).

В векторе же нагрузкитретьего конечного элемента интеграл в (12.20) должен быть вычислен по сумме поверхностей (S2+ S3), через которые происходит отвод тепла у этого элемента. Поскольку площадьS3= А, имеем:{FS2(3)} =hTOC PL211+ hTOCA01(12.21)51Пользуясь выражениями (12.18) и (12.19) построим глобальную матрицу теплопроводности стержня, а с помощью выражений (12.20) и (12.21) – глобальный вектор нагрузки всегостержня. Предварительно вычислим значения термов в этих выражениях: А=(см2),L=2,5(см), P=2(см):A/L= (см2)75[Вт/(смОС)]/2,5(см)=30(Вт/ОС);hPL/6= 10[Вт/(см2 ОС)]2(см)2,5(см)/6=8,3(Вт/ОС);hA= 10[Вт/(см2 ОС)] (см2)=10(Вт/ОС);hTocPL/6 = 10[Вт/(см2 ОС)] 40ОС2(см)2,5(см)/2=1000(Вт);Подставляя полученные значения в (12.18 – 12.21), последовательно находим:[KV(1)]=30[1-1-11] = [KV(2)]=[KS(1)]=8,3[2112] = [KS(2)](Вт/ОС)[KS(3)]=8,3[2112]+[F(1)]=1000[11] = [F(2)](Вт)[F(3)]=1000[11]+[10400[0001[KV(3)]01](Вт/ОС)(Вт/ОС)] (Вт)Объединяя матрицы по методу прямой жесткости, составляем систему (12.13):46,6-21,700T1-21,793,2-21,700-21,793,2-21,7T300-21,756,6T4T21000=200020001400Здесь проведено сокращение на множитель , так как он входит в обе части системыуравнений.

Значение Т1 известно (150оС), поэтому полученная система должна быть модифицирована перед решением. Подробно эта процедура изложена в разделе 14. После модификации система примет вид:46,6000093,2-21,700-21,793,2-21,7T300-21,756,6T41506990T25255=20001400После решения системы имеем:[T]T=[15067,3547,142,8]Результаты расчетов приведены в графическом виде на рисунке 12.3.

На том же рисунке цифрами в скобках отмечены теоретические значения температур, замеренные через каждые 1,5 см. Видно, что полученные в результате расчетов значения достаточно хорошо согласуются с истинными значениями на участке 2,5 см – 7,5 см, то есть ближе к правому концу52стержня. Решение по методу МКЭ можно было бы улучшить, если использовать более короткие элементы вблизи левого конца стержня.В рассмотренном примере площадь поперечного сечения стержня была постоянной. Рассмотрим элемент на рисунке 12.4. Площадь элемента А(х) и его периметр Р(х) меняются линейно по длинеот Аi и Рi на левом конце до Аj и Рj – на правом конце соответственно.

Рассмотрим методику вычисления температурного поля внутри этого элемента1. Записываем выражения для площади боковой поверхности А(х) и для площади периметра Р(х) стержня как функции его длины:Рис. 12.3Рис. 12.4A(х) = Ni Аi + Nj Аj(12.22)Р(х) = Ni Рi + Nj Рj(12.23)где Ni и Nj – линейные ФФ, определенные в примере 12.1.2. Составляем матрицу теплопроводности элемента, для чего заменяя в (12.8) dV наA(x) dx и учитывая (12.16) получим выражение для объемной составляющей матрицы теплопроводности KV(е):L2KV(е) =] 1 -1-1 1[A(x) dxL(здесь величина А(х) не постоянна вдоль оси ОХ, поэтому ее нельзя вынести за знак интеграла) .

Подставляя (12.22) в полученное выражение имеем:LKV(е) =[2L1-1-11]  [(1 -xLx AjL)Ai +] dx0выполняя интегрирование, получаем:KV(е) =2L[1-1-11=L]или:[(LAi +21 -1-11]L Aj2)=Ai + Aj2наконец, обозначая среднюю площадь элемента как Â =(Ai+Aj)/2, имеем:(е)KV =ÂL(е)[1-1-11]Формулы (12.23) = (12.6) с точностью замены площади ее средним значением.(12.23)53Матрица KS(е) с учетом (12.9) примет вид:LKS(е) = hTh[N] [N] dS =SNi2Ni NjNi NjNj2P(x)dx0Вычислим первый коэффициент, определяемый выражением:LLNi2 P(x)dx =0(Ni3 Pi + Ni2 NjPj) dx =L12(3 Pi +Pj)0Вычисляя остальные коэффициенты, получим окончательное выражение для поверхностнойсоставляющей матрицы теплопроводности элемента:KS(е) =hL(3 Pi +Pj)(Pi +Pj)12(Pi +Pj)(3 Pi +Pj)(12.24)Сумма (12.23) и (12.24) и определит выражение для искомой матрицы теплопроводности рассматриваемого конечного элемента.3. Составляем матрицу вектора сил элемента.

Согласно формуле (12.11), матрица вектора сил примет вид:LFS(е) = h TOC[N]T dS = h TOCSNi (N P + N N P )dxi ii j jNj0Откуда, перемножая матрицу-столбец на коэффициент, имеем:LFS(е) = h TOCNi(NiPi + NiNjPj) dxNj(NiPi + NiNjPj)0Вычислим интеграл для верхнего коэффициента матрицы-столбца:LL N P dx +  N N P dx =Pi2i0iij j(1 -3хLL)30LPj х2+2L20Pj х33L3L=00Подставляя пределы и записывая результат в матричном виде, получим:Pi3=+Pj6=1 [2 1] {6Pi }PjВычисляя остальные коэффициенты, получим окончательное выражение для векторанагрузки произвольного конического стержня:FS(е) = h TOChLTOC21{Pi}(12.25)54[N]T dS =612PjSПример 12.2.

Вычислить распределение температур в стержне из примера 12.1, имеющего коническую форму, если температура большего по диаметру основания конуса постоянна и равна 150оС.X4 =7,5см ( R4=0,5; A4=0,25; P4 = )X3 =5см( R3=0,83; A3=0,71; P3=1,66 )X2=2,5см ( R2=1,16; A2=1,35; P2=2,32 )X1=0( R1=1,5; A1=2,25; P1=3 )Рис. 12.5Решение.1. Разбиваем стержень на три конечных элемента длиной по L=2,5см.2. Рассчитываем геометрические характеристики (Â(e), Рq, Аq, Rq, где q=1,…,4) – результатырасчета приведены на рисунке 12.5.3. Рассчитываем значения коэффициентов, входящих в выражения для матриц выделенныхконечных элементов (12.23 – 12.25):/L = 75[Вт/(смОС)]/2,5(см)=30(Вт/см2 ОС);hL/12 = 10[Вт/(см2 ОС)]2,5(см)/12=2,1(Вт/см2 ОС);hTocL/6 = 10[Вт/(см2 ОС)] 40ОС2,5(см)/6=166,7(Вт/см);4.

Вычисляем согласно (12.23) объемную составляющую матрицы теплопроводности элементов:KV(1) =KV(2)KV(3)=1301,8-11301,0= 300,48-11-1-11= 54-54-5454-11= 30-30-3030-11= 14,4 -14,4-14,4 14,45. Складывая полученные матрицы по методу прямой жесткости, получаем объемную матрицу теплопроводности всего стержня:KV = 54-5400-5484-3000-3044,4-14,400-14,414,46. В соответствии с выражением (12.24) вычисляем поверхностную матрицу теплопроводности элементов:55KS(1) =2,1(9+2,32)(5,32)(5,32)(9+2,32)= 23,8 11,211,2 21KS(2) =2,18,66447,32= 18,2 8,48,4 15,4KS(3) =2,162,662,664,66= 12,65,65,69,87. Складывая полученные матрицы по методу прямой жесткости, получаем поверхностнуюматрицу теплопроводности всего стержня:23,811,200KS(1) +KS(2) +KS(3) = 11,239,28,4008,428,05,6005,69,88.

К полученной матрице необходимо добавить поверхностный интеграл, взятый по площадиА4,=0,25 (см2). Используя выражение (12.19), имеем:K(3)S=А4 =00100,250100=02,59. Складывая все матрицы, приходим к общей матрице теплопроводности стержня:77,8-42,800KV +KS(1) +KS(2) +KS(3) +K(3)S=А4 = -42,8123,2-21,600-21,672,4-8,800-8,826,710. По формуле (12.25) вычисляем вектор нагрузки для каждого элемента:FS(1) = 166,7FS(2) = 166,7FS(3) = 166,7211221122112{32,32{2,32{1,661,661,0}= 166,7 {6+2,32}= 166,7 {4,64+1,66}= 166,7 {3,32+1,03+4,642,32+3,321,66+2,0}={13761273}}= {1050940}}= {720610}11.

К полученной матрице необходимо добавить поверхностный интеграл, взятый по площади А4= 0,25 (см2). Чтобы воспользоваться выражением (12.21), вычислим произведение:(hTOCА4) = 10[Вт/(см2 ОС)]40(oC)0,25(см2)= 100100 F(3)S=А4 ={01}= {010012. Приходим к системе уравнений:77,8-42,800T11376}56-42,8 123,2 -21,60T22323=0-21,6 72,4-8,8T3166000-8,826,7T471013. Решать данную систему уравнений есть смысл только для того, чтобы проверить правильность ее получения.

Характеристики

Список файлов лекций