Файловый архив студентов. Файловый архив студентов.
1267 вузов, 4649 предметов.
/ Регистрация
FAQ Обратная связь Вопросы и предложения
Добавил:
Опубликованный материал нарушает ваши авторские права? Сообщите нам.
Вуз:
Пермский национальный исследовательский политехнический университет
Предмет:
[НЕСОРТИРОВАННОЕ]
Файл:

1323

.pdf
Скачиваний:
3
Добавлен:
15.11.2022
Размер:
16.27 Mб
Скачать
☆

Следующий шаг состоит в делении оси х на отдельные отрезки (числом от 50 до 100) и вычислении длин дуг, соответствующих каждому отрезку оси, по формуле

d% = У (Xj—Xt)2 + (Yj —Yf)2.

(16.3)

Сложение длин всех этих дуг дает длину дуги кривой в целом. Итак, в случае криволинейной границы возникает необходи­ мость определения координат одной или двух точек, которые де­ лят кривую на два или три участка равной длины. Расположение этих точек может быть определено в результате запоминания на­ капливаемых приращений длины дуги по х. Как только полная длина дуги станет известна, сравнительно просто пройти отрезок в обратном направлении и найти две точки на оси х, между ко­ торыми должен располагаться узел. Далее, считая изменение функции между этими двумя точками оси х линейным, определя­ ем х-координату узла. После определения х-координаты узла, его

у-координата вычисляется по формуле (16.2).

Задачи

169. Напишите программу для определения длины кривой, ко­ торая аппроксимируется кубичным одномерным элементом. Про­ верьте программу, рассматривая кривую у = 2 + х2 на отрезке от х = 1 до х = 2 . Сравните результаты вычислений по программе с результатами, полученными по формуле (16.1).

ВЫВОД УРАВНЕНИЙ ДЛЯ ЭЛЕМЕНТОВ С ПОМОЩЬЮ МЕТОДА ГАЛЁРКИНА

Уравнения для элементов, которые использовались в главах прикладного характера, были выведены в гл. 5 путем минимизации либо энергии деформации, либо функционала, связанного с рас­ сматриваемым дифференциальным уравнением. Существуют дру­ гие способы получения уравнений для элементов. Преимуществом этих способов является то, что отправной точкой для них служит непосредственно само дифференциальное уравнение и, кроме того, они исключают необходимость вариационной формулировки фи­ зической задачи. Один из способов, известный как метод Галёркина, был предложен в 1915 г. Галёркиным как приближенный метод решения краевых задач. В сочетании с интерполяционными соотношениями метода конечных элементов метод Галёркина весь­ ма эффективен при решении как краевых задач, так и задачи Коши.

В этой главе обсуждается применение метода Галёркина к решению дифференциальных уравнений первого и второго поряд­ ка. Эти области приложения достаточно хорошо иллюстрируют реализацию метода Галёркина, а также показывают возможность его применения в других прикладных областях.

17.1. Метод Галёркина

Метод Галёркина обсуждался многими авторами. Применения этого метода в сочетании с конечно-элементной моделью рассмат­ ривались в работах [2, 3, 6]1 Подробное изложение его содержит­ ся также в работе [4]. С помощью метода Галёркина получается приближенное решение дифференциального уравнения. При этом должно выполняться следующее условие: разность между прибли­ женным и точным решениями должна быть ортогональна функ­

циям, используемым при аппроксимации.

 

Если исходить

из

дифференциального уравнения Lu—f = О

(L — дифференциальный

оператор) и приближенное решение

ис­

кать в виде w =

то для него будем иметь Lu—/ = е, где

е —

остаток или ошибка, поскольку это решение только приближен­ ное. Мы хотим сделать е малой величиной, насколько это воз­ можно. При использовании одного из способов сделать е как мож-1

но меньшей величиной требуется выполнение равенства ] NiedR =

R

= 0 для каждой из базисных функций Ni. Это равенство матема­ тически означает, что базисные функции должны быть ортогональ­ ны ошибке по области R.

Применение метода Галёркина в сочетании с методом конеч­

ных элементов приводит к уравнениям

 

(17.1)

 

J адФ)^=о, р=/,и к

 

R

 

 

 

 

где ср — искомая величина,

которая аппроксимируется

соотноше­

нием

Ф= [ N l t N j , N k,

.]{Ф),

(17.2)

 

