Численный метод построения сглаживающих аппроксимирующих функций Текст научной статьи по специальности «Математика»

УДК 539.3

Ф. С. Хайруллин, О. М. Сахбиев

ЧИСЛЕННЫЙ МЕТОД ПОСТРОЕНИЯ СГЛАЖИВАЮЩИХ АППРОКСИМИРУЮЩИХ ФУНКЦИЙ

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

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

Key words: аpproximating function, numerical method, surface.

Article presents a numerical method for constructing the smoothing approximating functions. This function given by set of discrete points, which can be used to parameterize the middle surface and the boundary lines of thin shells. Also shows the results of calculations.

Параметризация срединной поверхности и граничных линий оболочки может быть произведена как аналитически, так и численно. Если оболочка имеет не каноническую форму, то ее срединную поверхность обычно задают сеткой дискретных точек, по которым строится интерполяционная или сглаживающая функция [1-3]. При использовании вариационного метода решения задачи в определяющем функционале имеются слагаемые, которые включают в себя частные производные второго порядка от искомой функции, которые с механической точки зрения обеспечивают минимальность «изгибания» функции. В данной работе предлагается численный метод построения сглаживающей функции, в основе условия «изгибания» которой лежит теория оболочек типа Тимошенко. Этот метод приводит к уменьшению порядка производных в функционале.

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

Предположим, что для искомой поверхности Q можно ввести поверхность отсчета Q0

таким образом, чтобы прямая, проведенная по нормали n 0 к поверхности Q0, пересекала Q

не более одного раза. В этом случае положение точки на поверхности Q можно представить в виде

r = Го +С По,

где r 0 = r 0(а 1, а2) - радиус-вектор точки М0 на поверхности отсчета Q0; r = r(а^ а2) -радиус-вектор точки М на поверхности Q, лежащей на нормали n0; £ = £(а.,, а2) -

расстояние между точками М и М0; а1, а2 - криволинейная гауссова система координат, заданная на поверхности Q0 .

Если система координат а1, а2 ортогональна на поверхности Q0 и координатные линии являются линиями главной кривизны, то коэффициенты первой и второй квадратичных форм искомой поверхности определяются по формулам [3]: aN = А* А* (бд +Л| Лj);

Ьц = _А^ kii = -Ai

А01 k0i

Л 2 5А01 dn i

V

А02 5а 2 5а 1

1 о 2; (1)

Ьі2 — — А\ А2 к 12 — — А-

Л і дА01 бЛ 2 ^ да 1

V А02 да 2

1 дС

Здесь АІ — Ао,- (і + С ко;) ; л і — — -; А0І, ко; - параметры Ляме и главные кривизны на

А і да,

поверхности О0; к11, к22, к12 - кривизны и кручение координатных линий на поверхности О; 5;у - символы Кронекера ( 5,, — 1; 5/у- — 0, если і ^ у ). Из формулы (1) видно, что в общем

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

Рассмотрим численный метод построения функции С — С(а1, а2) . Поверхность отсчета О0 разбивается на четырехугольные подобласти О0к с гладкими граничными линиями. Пусть в пределах подобласти О0к положение поверхности О задается совокупностью точек с координатами (а^, а£ у, С к), У — , где Ск - расстояние от точки с координатами (а к,-, а 2 у)

на поверхности отсчета О0к до поверхности О к по нормали п 0 .

Для построения сглаживающей функции предлагается использовать следующий функционал:

Ф(С, V, V 2) — Я<

1 д(А02^1)

А А

01 02

+

1 дС А02 да 1

да 1

ЙО0 + £ £р‘ [ , а2 у)—С к Г

2 1 д(А01¥ 2 ) 2 + р 7

1 2 С§ 2 0 % о А 1 V

¥1 +

_д^' А01 да2

+

(2)

к—1 у—1

где у1, у2 - функции, аппроксимирующие первые производные от искомой функции; р, рк -множители Лагранжа.

