Алгоритм расчета трехмерных электромагнитных полей электротехнических устройств методом конечных элементов Текст научной статьи по специальности «Математика»

где п. = 0,23 - номинальная процентная ставка (соответствующая депозитной процентной ставке банка); в = 0,08 - средний уровень инфляции.

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

1„

400

Полученные результаты по ЧДД > 0 и ИД > 1 свидетельствуют о том, что вложения дополнительных инвестиций в ПГУ экономически эффективны.

Проведенные практические технико-экономические расчеты по обоснованию теплофикационных ПГУ по сравнению с ПТУ на основе предлагаемых мер по совершенствованию методики сравнительной экономической эффективности таких расчетов позволяют сделать следующие выводы:

1) Известно, что оптимальные значения (Т4 /Т3 )опт для технико-экономических показателей ПГУ эг, эг т

и = ПТ4 /X,) не совпадают [1]. Поэтому расчет их

следует выполнять для ряда значений Т4/Т3 с целью получения минимального расхода топлива Вгт . Причем недовыработка электроэнергии у некоторых значений должна компенсироваться за счет замещаемой электрической мощности энергосистемы.

2) Сравнительные расчеты по энергетическим установкам должны выполняться в относительных единицах: при отпуске из отборов турбин 1 ГДж теплоты и выработке на ней соответствующего количества электроэнергии.

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

4) И, наконец, по выбранному варианту энергоустановки надо дать оценку экономической эффективности через интегральные показатели [2].

Нетрудно заметить, что несоблюдение этих правил может привести даже к обратному конечному результату.

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

Библиографический список

1. Соколов Е.Я. Теплофикация и тепловые сети: Учебник для вузов. — 7-е изд. — М.: Изд-во МЭИ, 2001.

2. Лисиенко В.Г., Щелоков Я.М., Ладыщев М.Г, Хрестоматия энергосбережения: справочное издание в 2-х книгах. - М,: Теплоэнергетик, 2003.

УСМАНСКИИ Юрий Тихонович, доцент кафедры «Теплоэнергетика».

Статья поступила в редакцию 24.08.06. © Усманский Ю. Т.

УДК ¿21.313.17

Ю. 3. КОВАЛЕВ Е. Г. АНДРЕЕВА Д. В. КОЛМОГОРОВ

Омский государственный технический университет

АЛГОРИТМ РАСЧЕТА ТРЕХМЕРНЫХ ЭЛЕКТРОМАГНИТНЫХ ПОЛЕЙ ЭЛЕКТРОТЕХНИЧЕСКИХ УСТРОЙСТВ МЕТОДОМ КОНЕЧНЫХ ЭЛЕМЕНТОВ

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

В декартовой системе координат уравнение Лапласа-Пуассона имеет вид:

У-А(х,у,2) = и(х,у,г),

(1)

где А(х,у,г) — вектор-потенциал магнитного поля, (х,у,2) е V (V - область моделирования), и(х,у,г)— правая часть, определяемая распределенной токовой нагрузкой. Кроме того, при решении задач магнитостатики на границе Г модели электротехнического устройства задаются однородные граничные

либо первого рода А(х,у,г) = 0, либо второго рода 5А(х,у,г)/ЭТТ = 0, п - нормаль к Г. Выражение (1) в декартовой системе координат (ДСК) для электромагнитных процессов будет иметь вид:

И,

'Э2А 8гА о А) Эу- ду-)

(2)

где ц, — относительная магнитная проницаемость материала модели, А — вектор-потенциал магнит-

Рис. 1. Трёхмерная регулярная сетка и трёхмерный элемент типа параллелепипед

ного поля, ] - вектор плотности тока. Уравнение (2) представляется тремя скалярными уравнениями относительно проекций вектор-потенциала А(х,у,г) на оси ДСК. Для их решения используется метод Галёркина в сочетании с методом конечных элементов (МКЭ). Суть метода состоит в нахождении аппроксимации функции Ах(х,у,г), Ау(х.у.г) и А2(х,у,г) значениями базисной функции внутри ограниченного объёма — конечного элемента, на которые разбивается весь объём моделирования (рис.1).