a L(f — дифференциальное

уравнение,

определяющее

ср. Пусть

£ (ф )=

- S -

+ 3 Т

+ 4 = 0 >

Ф '(°)=0,

тогда формула

(17.1)

принимает вид

 

 

 

L

 

 

 

 

 

 

 

 

dx —0.

(17.3)

Верхний предел интегрирования равен длине одномерной области, в которой ищется решение.

Высший порядок производных, которые могут содержаться в L(ф), не ограничен, он определяется физическим смыслом задачи. Однако высший порядок производных, который допустим в (17.1), на единицу больше порядка непрерывности интерполяционных со­ отношений. Так как все интерполяционные соотношения в этой книге нулевого порядка непрерывности (непрерывна ф, но не ее первая производная), в уравнения (17.1) могут быть включены производные порядка не выше первого. Это ограничение может вызвать некоторые сложности, но большинство из них можно пре­ одолеть сокращением порядка L(ф), используя интегрирование по частям.

17.2. Изгиб балки

Уравнение упругой линии балки

d 2y

_

М *

(17.4)

dx2

~

EI

 

служит удобной отправной точкой для данного обсуждения; оно описывает одномерную задачу, и, кроме того, это одно из уравне­ ний, которые наиболее часто рассматриваются инженерами. На­

помним,

что

М — изгибающий момент

в произвольной

точке х

балки, Н-см,

Е — модуль

упругости,

Н/см2, / — момент

инерции

сечения,

см4,

и

у — прогиб

балки, см. Будем предполагать, что

М — известная

функция

координаты

ху тем самым М/EI будем

считать заданной величиной.

 

 

 

Применение

метода

Галёркина к

уравнению (17.4)

приводит

к соотношению

 

 

 

 

 

 

 

 

 

| m

r ( - 0 - - - g - ) * =

O .

(17.5)

 

 

 

о

 

 

 

 

 

Интерполяционная функция для у определена на отдельном эле­ менте, поэтому уравнение (17.5) должно быть переписано в виде суммы

я

dV e)

М<е)

(17.6)

И

dxz

EI

 

 

 

 

где R — число элементов,

a L(e) — длина отдельного элемента.

Прежде чем начинать вычисления, необходимо выбрать функции

формы [ЛДе)] и преобразовать интеграл (17.6) таким образом, что­

бы он содержал производные порядка не выше первого.

 

Ограничимся линейной моделью для у:

 

 

y — NiYt -\-NjY j=

1___JL

L

-Ш = [Л ^ > ]{ У ).

(17.7)

 

1

L г

 

 

Момент М — функция

длины

элемента,

величина M/EI1'1 может

быть аппроксимирована также

с

помощью линейной модели:

-% -= [#<*>]

MJEI1

(17.8)

 

 

Mj/EI J

Интеграл

преобразуем путем интегрирования по частям, что дает

J [W(e)]r - 0 d x = [ N (e)] T

| 1—J d[N^ ]T-J - d x . (17.9)

J) Эта величина представляет собой кривизну балки. — Прим, ред.

Подставляя (17.9) в (17.5), получаем в качестве исходного урав­ нения соотношение

*/

— f

 

 