С механической точки зрения функционал (2) соответствует функционалу полной энергии изгиба оболочки с учетом поперечных сдвигов, в которой не происходит деформаций в срединной поверхности, коэффициент Пуассона у = 0. Функции у1, у2 являются углами сдвига нормали к срединной поверхности оболочки, множитель р определяет отношение жесткости на сдвиг к жесткости на изгиб оболочки. В этом функционале скомбинированы интерполяционные условия прохождения поверхности Ок вблизи заданных точек и условие минимального изгибания искомой поверхности. Так как в функционале (2) максимальными производными являются производные первого порядка, то для искомых функций требуется обеспечить только их непрерывность.

В качестве аппроксимирующих функций обычно используются кубические сплайн-функции [1-3]. Однако в этом случае требуется задание значений первых производных на границах области, что сделать довольно сложно, особенно для криволинейных границ. В работе [4] в качестве аппроксимирующих функций используются полиномы. Показывается, что в этом случае для получения удовлетворительных результатов не требуется выполнения граничных условий для первых производных.

В пределах подобласти О ok искомые функции представляются в виде [5]:

мм

С(оЦ, 0,2 )= X X Втп ^т() ^п(Р2 )’

тп 1т'

т—1 п—1 М N

¥1(а1, а2 )— X X Втп )т (Р1) )п ($2 );

т—1 п—1 М N

¥2 (а1, а2 )— X X Втп )т () )п ($2 );

т—1 п—1

2

О

0

2

где Б^п, Вт,, Втп - неизвестные постоянные; функции формы

)і(0і) — 1 — р1, )2(01) —р1, )т(01) — )Ж) [)2(01 )]т—2, т — 3М. (4)

Система координат 01, 02 вводится таким образом, что подобласть О0к отображается в единичный квадрат с равномерной координатной сеткой на границах. Представление решения в виде (3) позволяет выполнять условия стыковки подобластей и удовлетворять геометрическим граничным условиям. Для решения задачи используется вариационный принцип Лагранжа.

Весовые коэффициенты рк определяются из условия выполнения неравенств [2]:

^ ІІСК ■ а і)—с; Г <*«.

где вк - среднеквадратичное отклонение значений сглаживающей функции в узловых точках от узловых значений Ск для подобласти Ок .

Аналогичным образом строится сглаживающая аппроксимирующая функция для одномерного случая. Предположим, что на некоторой поверхности О задана кусочно гладкая кривая Х0 . Разобьем кривую Х0 на гладкие элементы Х0к, к — 1,К, представимые в явном виде: а 2 — ^ (а;).

Требуется построить кривую X, которая лежит на поверхности О и проходит через заданные точки (а^-, Ск), у — 1, Л, где Ск - расстояние от точки с координатами (а^, а2у) на линии Х0к до искомой линии X, отсчитываемое по нормали т 0 к линии Х0к. Таким образом, требуется построить одномерную функцию С — С(а;), которая удовлетворяет условиям С к — С(аку ), у — 1, Л, к — 1, К. Если кривая Х0 совпадает с координатной линией а;, то искомая функция С = а2 — С(а;) определяет уравнение кривой X.

Для определения искомой функции предлагается использовать следующий функционал:

Ф(С, у) = J

+ p

у + .

ds + £jrp; [)-Сj<, (5)

к=1 j=1

^ А01А02 1 / V А02 1

где у - функция, аппроксимирующие первую производную от искомой функции; А01, А02 -параметры Ляме, определенные на линиях Х0к ; р, р1к - множители Лагранжа.

В пределах элемента Х0к искомые функции представляются в виде:

м м

с(о,) = х ет1 ^(р,) у(о,) = хвт ^«в,); (6)

где Bm, Bm - неизвестные постоянные; функции формы tm имеют вид (4): З1 - безразмерная дуговая координата на кривой Х0к .

В таблице 1 представлены результаты аппроксимации функции y = sinx на отрезке [о, я/2], полученные на основании функционала (5) с использованием аппроксимирующих функций (б) и другими авторами.

Значения функции задавались с равномерным шагом, за линию отсчета принята линия х=0, полагалось, что К=1, J=5, M=J, s = 10-8. Чтобы показать влияние на результаты расчетов погрешностей, связанных с заданием значений f'(a) и f'(b), в таблице приводятся результаты, которые получены на основании сглаживающего сплайна и программы, опубликованной в работе Р]. При этом предполагалось, что на левом конце отрезка задается

