Добавил:
Опубликованный материал нарушает ваши авторские права? Сообщите нам.
Вуз: Предмет: Файл:

книги / Метод конечных элементов для решения электротехнических задач

..pdf
Скачиваний:
17
Добавлен:
12.11.2023
Размер:
877.09 Кб
Скачать

Отношения

x b и y a

 

позволяют записать функции формы

в безразмерном виде:

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

N1

 

1

 

 

(b x)(a y)

1

1

 

x

 

1

y

 

 

 

;

4ab

4

b

a

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

N2

 

1

 

 

(b x)(a y)

1

 

 

x

 

y

;

 

 

 

 

 

1

 

 

1

 

 

 

 

4ab

 

4

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

b

 

a

 

N3

 

1

 

 

(b x)(a y)

1

 

 

x

 

y

;

 

 

 

 

 

1

 

 

1

 

 

 

 

4ab

 

4

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

b

 

a

 

N4

 

1

 

 

(b x)(a y)

1

 

 

x

 

y

 

 

 

 

 

 

 

1

 

 

 

1

 

 

 

 

.

4ab

 

 

4

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

b

 

a

 

Здесь: 1

x

1 и

 

 

1

y

1.

 

 

 

 

 

 

 

 

 

 

 

 

 

 

b

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

a

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Если ввести следующие обозначения:

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

x

и

y

,

 

 

 

 

 

 

 

 

 

(2.37)

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

b

 

a

 

 

 

 

 

 

 

 

 

 

то функции формы можно записать в виде

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

N 1

(1 )(1 );

 

 

 

 

 

 

 

 

 

 

 

 

 

1

 

4

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

N2

1

(1 )(1 );

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

4

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

(2.38)

 

 

 

 

 

 

 

 

1

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

N3

(1 )(1 );

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

4

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

N4

1

(1 )(1 ),

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

4

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

где и являются локальными координатами, которые изменяют-

ся в диапазоне от 1 до 1.

Элемент слокальными координатами и показаннарис. 2.11.

21

((–1,1,11))