Разбиение области моделирования на КЭ по своему виду может быть регулярным и нерегулярным. При регулярном разбиении структура связей каждого конечного элемента (КЭ) с соседними КЭ является постоянной. Использование регулярного разбиения позволяет получить несколько преимуществ: упростить процесс разбиения трёхмерной области моделирования на КЭ, сделать этот процесс полностью автоматическим, упростить структуру данных для хранения параметров КЭ в памяти компьютера.

В методе Галёркина решение задачи находится исходя из условия ортогональности невязки к базисной функции Ы Дх.у.г) по области моделирования:

г1 (а2Ах | а'Я | а2А, м Эх2 ду2 ду2

Нс1У=о.

(3)

В итоге преобразования Галёркина уравнения Пуассона (3) приводят к системе интегральных уравнений:

X КГ-1

м

е-А "ах

V

5 А,

+----11

2 ^

с!У +

Базисная функция для узла КЭ находится из условия: = 1 в узле ш, в других узлах Ыт = 0 .

После преобразований последнего выражения получим:

V Ь /дКГ ФЦ , ¿КГ ,

^ м ] гЫ г* гкт

г=1 I Нг уД

эКГ

&

аг

дх дх ду ду

{Ав)с1У-ц0 рг[МтГ(1У = 0

Для элемента г типа параллелепипед, содержащего 8 узлов: ¡, к, 1, т, п, о, р, (рис. 1), можно записать:

где г — номер КЭ, Р - общее количество КЭ, на которые разбивается область моделирования; V,. — объём элемента г, [Ыт]г — базисная функция конечного элемента г, Ап — проекция вектор-потенциала магнитного поля на ось г), )п = ]л(х,у,г) — проекция вектора плотности тока на ось.

Приближенное решение на элементе г по МКЭ:

А=К]>Д,

где [Ыт]г — матрица-строка базисных функций, {Ат}г- матрица-столбец значений вектор-потенциала в узлах расчетной сети.

II, и, л,

П О Н-Г

агм. ЯМ,

дх ду дг

ЭЫ, дN¡

~дх

ЭЫ, ЭЫ, ЗЫ,

"Эх 9у дг

аг^

дх ду дг

эы0

дх ду дг

ЭМ„

дх ду дг

9Н„,

дх ду дг

тл ЗМ,,

дх ду дг

аы, агм. аы„ агч, ам„ ЯМ,,,

дх "аГ ах дх ах дх ах ах

згм, ам, амк ам„ аы„

'ду ~ду Эу ду ау ду ау ду

аы, аы. 314, аы^ аг^

. дг дг дг. СП. аг дг аг.

с1хс1ус1г

А,'

А) N.

А, N.

Ак Г Г Г, N.

А„ ООО N0

А,

Ат

V N.

Зхс1ус1г

или

li^l'.h,

Jjj—[B]rdxdydz-

П I» I) ^r

A, A| A, Ak A„

A, Am A,

-Hi

N,

Ni N, Nk N0

Np

Nm N„

dxdydz

(5)

Элементная матрица имеет размер 8x8 (64 элемента):

а12 а13

°21 а22 °23

«31 а32 а33

°41 а42 °43

<»51 а52 а53

Об1 °62 °63

а71 а72 ^73

,с81 а82 а63

«11 й2, а31 «41 а51 «61 «71 «81

а12 а22 а32 а42 а52 «62 а72 «82

а!3 а23 «33 а43 °53 РбЗ «73 «83

а, ,а,, +а12а12+а13а,3 а, ,а21 +а12а22+а13а23 а2 ,а,, + а22а, 2 + а23а, 3 а2 ,а2, + а22а22+а^а^

(6)

Матрица коэффициентов элементной СЛАУ выглядит следующим образом:.

