CC BY
36
УДК 62-278:621.184.64 Доцент С.А. Подгорный
(Кубанский гос. технол. ун-т) кафедра автоматизации производственных процессов
Оценка точности аппроксимации, полученных методом конечных элементов
В работе дается оценка точности аппроксимаций, полученных методом конечных элементов, и показано, что в случае бесконечного числа конечных элементов можно получить аб -солютную точность. Найденная зависимость позволяет для заданной точности определить необходимое число конечных элементов.
The paper assesses the accuracy of the approximations obtained by finite element method, and it is shown that in the case of an infinite number of finite elements can be absolute accuracy. The obtained dependence allows determining the required number of finite elements for a specified accuracy.
Ключевые слова: комбинационные составляющие, цепные дроби, ряд Фарея
В настоящее время метод конечных элементов и получил широкое применение для решения разнообразных задач математической физики [1,2]. Одной из важных задач в анализе континуальных систем является преобразование непрерывных функций в набор конечных элементов [3]. Метод конечных элементов позволяет строить удобную схему алгебраических уравнений относительно узловых значений искомой функции. При этом приближенная аппроксимация решения осуществляется на типовом полиномиальном элементе. В данной работе предлагается дополнить используемые полиномиальные элементы ступенчатой функцией Хэвисайда и функциями округления, что позволяет не только упростить и формализовать запись пробной кусочно-непрерывной функции, но и избежать интегрирования при ее использовании для решения континуальных задач методом Галеркина [4]. Специфическим в методе конечных элементов является построение семейства функций, определяемых ограничен -ным числом параметров. Например, семейство функций и(х) при Хтт < х < Хтах. Интервал Хтт ... Хтах представляет собой одномерную область существования решения задачи, которая разбивается на конечное число частей (элементов), соединяющихся между собой и с концами интервала в узловых точках (узлах) Хг В пределах каждого элемента задается функция виде линейного полинома. Она определяется своими значениями и(Х1) в узлах и на концах элемента. Учитывая, что в континуальной задаче функция является непрерывной, ее значения в каждом узле для соседних элементов
© Подгорный С.А., 2013
должны совпадать. Для этого введем функции округления: - функция пола, которая определяется как наибольшее целое, меньшее или равное х, а именно |_х] = п ^ х - 1 < п < х;
Гх!
- функция потолка, которая определяется как наименьшее целое, большее или равное х, а именно |_х] = п ^ х < п < х +1. Используя эти функции округления, введем семейство кусочно-линейных непрерывных функций следующего вида:
u(x) = {[ф(х - a)J - [ф(х - Ь)]}-
-(d -<
Ь - a
(1)
где Ф(х) - функция Хэвисайда, единичная ступенчатая функция, чьё значение равно нулю для отрицательных аргументов, единице для положительных аргументов и 1/2 в нулевой точке; а, Ь - интервал окна функции и(х), на котором она отлична от нуля (а < Ь); с, d - параметры уравнения прямой на интервале а, Ь.
Функция и(х) на границах интервала (а, Ь) обращается в ноль. Используем (1) для описания пробной функции на произвольном интервале, задаваемом полушириной Ь = (Ь - а)/2 одномерного единичного конечного элемента относи -тельно узла Хг В этом случае параметры кусоч-но-линейной непрерывной функции (1) соответственно равны а = Х1 - Ь, Ь = Х1 + Ь, с = 0, d = 1, а функция приобретает вид:
я (x )=Mx-X - h)]J-N* - X )!}• +
+ {ф(х - Х,)~|-|_Ф[х-(X, + h
h
-x + X, h
(2)
Функции ф1(х) изображаются в виде ломаных и определяются конечным числом параметров - своими узловыми значениями, отличными от нуля в пределах небольшого числа элементов вблизи узла Хг Поскольку
и(х) по своему физическому смыслу должна быть непрерывной функцией, выберем ф1(х) в виде кусочно-линейных функций, отличных от нуля на двух соседних элементах относительно выбранного. Каждая такая функция ф1(х), 1 = 0, 1, 2, ..., п-1, равна единице в XI и нулю во всех остальных узлах. При этом набор функций и(х) будет состоять из непрерывных функций, линейных в пределах элементов с изломами в узлах и определяемых своими узловыми значениями и;. Каждую из таких функций можно изобразить в виде ломаной линии. Метод конечных элементов заменяет задачу отыскания функции на задачу отыскания конечного числа ее приближенных значений в отдельных точках-узлах. При этом, если исходная задача относительно функции состоит из дифференциального уравнения с соответствующими граничными условиями, то задача метода конечных элементов относи-тельно ее значений в узлах представляет собой систему алгебраических уравнений. С уменьшением максимального размера элементов увеличивается число узлов и неизвестных узловых параметров. Вместе с этим повышается возможность более точно удовлетворить урав-нениям задачи и тем самым приблизиться к искомой функции. Рассмотрим применение метода для описания одномерного профиля. В этом случае пробная функция, определяемая уравнением (3), представлена линейной комбинацией функций (2) с коэффициентами
и — 1
и (X )=£ иг (X)
(3)
где п - число конечных элементов, аппроксимирующих континуальную зависимость на заданном интервале. Пусть на концах интервала Хшт...Хшах функция (3) обращается в нуль, а на рассматриваемом промежутке равна единице. В этом случае приравниваем функцию (3) этой величине:
и-1
X иг -?г(Х ) = 1
г=0
(4)
и, умножая левую и правую часть равенства (4) на ф^х), где} = 0, 1, 2, ., п-1:
п
Фу(хиг ) = 1 'Фу(х)
(5)
интегрируем систему уравнений (5) по области существования решения Хшт ... Хшах , согласно методу Галеркина получаем систему алгебраических уравнений, используемую для определения коэффициентов и;:
Г
Фу (хи1 х)
¿X = | ср у (х )&
(6)
При этом определенные интегралы имеют простые алгебраические выражения во всех узлах Х; на интервале существования решения:
0, г < j -1
Фу (хи -ФХх)
¿X =
-, г = у -1 6
2 и ■ ■ — ■ И, г = у
3
П . . ,
—, г = у +1 6
0, г > ] +1
| (Ру (х ^х = И
(7)
(8)
где Ь = (Хшт - Хшах)/(п + 1) - полуширина конечного элемента. Используя соотношения (7) и (8) в системе уравнений (6), получаем слева трехдиагональную разреженную матрицу коэф-фициентов, умноженную на вектор неизвестных {и;}, а справа - вектор из элементов Ь:
(9)
■2 И 3 " И 6 0 0 0
И 6 2 И И 6 0 0
0 И 6 2- И 3 И 6 0
0 0 И 6 2 И 3 И 6
0 0 0 И 6 2 И 3
и 0 И
и 1 И
и 2 = И
и„ - 2 И
и„-1 И
Умножая левую и правую часть матричного уравнения (9) на (1/Ь), получаем матрич-
ное уравнение, не использующее па]
2 3 1 6 0 0 0
1 6 2 3 1 6 0 0
0 1 6 2 3 1 6 0
0 0 1 6 2 3 1 6
0 0 0 1 6 2 3
и 1
и 2
и3 =
ип-1
ип
эаметр (Ь):
(10)
2 3 6 0 0 0
6 2 3 б" 0 0
0 т 2 3 т 0
0 0 6 2 3 6
0 0 0 2 3
и 1 и {X 0)
и 2 и (X„)
и 3 = и (X 2 )
и п _1 и (Хп - 2 )
и п и (Хп)
Очевидно, что вектор в правой части представляет собой значения аппроксимируемой функции, следовательно, в общем виде коэффициенты (ш) определяются следующим матричным уравнением:
(11)
Здесь матрица (11) является трехдиаго-нальной, и решение может быть получено методом Гаусса, матричным методом (умножением на обратную матрицу) или методом прогонки. Для систем малой размерности метод решения не существенен, но с увеличением точности решения необходимо увеличивать число конечных элементов на области существования решения, а это приводит к увеличе-
нию числа неизвестных в уравнении (11). Если требуется большая точность решения, необходимо использовать метод прогонки, позволяющий решать системы большой размерности без существенных ошибок округления. Рассмотрим применение этого метода в данном случае. Метод прогонки является двух шаговым. Вначале вычисляем вспомогательные величины ш, зависящие от элементов трехдиа-гональной матрицы из правой части (11):
1 1
а0 =--, а. =--
0 4 1 4 + а0
,«._, = - -
1
4 + а.
(12)
Рекуррентные соотношения (12) фактически представляют собой цепную дробь следующего вида (ряд Фарея):
4 -
(13)
4 -
Следовательно, (13) может быть представлена рекуррентной формулой, образованной следующей числовой последовательностью:
* д-' )2 -1 *, (14)
где д0 = 4; = 15. В этом случае элементы вектора (12) задаются следующими соотношениями:
1 * 0 * ¡-1 «0 =--, «1 = -—,- = (15)
Следующим этапом первого шага метода прогонки является расчет коэффициентов р1, которые для уравнения (10) (единичный прямоугольный профиль) определяются следующими соотношениями:
Р. = т, Р1 =
6 "^0
, ••• , Р1 =
6 - А-1
(16)
Используя соотношения Эйлера, преобразуем (16) к удобному для вычислений рекуррентному виду, образованному следующими числовыми рядами:
= 3, а1 = 6,
Ь0 = 2, Ь1 = 5,
Ь =
3 -(-1)' 2
3 -(-1)' 2
• Ь + Ь ¡_
(17)
(18)
к
1
4
а'-1 + а'- 2
а
а
0
Используя (17) и (18), получаем рекуррентную формулу для расчета элементов вектора (16):
^ = Ьт. л = ^ ■■■, = V (19)
и0 и1 и'
Определив из (15) и (19) значения вспо- шагу метода прогонки - расчету весовых ко-
могательных векторов, переходим ко второму эффициентов ш по формулам:
ип -1 = -4-^^ , ип - 2 = « п - 2 • ип -1 + Р п - 2 , - , и 0 = « 0 • и 1 + (20)
4 + а п - 2
Используя значения весовых коэффици- аппроксимацию этой функцией исходного
ентов и;, для пробной функции (3) получаем единичного прямоугольного профиля:
п-1
)=Т\ик -Фк (х п)]
к=0
где фк(х, п) определяется следующей формулой:
Фк (X, п)
2 • к + 2 1
Ф| X--+ 1
п +1
2 • к + 4 , Ф1 X--+ 1
2 • к 1
Ф| X--+ 1
п + 1
2 • к + 2 1
Ф| X--+ 1
п + 1
п + 1 п
, п X п • X
к------
X
- к + — + 2 2 2
2 2 п • X 1
(22)
Интегрирование функции (21) на интервале прямоугольного профиля дает следующую числовую последовательность:
4 6 _ 84
2 5 3 7 4 95
f _ 47 f _ 156
5 52 6 497
f_:
90
Таким образом, предложенный конечный элемент (22) для аппроксимации континуальных функций позволяет использовать числовые ряды для оценки точности описания континуальной функции, не прибегая к численному интегрированию. Рассмотрим данную числовую последовательность (23) в виде графика в осях {1/п - 1-Г(п)} (рисунок 1).
) с 1 0 2 ( ,3 ( ,4 (
Рисунок 1 - Экстраполяция числового ряда невязок -[ 1 -Т(п)] от обратного числа конечных элементов - (1/п)
Как видно из представленного графика (рисунок 1), уравнение линии тренда (Я=1) в точке ноль дает невязку, равную нулю, что соответствует числу конечных элементов стремящихся к бесконечности:
ЫЕ(п) = 1,2749(1/п)4 - 0,603(1/п)3 -0,347(1/п)2 + 0,5649(1/п)
Следовательно, описание континуальных функций конечными элементами может дать большую точность, а уравнение (24) является инфимумом для описания одномерных континуальных функций методом конечных элементов.
(24)
248
f. :
3411
10308
f.. :=
1286
7 97 8 265 9' 3620 10 ' 10879 11 ' 1351
(23)
ЛИТЕРАТУРА
1 Подгорный, С.А. Метод конечных элементов в решении задач теплопроводности [Текст] / С.А. Подгорный, З.А. Меретуков, Е.П. Кошевой, B.C. Косачев // Вестник ВГУ-ИТ. - 2013. - №2. - С. 10-15.
2 Kosachev, V.S. Using rounding function in problems of finite-element analysis [Text] / V.S. Kosachev, E.P. Koshevoy, S.A. Podgorny // Studies in mathematical science. - 2012. - V. 4. -№2. - P. 17-24.
3 Gandhi, D. Mixer spur analysis with concurrently swept LO, RF and IF [Text] / D. Gandhi, C. Lyons // Tools and techniques. -2003. - V. 46. - № 5. - P. 212.
4 Thomee, V. Galerkin finite element methods for parabolic problems [Text] / V. Thomee // Springer series in computational mathematics. - 2006. - P. 370.
REFERENCES
1 Podgorny, S.A. The finite element method in solving problems of heat conduction [Text] / S.A. Podgorny, Z.A. Meretukov, E.P. Koshevoy, V.S. Kosachev // Bulletin of VSUET. -2013. - № 2. - P. 10-15.
2 Kosachev, V.S. Using rounding function in problems of finite-element analysis [Text] / V.S. Kosachev, E.P. Koshevoy, S.A. Podgorny // Studies in mathematical science. - 2012. - V. 4. -№2. - P. 17-24.
3 Gandhi, D. Mixer spur analysis with concurrently swept LO, RF and IF [Text] / D. Gandhi, C. Lyons // Tools and techniques. - 2003. - V. 46. -№ 5. - P. 212.
4 Thomee, V. Galerkin finite element methods for parabolic problems [Text] / V. Thomee // Springer series in computational mathematics. - 2006. - P. 370.