(1,(1,11)

4

3

 

 

 

 

1

2

((–1,1,–1)1)

(1,(1,–1)1)

Рис. 2.11. Локальная система координат

для прямоугольного элемента

Локальную систему координат можно использовать не только для прямоугольных конечных элементов, но и для непрямоугольных четырехугольных элементов (рис. 2.12). При этом ориентация сторон таких элементов относительно системы координат xy может

быть произвольной. На рис. 2.12 показан четырехугольный элемент общего вида и его локальная система координат и . При этом

формулы (2.38) идентичны функциям формы для указанного конечного элемента.

Для того чтобы найти частные производные по координатам x и y от функций формы, необходимо связать между собой системы

координат и xy .

 

 

 

0

 

0,5 (X3,Y3)

 

 

 

 

 

 

(X4,Y4)

 

3

 

 

 

 

 

 

4

 

 

 

 

 

 

 

0

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Y1

 

1 (X1,Y1)

 

 

2

 

 

 

 

(X2,Y2)

0

 

X1

 

 

x

 

 

 

 

Рис. 2.12. Четырехугольный элемент

22

Поскольку элементные вклады содержат, как убедимся позже, декартовы производные, необходимо преобразовать производные по и в производные x и y . Общее преобразование имеет вид

x y ;x y

x y .x y

Или в матричной форме

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

1

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

x

,

 

x

J

 

 

,

 

 

 

J

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

y

 

 

y

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

где J – якобиан преобразования координат,

(2.39)

(2.40)

x

 

y

 

 

 

 

 

 

 

 

 

J

 

 

,

(2.41)

 

x

 

y

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

а обратный якобиан определяется по формуле

 

 

 

 

 

 

 

 

 

 

y

y

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

1

 

 

 

 

 

 

 

 

 

 

 

J 1

 

 

 

 

,

(2.42)

 

 

 

 

det

 

 

 

 

x

 

 

 

 

 

 

 

J x

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

где det J

x y

 

x y

– определитель якобиана J .

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Соотношения (2.37) устанавливают связь между координатами x , y и , с началом координат в центре прямоугольного эле-

мента (см. рис. 2.10).

Для произвольного четырехугольного элемента, показанного на рис. 2.12, связь между координатами x , y и , определятся

уравнениями

23

x N1 X1 N2 X2 N3 X3 N4 X4 N

X ;

(2.43)

y N1Y1 X2Y2 N3Y3 N4Y4 N Y ,

 

 

 

где N1 , N2 , N3 и N4 – функции формы, определяемые выражения-

ми (2.38), описывают преобразование координат.

После подстановки в якобиан (2.41) уравнений (2.43) с учетом выражений для функции формы (2.38) получим

J1 (1 ) 4 (1 )

(1 )

(1 )

 

 

X1

Y1

 

 

(1 )

(1 )

X

 

Y

 

. (2.44)

(1 )

 

 

2

2

 

(1 )

X3

Y3

 

 

 

 

 

 

Y4

 

 

 

 

X4

 

 

Необходимо отметить, что функция, описывающая преобразование координат, может не совпадать по порядку с интерполяционной функцией, что позволяет описывать геометрию элемента независимо от аппроксимации неизвестной величины. Это позволяет сочетать как интерполяционные полиномы высокого порядка с элементами простой геометрии, так и элементы сложной формы с простыми интерполяционными полиномами. Если форма конечного элемента описывается полиномами той же степени, что и базисные функции, то такой конечный элемент называется изопараметрическим. Если степень этих полиномов меньше, чем у базисных функций, то элемент называется субпараметрическим, в противоположном случае – суперпараметрическим [3].

2.5. Интерполяционные полиномы для дискретизованной области

На рис. 2.13 представлена область исследования, составленная из пяти конечных симплекс-элементов. Нумерация элементов приведена в круглых скобках. Узлы конечно-элементной сетки также пронумерованы. Число узлов равняется шести. На рис. 2.13 величины U1 , U2 , U3 , U4 , U5 , U6 представляют собой значение функции

в узлах, а Xq , Yq – координаты узлов сетки ( q 1, 2, ..., 6 ).

24

Для локальной нумерации узлов элемента используются индексы i , j и k . На рис. 2.13 для второго элемента приведена локаль-

ная нумерация узлов. Звездочкой в каждом элементе обозначен локальный узел i . При этом локальная нумерация производится против часовой стрелки.

 

 

 

(

 

, X , Y

)

 

 

 

 

 

 

2

(U22,X2,2Y2)2

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

k *

 

 

 

 

 

 

 

(U4, X4, Y4)

 

 

 

 

 

 

 

 

 

 

(U4,X4,Y4)

* i

(2)

(1)

 

 

 

 

 

4

 

 

 

 

 

 

 

 

((U11,,XX1,1Y, 1Y)1)

 

j

 

((U3,,XX,3Y, Y)3)

 

 

 

 

 

 

 

1

 

 

 

 

3 3

3

 

 

 

(3)

 

3 *

(5)

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

(4)

 

 

 

 

 

 

 

5*

 

 

 

*

6

 

 

 

((UU5,5,X5,YY55))

 

 

 

((UU66,,X6,YY66))

Рис. 2.13. Пятиэлементная область

Связь между глобальной и локальной нумерацией узлов конечных элементов приведена в табл. 2.1.

 

 

 

 

 

Таблица 2 . 1

Связь между глобальной и локальной нумерацией

 

узлов конечных элементов

 

 

 

 

 

 

k

Номер элемента

 

i

 

j

 

 

 

 

 

1

1

 

2

 

3

2

 

4

 

3

2

3

 

5

 

3

4

4

 

6

 

3

5

5

 

3

 

6

1

Если для каждого конечного элемента в формулу (2.26) вместо

локальных индексов i ,

j и k

подставить соответствующие им гло-

бальные номера узлов из табл. 2.1, то получим следующую систему уравнений [3]:

25

u(1) N2(1)U2 N3(1)U3 N1(1)U1 ;

u(2)

N4(2)U4 N3(2)U3 N2(2)U2 ;

 

u(3)

N5(3)U5

N3(3)U3

N4(3)U4

;

(2.45)

u(4)

N6(4)U6 N3(4)U3 N5(4)U5 ;

 

u(5)

N (5)U

3

N (5)U

6

N (5)U

.

 

 

3

6

1

1

 

 

В уравнениях (2.45) интерполяционные функции определяются глобальными координатами и глобальными значениями функции в узлах. При этом конечные элементы объединяются в ансамбль.

Расширенная форма записи уравнений (2.45) имеет вид

u(1)

N (1)U

1

N (1)U

2

N (1)U

3

0U

4

0U

5

0U

;

 

 

1

 

 

2

 

 

3

 

 

 

 

 

 

 

 

6

 

 

 

u(2)

0U1 N2(2)U2 N3(2)U3 N4(2)U4 0U5 0U6 ;

 

u(3)

0U1

0U2

N3(3)U3 N4(3)U4

N5(3)U5

0U6

;

(2.46)

u(4)

0U1

0U2

N3(4)U3 0U4 N5(4)U5 N6(4)U6 ;

 

u(5)

N (5)U

1

0U

2

N (5)U

3

 

0U

4

0U

5

N (5)U

6

.

 

 

1

 

 

 

 

 

 

3

 

 

 

 

 

 

 

 

6

 

 

Уравнения (2.46) можно записать в матричной форме

u(1)

 

N (1)

N (1)

N (1)

0

0

0

U

 

 

u

 

 

 

1

2

3

N4

0

0

U

1

 

 

(2)

0

N2

N3

2

 

 

 

 

 

(2)

(2)

(2)

 

 

 

 

 

 

 

 

 

 

0

0

N3(3)

N4(3)

N5(3)

0

 

 

 

. (2.47)

u(3)

 

U3

 

u

(4)

 

0

0

N3

0

N5

N6

U

4

 

 

 

 

 

 

 

(4)

 

(4)

(4)

 

 

 

 

u(5)

 

N (5)

0

N (5)

0

0

N (5)

U

5

 

 

 

 

 

 

1

 

3

 

 

6

 

 

 

Таким образом, продемонстрировано, как с помощью заданного множества кусочно-непрерывных функций отдельных конечных элементов осуществляется аппроксимация искомой функции на всей области исследования. При этом для определения узловых значений функции необходимо применение других методов, таких как метод Ритца, методнаименьших квадратов, методГалёркинаи другие[1–5].

26

3.ПРИМЕНЕНИЕ МЕТОДА ГАЛЁРКИНА

ВМЕТОДЕ КОНЕЧНЫХ ЭЛЕМЕНТОВ

Внастоящее время галёркинский метод конечных элементов находит широкое применение при решении дифференциальных уравнений математической физики [1–5].

Рассмотрим дифференциальное уравнение

L f 0 ,

 

(3.1)

где L – дифференциальный оператор;

– искомая функция;

f

заданная функция.

 

 

 

Допустим, что приближенное решение задачи L f

0

запи-

сывается в виде

 

 

 

u NiUi

,

 

(3.2)

где Ni – базисная функция (функция формы); Ui

– узловые значе-

ния функции u .

 

Поскольку искомое решение приближенное, то для него будем

иметь

 

Lu f ,

(3.3)

где – невязка (ошибка) решения.

Чтобы невязка была минимальной, необходимо обеспечить ра-

венство нулю интеграла по объему V :

 

Wi dV 0 ,

(3.4)

V

 

где Wi – весовая функция.

 

В результате подстановки (3.3) в (3.4) получим

 

Wi (Lu f )dV 0 .

(3.5)

V

 

Условие (3.4) или (3.5) означает условие ортогональности невязки к выбранному базису.

27

В методе Галёркина весовые функции должны быть равны базисным функциям (функциям формы)

Wi Ni .

(3.6)

Тогда уравнение (3.5) запишется как

 

Ni (Lu f )dV 0 .

(3.7)

V

 

Ниже будут рассматриваться задачи, в которых высший порядок производных дифференциальных уравнений, описывающих тот или иной физический процесс, равен двум. При этом для аппроксимации искомой функции будут использоваться линейные интерполяционные полиномы (полиномы первой степени). Поэтому в уравнении (3.7) допускаются производные порядка не выше первого. Данная проблема решается снижением старшей производной выражения Lu путем интегрирования по частям.

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

 

 

i

 

 

 

 

 

0 .

(3.8)

 

 

N

(Lu f )dV

(e) V ( e )

 

 

 

 

Матричную форму записи полученной системы алгебраических уравнений относительно узловых значений искомой функции можно представить в виде

 

 

K U G F ,

 

(3.9)

где

[K ] [k](e)

– глобальная

матрица коэффициентов;

 

(e)

 

 

 

 

{F} { f }(e) – глобальный вектор-столбец

свободных

членов;

 

(e)

 

 

 

 

U G

– глобальный вектор-столбец узловых неизвестных;

k (e)

локальная матрица

коэффициентов;

{ f }(e)

локальный

вектор-

столбец свободных членов.

В результате решения системы уравнений (3.9) определяются узловые значения искомой функции.

28

4. РЕШЕНИЕ ОДНОМЕРНЫХ СТАЦИОНАРНЫХ ЗАДАЧ МКЭ

4.1. Теплопроводность плоской стенки

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

Рассмотрим процесс стационарной теплопроводности плоской бесконечной стенки толщиной (рис. 4.1) с заданными условиями теплообмена с окружающей средой на ее границах. Дифференциальное уравнение стационарной теплопроводности в трехмерной постановке в декартовой системе координат запишется как [6]

 

 

t

 

 

t

 

 

t

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

qV 0

,

 

 

x

 

 

 

 

 

x

 

y

 

y

 

z

z

 

где x ,

y ,

z – декартовыекоординаты, м; t

температура,

°С;

– коэффициент тепло-

проводности, Вт/(м·°С);

qV – удельная мощ-

ностьвнутреннегоисточникатепла, Вт/м3.

 

Будем считать, что в градиенты температуры по координатам y и z равны

нулю. Тогда

d

 

dt

qV

0 ,

(4.1)

 

 

 

 

 

 

dx

 

dx

 

 

 

t

0

x

 

 

Рис. 4.1. Плоская стенка

В результате применения метода Галёркина к уравнению (4.1) получим

T

d

 

 

du

 

 

N

 

 

 

 

 

qV dV 0 .

(4.2)

 

V

dx

 

dx

 

 

29

где N T – транспонированная матрица функции формы; u – при-

ближенное решение.

Уравнение (4.2) преобразуем в уравнение, порядок производной которого не выше первого.

Так как

d

 

T

 

du

T

d

 

du

 

d N T

 

du

 

 

 

N

 

 

 

N

 

 

 

 

 

 

 

 

 

,

 

 

dx

dx

 

 

 

dx

 

dx

 

dx

 

 

 

dx

 

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

T

d

 

du

 

d

 

T

 

du

 

d N T

 

du

 

N

 

 

 

 

 

 

 

N

 

 

 

 

 

 

 

.

(4.3)

 

 

dx

 

dx

 

dx

 

dx

 

 

 

dx

 

 

 

dx

 

Интегралотлевой иправойчастей выражения(4.3) запишетсякак

T

d

 

du

 

 

 

 

d

 

T

 

du

N

 

 

 

dx

dV

 

 

 

 

N

 

 

dV

 

 

 

V

dx

 

 

V dx

 

 

 

dx

 

 

 

 

 

 

d N T

 

 

 

du

 

 

 

 

 

 

 

 

 

 

 

 

 

dV .

 

 

 

 

 

 

 

dx

 

 

 

 

 

 

 

 

V

 

 

 

 

dx

 

 

 

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

 

d

 

T

 

du

 

T

 

du

 

 

 

N

 

 

dV

N

 

 

dS .

 

V

dx

 

 

 

dx

S

 

 

 

dx

(4.4)

(4.5)

Здесь

T

 

du

– интеграл по замкнутой поверхности

N

 

dS

S

 

 

 

dx

 

тела.

Тогда

 

 

 

 

 

 

 

 

 

 

T

 

du

T

 

du x

N

 

dS S0

N

 

 

,

S

 

 

 

dx

 

 

 

dx x 0

30

Соседние файлы в папке книги