[K]r= }}f-[B]rdxdydz.

О (I О

(7)

В результате объединения матриц [к], для отдельных КЭ получаем ленточную матрицу глобальной СЛАУ. Результирующая матрица имеет разреженный характер.

При подстановке в (4) конкретных базисных функций вывод коэффициентов [К]г представляет собой нетривиальную задачу. Базисные функции для КЭ типа параллелепипед будут иметь вид [2]:

N =

-х hv-y h -5

hx hy hz

;Nj =

hK-x _y_ h7 ^z

h„ hy "hT"!

\ -У x h,-z x у hz-z

N' ~' h. hv h,

N> = h„ h.

hx -x h„ -y z

N = --X— vT ; N = 7" hv nv h., P h

x h -y z

hy h,

JX_ у z

hx hy К

h,-x у z

hy К

Рис. 2. Конечный элемент типа «куб»

Программа для вывода коэффициентов СЛАУ типа параллелепипед, содержащего 8 узлов (рис.1), в Maple 9 выглядит следующим образом:

1. > with(LinearAlgebra):

2. Ni: = (hx-x) /hx"(hy-y)/hy* (hz-z)/hz;

3. Nj: = (hx-x)/hx*(y)/hy(hz-z)/hz;

4. N1: = (x)/hx*(hy-y)/hy*(hz-z)/hz;

5. ...

6. M: = Matrix(8,3,symbol = a);

7. a[l,l): = diff(Ni,x);

8. a[l,2]: = diff(Ni,y);

9. ...

10. a[8,3): = diff(Nn,z);

11. B: = M.Transpose (M);

12. MM: = Matrix(8,3,symbol = aa);

13. BB: = MM.Transpose(MM);

14. BB [ 1,1];

15. B[l,l];

16. В[8,8];

17. > for i from 1 to 8 do

18. for j from 1 to 8 do

19. printf("K[%d,%d] = ",i,j);

20. m(x,y,z): = B[i,j];

21. print(expand(int(int(int(m(x,y,z),x = 0.,hx), y = = 0..hy),z=0..hz)));

22. end do;

23. end do{

hx hz , hyhz , hx hy

24. K[l,l] =

25. K[l,2] =

26. K[8,8| =

hyhz

9 hx hxhz IS hx 9 hy hx hx

Thy

9 hz hx hy 18 hz

hy hz hx hy 9 hx 9hz

Поэтому для ускорения получения решения и исключения ошибок в преобразованиях для вывода коэффициентов элементной СЛАУ использовался пакет математических исследований Maple 9. Была разработана программа на встроенном языке описания формул.

В строке 1 подключается библиотека линейной алгебры (ЦпеагА1деЬга). В строках 2-4 задаются базисные функции для узлов КЭ. В строке 7 создаётся матрица М 8x3. В строках 7-10 — присвоение элементам матрицы М значений производных БФ. Строка 11 — присвоение матрице В значения транспонированной матрицы М. Строка 12 — создание матрицы ММ, строка 13 — создание матрицы ВВ из произведения матрицы ММ и транспонированной матрицы ММ. Матрицы ВВ и ММ служат для визуального контроля правильности умножения матриц. Строка 14 — вывод на экран элемента (1,1) матрицы ВВ. Строки 15-16 — вывод на экран элементов матрицы В.

Поскольку матрица В имеет 64 элемента, каждый из которых должен быть проинтегрирован по объёму, в строках 17-23 реализовано вычисление и вывод на экран значений матрицы В, проинтегрированных по х, у и г с использованием циклов для перебора

Таблица 1

1 2 3 4 5 6 7 8

1 h/3 0 0 -h/12 0 h/3 -h/12 -h/12

2 0 h/3 -h/12 0 -h/12 0 -h/12 0

3 0 -h/12 h/3 0 -h/12 0 -h/12 -h/12

4 -h/12 0 0 h/3 -h/12 -h/12 0 -h/12

5 0 ■h/12 -h/12 -h/12 h/3 0 -h/12 0

6 h/3 0 0 -h/12 0 h/3 -h/12 -h/12

7 -h/12 -h/12 -h/12 0 -h/12 -h/12 h/3 0

8 -h/12 0 -h/12 -h/12 0 -h/12 0 h/3

элементов. Последний приём позволяет избежать рутинного ввода формул коэффициентов искомой матрицы вручную (всего 64 формулы), кроме того, даёт свободу при необходимости модификации формул, В строках 24-26 показаны значения коэффициентов результирующей матрицы К (8x8) — их всего 64.

После проведенных расчетов можно получить элементные уравнения (матрицы) для конечного элемента типа «куб» (рис.2) с повторяющимися коэффициентами, довольно просто вычисляющимися:

к.11 к12 к13 к! 4 к!5 к16 к17 к18~ к21 к22 к23 к24 к25 к26 к27 к28

к71 к72 к73 к74 к75 к76 к77 к78 к31 к82 к83 к84 к85 к86 к87 к88

A1 Fl'

A2 F2

•

A7 F7

A8 F8

Ниже приведена таблица коэффициентов для элементной СЛАУ (табл.1).

Свободные составляющие вычисляются следующим образом:

= = Р6 = Я8 = - —-;

8

Р2 = РЗ = Р5 = ¥7 = .

8

Далее берется «регулярный элемент» из 8 кубов с 27 узлами и для него составляется уравнение из гло-бальной СЛАУ. Оно пишется относительно центрального узла и будет содержать 27 слагаемых. Для следующего «регулярного элемента» процедура по-в-торяется, причем по прежним правилам и рекуррентным соотношениям. Например, глобальный коэффициент (рис.1) у202=к26 = 0; у 102 = к25 + к46 =

= -h/12-h/12 = -h/6; у002 = к45 = -h/12; у201 = к23 + к86 = -h/12-h/12 = -h/6; уОЮ = k21 + k43 + + k85 + k76 = 0 + 0 + О-h/12 =-h/12; yl 11 = k22 + + к 11 + кЗЗ + k44 + k55 -I- кбб + k77 + k88 = 8*h/3 и т.д.

Кроме всего прочего пакет Maple 9 позволяет копировать полученные выражения для коэффициентов в текстовом виде прямо в редактор программы Borland С 4- -I- Builder, на котором написан программный пакет (ПП) АПЭМ.

Результаты вычислений рассмотренной программы могут быть использованы в коде программного пакета An3M-3D расчёта трёхмерных магни-тостатических полей.

Предлагаемая схема получения коэффициентов элементной СЛАУ для решения квазистационарных задач исследования электромагнитных процессов ЭУ на основе введенного «регулярного элемента» позволяет по рекуррентным соотношениям формировать уравнения глобальной СЛАУ, тем самым сокращая ресурсоемкость программного алгоритма.

Библиографический список

1. Андреева Е.Г., Ковалев В.З, Математическое моделирование электротехнических комплексов: Монография/ Под общ. ред. Ю.З. Ковалева. - Омск: Изд-во ОмГТУ, 1999. -172 с.

2. Зенкевич О., Морган К. Конечные элементы и аппрок-с-мация: Пер. с англ. - М.: Мир,1986. - 318 с.

3. Сегерлинд Л. Применение метода конечных элементов: Пер. с англ.. - М.:Мир, 1979.-392 с.

КОВАЛЕВ Юрий Захарович, доктор технических наук, профессор, заведующий кафедрой «Электрическая техника».

АНДРЕЕВА Елена Григорьевна, доктор технических наук, профессор кафедры «Электрическая техника».

КОЛМОГОРОВ Дмитрий Викторович, инженер ЦИТ и ДО.

Статья поступила в редакцию 17.08.06. © Ковалев Ю. 3., Андреева Е. Г., Колмогоров Д. В.