[NM\T JlL

J L

+ [ N M\T J L \ d x = 0. (17.10)

d x

J V d x

dx

E l I

x i

L (e )

 

 

Первое слагаемое в интеграле (17.10) представляет собой мат­ рицу коэффициентов элемента [&(е)] в уравнении

[£<*>] {У} = {р'}.

(17.11)

После суммирования по элементам второе слагаемое в интеграле (17.10) будет соответствовать вектор-столбцу {F}. Член вне ин­ теграла в (17.10) вносит вклад в вектор-столбец {F} при усло­ вии, что производная dy/dx (наклон) определена на каждом кон­ це элемента. Этот член не учитывается, если ничего не известно о наклоне в узловых точках.

Вычислим интегралы в (17.10):

или

(17.12)

Верхний индекс (е) отброшен, так как рассматривается отдель­ ный элемент. Интегрирование выполняется в пределах от нуля до L, так как функции формы в (17.7) выражены в местной системе координат с центром в i-ы узле. L — длина элемента. Для второ­ го интеграла имеем

L

MJEI 1

j W I N ) (17.13)

Mj/EI J d x = J t

0

Применение полученных соотношений иллюстрируется на следую­ щем примере.

Пример

170.

Консоль длиной 150 см жестко

закреплена в точке

х = 0 и

подвержена действию сосредоточенной

силы величиной

4500 Н на свободном конце. Изгибная жесткость сечения EI рав­ на 8,5-10е Н/см2. Требуется определить прогиб консоли в шести узловых точках, расположенных на равном расстоянии друг от друга.

Разбиение области на элементы, а также эпюра изгибающего

момента показаны

ниже. Значения

величины М/EI в узлах даны

в следующей таблице:

 

 

Узел

M /EI

Узел

M /EI

1

—0,000794

4

—0,000318

2

—0,000635

5

—0,000159

3

—0,000476

6

0 , 0

К задаче 170. Е/ = 8,5-10® Н-см2; длина каждого элемента 30 см.

Запишем уравнения, соответствующие первому элементу:

—

1

— Г

1

30

' 2

Г

MJEI

Г

i

Yi ]

 

 

 

 

30

1

1 . \ у , \ 1

 

 

 

Mj/EI

 

— 1

6

. 1

2

 

 

>—

г

dy_

 

 

 

 

X

dx

Iх= 0 = { o l -

 

 

L

 

 

 

 

Последний член в левой части уравнений исчезает, поскольку в точке х = 0 угол поворота равен нулю. Этот член не появляется в других уравнениях, потому что ничего не известно об угле по­

ворота ни в одной из других узловых точек. Коэффициенты в уравнениях для остальных элементов остаются теми же самыми, так как все элементы имеют одинаковую длину. В результате суммирования матриц элементов получаем следующую систему уравнений:

' 1

—1

' Ч

—1

2 — 1

У ,

 

—1 2 —1

У *

 

—1 2 — 1

у *

 

—1 2 —1

Уь

 

— 1 1

У .

 

 

 

"2

1

 

 

 

 

—0,000794'

0

 

 

 

1

4

1

 

 

 

—0,000635

0

 

 

-j“ 150

1

4

1

 

 

—0,000476

0

 

 

 

1

4

1

 

—0,000318

—

 

 

 

 

 

 

0

 

 

 

 

 

 

1

4

1

—0,000159

0

или

 

 

 

 

 

 

1

2

0,0

л

 

 

 

 

 

 

 

 

 

 

- 1

—1

 

 

 

 

 

fKi]

 

0,33345

 

—1

2

— 1

 

 

 

 

 

0,57155

 

 

 

 

 

у *

 

 

 

—1

2

— 1

 

 

 

У 3

 

0,42855

 

 

 

— 1

2 — 1

 

 

у *

 

0,28605

 

 

 

 

—1

2

—1

 

У 5

 

0,14310

 

 

 

 

— 1

1

 

 

 

0,02385

 

Последняя система уравнений может быть решена без обра­ щения матрицы, потому что Yi=0. Из первого уравнения

Кл—К2= 0,33345

получаем

У 2= — 0,33345.

Второе уравнение теперь может быть решено относительно Уз, третье — относительно Уд и т. д. Таким образом, для У могут быть записаны рекуррентные соотношения. Определяя п как глобаль­

ную степень свободы, будем иметь

- Y n_ 1 + 2 Y a - Y n+1= F n, п > 2,

ИЛИ

Y n+i = 2 Y n — —F „, п > 2.

Окончательные значения прогибов в узловых точках даны н и ж е. Эти значения очень хорошо согласуются с теоретическими резуль ­ татами.

У зе л

М е то д

Т е о р е т и ч е с к и е

 

М е т о д

Т е о р е т и ч е с к и е

кон ечн ы х

У з е л

к о н е ч н ы х

зн а ч е н и я

з н а ч е н и я

 

эл ем ен то в

 

э л е м е н т о в

 

 

 

 

1

0

0

4

—2,5719

—2,5729

2

—0,3334

—0,3335

5

—4,1929

—4,1929

3

—1,2385

—1,2388

6

—5,9550

—5,9559

В рассмотренной выше задаче были получены рекуррентны е

уравнения благодаря

наличию

граничных условий

у ( 0 ) — 0 и

у ' (0) = 0 того ж е типа,

что и в

задаче Коши. Если

граничные

условия состоят в задании прогиба в двух или более точках, та линейная система уравнений должна быть решена с помощьк>

процедуры, описанной в гл. 7. Последняя

компонента

вектора

{Г}

не существенна для данного примера,

потому что

значение

Ув

может быть

определено из пятого уравнения. Ш естое

уравне­

ние

системы не

будет выполняться, если в

него подставить най­

денные числовые значения У5 и Уб. Это объясняется тем, что, ре­ шая задачу Коши, мы произвольно урезали область задания ис­ ходного дифференциального уравнения.

17.3. Д вум ерны е уравнения теории поля

Ряд задач физики и техники может быть описан уравнением вида

L(q>) =

<?2Ф

д2<р

( 1 7 . 1 4 )

А*8

ду* + Q — 0 .

В предыдущих главах это уравнение было использовано при ре­ шении задач о течении жидкости, переносе тепла и о кручении упругого стержня. Здесь уместно обсудить применение метода Галёркина к решению уравнения (17.14).

Подстановка (17.14) в (17.1) дает

J № ( - 0 - + - 0 - + Q) d7==O-

( 1 7 1 5 >

V

Прежде всего необходимо преобразовать (17.15) в уравнение, со­ держащее только первые производные по * и у . Замечая, что

д

д2Ч> . d[N]T dq>

можно записать

 

д2<р

 

 

 

 

 

 

(17.16)

[NV дх2

 

 

 

 

 

 

Теперь первое

слагаемое

в

объемном

 

интеграле преобразуется

к виду

 

 

 

 

 

 

 

 

f № - 0 -

d v = \ -Л - П а/г ■ £ ) d

V

-

\ H

f f L - ^ - d V .

(17.17,

V

V

 

 

 

V

 

 

 

По теореме Остроградского — Гаусса имеем

 

 

 

 

f - J r ( W f ^

)

d V = U

N

V

^

l xd S .

(17.18,

V

 

 

s

 

 

 

 

 

Точно так же можно преобразовать интеграл

 

 

Объединяя соотношения (17.17) и (17.18) с аналогичными соот­ ношениями для приведенного выше интеграла и учитывая, что

dV =tdA и dS =td%6,

уравнение (17.15) можно переписать в виде

г ( 1 К + | г ' « ) ‘ й 5 -

m r Q ) dA = 4 . (17.19,

А

Здесь предполагается, что толщина элемента t равна единице.

Поверхностный интеграл в (17.19) может

быть

выражен че­

рез

величину

ду/дп,

где п — внешняя нормаль к

поверхности.

В результате имеем

 

 

 

 

 

 

[ ( д т т i!L

1 Ш Ш * L \ d A -[{N \T QdA —

 

 

 

J V

дх

дх

+

ду

ду I

J

 

 

 

 

А

 

 

 

 

 

А

 

 

 

 

 

 

 

 

 

 

 

—§1М\т^ г № = 0 .

(17.20)

 

 

 

 

 

 

 

%

 

 

 

Первый

интеграл

в

(17.20)

вносит

вклад в

матрицу [&(е)], вто­

рой— в

вектор {/<е)}

из уравнения

[•#*>] {ф} =

{/<с>}. Третий

интег­

рал участвует в образовании обеих величин [№е)] и {/(с)}. Очевид­ но, если ду/дп обращается в нуль на границе, то третий интеграл исчезает.

Неизвестная функция <р в уравнении (17.20) определяется со­ отношением

ф= [# ]{(£ },

так что

£ - [« 1 |Ф | "

5 - |Л,|1Ф>-

Подставляя полученные формулы в первый интеграл (17.20), имеем

что соответствует выражению

 

J [В]т[В]А4{Ф}

(17.21)

А

 

для задач теории поля.

 

Третий интеграл в уравнении (17.20)

заслуживает более де­

тального рассмотрения, поскольку образует конвективную матри­ цу >в задаче о переносе тепла. Предположим, что мы имеем ряд треугольных элементов вдоль вертикальной границы, как показа­ но на фиг. 17.1 (эта граница выбрана только из соображений простоты), и хотим вычислить интеграл

\N]T - ^ - d ^

(17.22)

25

вдоль этой границы. Поток тепла вдоль границы соответствует теп­ ловым потерям, вызванным конвективным теплообменом, и пред­

ставляется величиной

 

^ T = h ( <Р5-Фсе),

(17.23)

где фs — температура границы тела, а фоо — температура

окру­

жающей среды.

 

Температура внутри элемента дается соотношением

 

Ф= В Д + N p j + NkOk,

(17.24)

откуда имеем для точек поверхности

 
Соседние файлы в предмете [НЕСОРТИРОВАННОЕ]
Помощь Обратная связь Вопросы и предложения Пользовательское соглашение Политика конфиденциальности