З41

точное значение, т.е. f'(0) = 1, на правом конце может быть погрешность: f'(л /2) = '’ь. В четвертой строке представлены результаты, когда значения f'(а) и f'(Ь) вычислялись

1 f'(b) = —(3 fK - 4 fK- + fK_2)

Таблица 1

x=1.178 x= л / 2

Работа f'b f'( x) f''( x) f'( x) f''( x)

[2] 0 0.383 -0.925 0 -1.14

0.01 0.38 -0.955 0.01 -0.929

0.1 0.356 -1.16 0.1 -1.41

0(h2) 0.378 -0.961 0.015 -0.89

[4] - 0.384 -0.924 -0.003 -1.03

- 0.383 -0.92 -0.002 -1.02

Точное - 0.383 -0.924 0 -1

Из таблицы 1 видно, что при использовании кубических сплайнов результаты расчетов существенно зависят от граничных значений производных. Даже при точном задании значения f '(b) =0 результаты, полученные с помощью сплайнов, менее точны, чем с помощью полиномов.

При аппроксимации сплайнами в точках, находящихся между узлами, погрешность в определении производных функции может быть значительно большей, чем в самих узлах. При использовании же полиномов в качестве аппроксимирующих функций этого не наблюдается [4]. Для иллюстрации этого факта в таблице 2 приводятся значения производных f'(x), f''(x)

и их проценты расхождения от точных значений для функции y = sin x на отрезке [0, 2л],

вычисленные не в узловой точке. Предполагалось, что J=9, f0(x) = 0.8 sin x . Результаты

приведены для точки xc = 0.942 .

Таблица 2

Точное I II III IV

f'( x) 0.588 0.588 0.559 0.604 0.422 0.608

0% 6% 2.7% 28% 3.4%

f''( x) -0.809 -0.808 -0.770 -0.925 -1.267 -0.951

0.1% 4.8% 14% 57% 18%

В столбцах таблицы I, II, III, IV представлены значения производных, подсчитанные с использованием различных интерполяционных функций по точным узловым значениям. Длина отрезка разбивалась на 10 равных интервалов. Значения вычислялись для середины второго интервала в точке хс . Столбец I соответствует линейной интерполяции

f ' = ^2 +'з f '' = .

'с ~ 2 , 'с ~ 2 ’

столбцы II и III - интерполяции с использованием функций формы метода конечных элементов в виде одномерных функций Эрмита [6]:

242

f'(x) = f2'H,(x') + f2"H2(x) + /3'Нз(х) + fj"H4(x), f''(x) = [f'(x)] . f'(x) = [f2 H,(x) + f2'H2(x) + f, H3(x) + f3'H4(x)]', ](x) = [f'(x)]',

где

—2 —3 — —2 —3

H1(x ) = 1 - 3 x + 2 x , H 2 (x ) = x - 2 x + x , h3(x) = 3 x2 - 2 x3, h4(x) = x3 -x2, x = -x—^. x3 — x2

Здесь x2 = 0.628, x3 = 1.256 - граничные точки второго интервала,

f2' = 0.809, f2''=-0.588, f3' = 0.309, f3''=-0.951. Столбец IV соответствует интерполяции

кубическим сплайном [2].

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

В таблице 3 приводятся результаты аппроксимации цилиндрической поверхности,

заданной уравнением (z + 120)2 + x2 = 1302 в области 0 < x, у < 50. Поверхностью отсчета

являлась плоскость. Здесь 3 х 3, 5 х 5 - количество точек аппроксимации соответственно по

осям x и у; Д1 , R1 - коэффициент Ляме и радиус кривизны в окружном направлении цилиндра. При аппроксимации 3 х 3 значения функции задавались в точках с координатами x = 12,5; 25; 37,5; у = 0; 25; 50. При аппроксимации 5 х 5 значения функции дополнительно задавались в точках с координатами x = 5;45; y = 12,5; 37,5. Кроме этого выполнялись граничные условия:

V2 = 0 при x = 0, £ = 0 при x = 50, у1 = 0 при у = 0, у = 50.

При определении множителей Лагранжа полагалось s = 10-5. Вследствие симметрии в точках с одинаковыми координатами по оси x множители Лагранжа выбирались равными. При аппроксимации 5 х 5 в направлении оси x приняты следующие значения множителей Лагранжа: p12 3 = 4; р4 = 80; р5 = 160; p = 4,8 .

Таблица З

y x A R1

Точное 3 x 3 5 x 5 Точное 3 x 3 5 x 5

25 0 1 1 1 130 152 138

25 1,02 1,02 1,02 131 127

50 1,08 1,08 1,08 229 130

18,8 0 1 1 1 130 146 138

25 1,02 1,02 1,02 127 127

50 1,08 1,08 1,08 536 129

Как видно из таблицы, уже пятиточечная аппроксимация по каждой из координат позволяет получать удовлетворительные результаты не только по коэффициентам первой квадратичной формы, но и по коэффициентам второй квадратичной формы. Также следует, что погрешность определения значений не в узловых точках (строки у = 18,8 ) такая же, как и в узловых.

Результаты аппроксимации эллипсоидальной поверхности, заданной в параметрическом виде

x = 50 • sin p1 • cos p2;

<у = 60• sinp1 • sinp2; 0 <p1, p2 < 2,

z = 60 • cos p1;

представлены в таблице 4.

За поверхность отсчета принята сфера x = 50 • sin a1 • cos a2;

<у = 50 • sina1 • sina2; 0 <a1, a2 < K.

z = 50 • cos a1;

Значения функции задавались на равномерной сетке 5 х 5 в точках с координатами по оси a1 = 0,262; 0,524; 0,785; 1,05; 1,31; по оси a2 = 0; 0,393; 0,785; 1,18; 1,57. Выполнялись граничные условия:

V2 = 0 при a1 = 1,57 ; y1 = 0 при a2 = 0, a2 = 1,57 .

Для множителей Лагранжа приняты значения Pj = 12, р = 6, при этом s = 0,15 • 10-3 .

Таблица 4

a2 a1 A А r2

Точное Точное Точное Точное

0,785 0,02 60 60 1,2 1,2 50 49,8 53,3 50,7

0,408 59,2 59,2 23,5 23,5 53,1 51 54,4 50

0,795 57,2 57,2 41 41,2 59,5 59,9 56,5 55,2

1,18 55,2 55,2 51,7 51,9 65,2 64,9 58,3 60,1

1,57 54,3 54,3 55,2 55,3 67,3 69,9 58,9 61,7

1,18 0,02 60 60 1,2 1,2 57 56,4 45,8 46,7

0,408 59,7 59,7 23,7 23,7 57,9 56,6 46 42,5

0,795 59,1 59,1 42,4 42,5 59,9 59,6 46,5 44,1

1,18 58,4 58,4 54,6 54,8 61,8 61,2 47 45,8

1,57 58,2 58,1 58,8 59 62,2 64,5 47,2 46,7

Литература

1. Марчук, Г.И. Методы вычислительной математики. / Г.И. Марчук. - М.: Наука, 1980. - 536 с.

2. Григолюк, Э.И. Многослойные армированные оболочки. / Э.И. Григолюк, Г.М. Куликов. - М.: Машиностроение, 1988. - 288 с.

3. Корнишин, М.С. Вычислительная геометрия в задаче механики оболочек. / М.С. Корнишин, В.Н. Паймушин, В.Ф. Снигирев. - М.: Наука, 1989. - 208 с.

4. Якупов, Н.М. Расчет упругих тонкостенных конструкций сложной геометрии. / Н.М. Якупов, М.Н Серазутдинов. - Казань, 1993. - 206 с.

5. Хайруллин, Ф.С. Метод расчета тонких оболочек сложной формы / Ф.С. Хайруллин. // Известия РАН. Механика твердого тела. - 1998. - №3. - С. 30-33.

6. Норри, Д. Введение в метод конечных элементов. / Д. Норри, Ж. де Фриз. - М.: Мир, 1981. - 304 с.

© Ф. С. Хайруллин - д-р физ.-мат. наук, проф. каф. теоретической механики и сопротивления материалов КГТУ, [email protected]; О. М. Сахбиев - соиск. той же кафедры.