3/2006
ДИСКРЕТНО-КОНТИНУАЛЬНЫИ МЕТОД ГРАНИЧНЫХ ЭЛЕМЕНТОВ ДЛЯ РАСЧЕТА СТРОИТЕЛЬНЫХ КОНСТРУКЦИИ
П.А. Акимов, А.Б. Золотов
Введение. Дискретно-континуальный метод граничных элементов (ДКМГЭ) применяется для определения напряженно-деформированного состояния (НДС) конструкций, зданий и сооружений с постоянством физико-геометрических характеристик по некоторому координатному направлению, называемому далее основным направлением. Примерами таких объектов являются балки, балки-стенки, тонкостенные стержни, плиты, оболочки, длинные фундаменты, высотные и протяженные здания, плотины, трубопроводы, рельсы и т.д. [1-4]. ДКМГЭ относится к классу полуаналитических методов. Полуаналитические постановки являются современными математическими моделями, реализация которых в наши дни стала возможной вследствие существенного роста производительности компьютерной техники. В основе ДКМГЭ лежат граничные псевдодифференциальные уравнения. Соответствующие псевдодифференциальные операторы аппроксимируются дискретно с привлечением анализа Фурье или вейв-лет-анализа. Преимущества ДКМГЭ перед другими методами численного моделирования заключаются в двукратном понижении размерности задачи (дискретизации подвергается не вся расчетная область, а только граница ее поперечного сечения, т.е. решается, по сути, одномерная задача и задается лишь шаг по контуру), возможности проведения детального анализа отдельных напряженных зон, упрощенном этапе подготовки данных, алгоритмиче-
ской простоте и высокой степени универсальности. Разработаны два варианта ДКМГЭ - непрямой (НДКМГЭ) и прямой (ПДКМГЭ), при этом непрямой, как и в стандартном МГЭ, применяется несколько шире прямого. В статье будут кратко рассмотрены примеры постановок некоторых краевых задач в рамках ДКМГЭ, изложены различные способы их численной реализации и приведены актуальные примеры практических приложений.
1. Примеры постановок краевых задач расчета конструкций в рамках ДКМГЭ. На этапах разработки и исследования ДКМГЭ авторами сформулированы постановки целого ряда актуального краевых задач расчета конструкций [4,11]. Построены операторные формулировки граничных псевдодифференциальных уравнений краевых задач для уравнения Пуассона, двумерной и трехмерной теории упругости, даны их регуляризации.
Пусть О - область, занимаемая конструкцией с границей 2 = дО; -дельта-функция границы 2; основное направление вдоль Х3. Используя технику метода расширенной области [1], область О дополняется до стандартной ш с продолжением реальных физических характеристик, решение на области ш ищется в классе функций, допускающих разрывы при переходе через 2. Тогда постановку второй краевой задачи трехмерной теории упругости можно представить в виде системы дифференциальных (по Х3) уравнений 1-го
порядка с операторными коэффициентами:
И' = ЬоИ + +5Н%, X еО ; /оИ = £ х е5 ,(1)
где и = [ йт Vх ]т ; и' = дзи ; V = й' = д3й ; 5; = д/дх1 ; (2)
" 0 Е " " 0 " " 0 "
Ь« = Т Ь 2 Т ¡1 ; г« =- Ь-1Ё ; Чо =- _ Ь^ J
Ч = -Ь2Ч ; (3)
Е - тождественный оператор; Ч - скачок естественных краевых условий при переходе через 2; Б = [Б1 Б2 Б3]т, £ = [£ Г2 £,]т - векторы составляющих нагрузок внутри О и на 2 соответственно; /« - оператор, задающий условия на 2;
V = V v2 У3]т - вектор составляющих внутренней нормали к границе (Уз=0);
5; =д/дх1, 1 =1, 2, 3 ; х = [х1 х2 х3]т - вектор координат точки;
Ь2 = к
0
1 0
0 0
у + 2
(у+2)д2 +д 2
(У + 1)д1д 2 0
Т =-Кг +1)
0 0 0 0 д, д 2
д,
д.
Я
; У = ~ К
(4)
Ц0 = -К (у + 2)д1 уд2 0 0 0 у
(у + 1)д:д 2 0
д? + (у + 2)д 2 0 0 д2 +д
1о = ^0,1 + ^/а,2 ; (5)
1а,1 =
д2 0
д1 0 0 0 0
; /0,2 = К
д
д.
0 0 0 0
2
уд1 (у + 2)д2 0 0 0 х 0 0 д 2 0 1 0
, (6)
0 д1 1 0 0 J _ „ „2
где Я, к - коэффициенты Ламе материала конструкции.
В соответствии с теорией операторов и методиками, описанными в [1]
строится фундаментальная матрица-функция системы И' = ЬоИ, далее производится ее свертка по х3 с обеими частями (1). После редукции имеем: И = ехр(- | V2 II х3 |)[Р1Г0 + ^п^)^ + х3Р2,0 + | х3 | Р^Б« + + ехр(- | V2 || х3 |)[Р1^0 + ^п^)^ + х3Р2,0 + | х3 | Р2,1]*(5Н~),
где Р1г0, Р1Г1, Р2Г0, Р2Г1 - редуцированные псевдодифференциальные операторы по хь х2, каждый представляется в виде суммы операторов вида дрд2 | V2 Р
( р,е2 ) с числовыми матричными коэффициентами; | V2|= д/-(д^ +д2) .
Подставляя (7) в граничные условия (1), получаем редуцированную разрешающую систему псевдодифференциальных уравнений относительно Ч :
¿V; ехр(-1V21| х3 »Ю^ + 81Еп(х3>д^11 + х3РГ2 0, +1 х31 Р^]*^~) =
2 _ (8) где = £-¿V;ехр(-1V2 ||х3 |)ЮГД1 + 81Еп(х3)д1>ц + х3дГ2,0Д +|х3|дГ2Лд]*Б0Г, х^2 + 0,
=1
Рм,к = «г,, 1 = 1,2; , = 0, 1; к = 1,2 .
(9)
Для дополнительной регуляризации ядер псевдодифференциальных операторов можно, воспользовавшись свойствами свертки, дважды перебросить произ-
водную на неизвестную функцию и перейти к отысканию = сЭ^ . Альтернативой является двукратное интегрирование граничных уравнений по х3.
Примеры постановок (в том числе регуляризованных) краевых задач расчета конструкций в рамках ПДКМГЭ имеются, например, в [4].
2. Исследование и регуляризация ядер основных псевдодифференциальных операторов ДКМГЭ. Вообще, проблема решения интегральных и интегро-диффе-ренциальных уравнений с ядрами вида х~к, | х |~к, к > 0 возникает при рассмотрении разнообразных технических задач. Указанные ядра не всегда могут быть вычислены в смысле Коши. Ядра типа 1п | х | и х-1 хоть и являются интегрируемыми в каком-то смысле, в некоторых точках принимают бесконечные значения, что приводит к достаточно сложным формулам численного интегрирования. Перечисленные функции правильнее понимать в обобщенном смысле, т.е. в виде их регу-ляризаций. Многие из существующих формул регуляризации (канонические, неканонические и т. д.) являются неоднозначными с точки зрения численной реализации, оказываясь полезными, главным образом, для теоретических изысканий. Естественной представляется регуляризация Ур(Г(х)) функции Ях) как производной соответствующего порядка от некоторой непрерывной, например,
Ур(1/хк) = (-1)к1[х(1п I X I -1)](к+" /(к -1)! ; Ур(1/ | X I) = [| X I (1п I X I -1)](2). (10)
Таким образом, регуляризованная обобщенная функция может быть представлена как последовательность конечно-разностньгх производных соответствующего порядка с параметром И от некоторой непрерывной функции, т.е.
Ур(Г(х)) = ИшБ8Р(х) , (11)
где Б(х) - непрерывная первообразная функции Д(х) Б-го порядка, Р(х) = Г( 5)(х); О8 - дифференциальный конечноразностный оператор Б-го порядка с шагом И.
Значения Ур(Д(х)) в достаточно удаленной от начала координат точке, должно практически совпасть с соответствующим значением функции Д(х), там можно использовать и саму Д(х). Альтернативным способом регуляризации сингулярных ядер является их аппроксимация рядами Фурье.
Ядра псевдодифференциальных операторов ДКМГЭ строятся с использованием прямого и обратного преобразований Фурье обобщенных функций. Например, в трехмерном случае для произвольной функции Д имеем:
| V21-1 ехр(-1 V2||хз |)Г(х1,х2) = (1/2п)[1/Я] * ^х^); к =1 х2 + х2 + х? . (12)
х1,х2 V 1 2 3
Авторами исследованы предельные свойства операторных ядер, даны их осе-симметричные интегро-дифференциальные представления [7].
3. Численные реализации ДКМГЭ.
3.1. Численные реализации ДКМГЭ с использованием рядов Фурье. Применение аппарата рядов Фурье при численной реализации ДКМГЭ ведет к сравнительной алгоритмической простоте и наглядности метода, а также обеспечивает «распад» разрешающей системы уравнений на набор подсистем меньшего порядка. Кратко рассмотрим без ограничения общности непрямой вариант ДКМГЭ применительно к трехмерной задачи теории упругости.
Для трехмерных задач принимается дискретно-континуальная модель границы нпоперечном сечении - сеточная аппроксимаиия гранииь,, а . продольном 2 Л
направлении (ось 0х3) задача остается континуальной и 2 разбивается на дискретно-континуальные граничные элементы (ДКГЭ) 2
Введем вдоль Г локальную координату X е [0, 1] (рис. 1, 2) и произведем на элементе локальную перенумерацию узлов: 1 ^ 1 ; 1 +1 ^ 2. В качестве основных неизвестных на ДКГЭ принимаются компоненты вектор-функции
~(х3) = Ч2 Ч3]т , обозначаемые ~ = [~11 Ч21 Ч31]т, 1 = 1,2, ..., ^е1 . В простейшем случае вдоль Г используется кусочно-постоянная аппроксимация и в пределах Г неизвестные считаются постоянными. Определяем:
со = { (х1, х2, х3) : - /1 < х1 < /1; - /2 < х2 < /2; - /3 < х3 < /3 } .
(13)
В качестве ортонормированного базиса в Ь3(ш) используется система функ-
N < к1 < N - N. < к2 < N. - N < к3 < N к = (к1, к2, к3)
ций ( 1 1 1 ; 2 2 2 ; 3 3 3 , ):
Фк(х) = ^к1(х1)^к2(х2)^к3(х3) ; Рк^) = ехр^х^/^ ; Як1 = /11 . (14) Для улучшения сходимости ряда Фурье используются а - множители Ланцоша.
Рис. 1. Примеры дискретно-континуальных моделей границ рассматриваемых областей в
рамках ДКМГЭ.
22
Рис. 2. ДКГЭ и его характеристики
Псевдодифференциальные операторы вида Р^, 1 = 1, 2; ] = 0, 1 и Qг,j,k, 1 _ 1,2;
J'
] = 0, 1; к = 1,2 аппроксимируются с привлечением рядов Фурье:
Р№,х2) - X XР1^к1к2^к1(х1)Фк2(х2>, 1 = 1, 2; ] = 0, 1; (15)
2' ^1Цк1к/к1к
к1 =- к2 =-Ы2
N ы2
Q[,j,kf(Xl,x2) - X XQLk,k1k2fk1k2Фk1(xlMc2(x2), 1 = 1, 2; j _ 0, 1; к = 1, 2, (16)
к1 =-N к2 =-Ы2
где ^ РГЛк2, QГ,j,k,klk2 - соответствующие коэффициенты Фурье.
Дискретно-континуальные аналоги уравнений (8) записываются в точках
х+_ (х+1,х+1,х3) е Г + 0, 1 = 1, 2, ..., ЫЬе1 с координатами
х
11 = х11 + 0.5 • Ь11 +У11ё1 ; х+1 = х21 + 0.5 • Ь21 +v2ld1 .
(17)
С позиций МГЭ, рекомендуемое значение dl = 0.01^, но di должна задаваться и так, чтобы избежать осложняющих факторов явления Гиббса, возникающего вблизи Г, и потенциально возможных паразитических эффектов, связанных с аппроксимациями дельта-функций, сосредоточенных на границе.
Глобальный вектор неизвестных формируем в виде
Г(хз) = [ ~1Т(хз) ... Щм(хз) ]Т; Г(хз) _ ¿^О^^з) . (18)
к3 =-Ы3
Набор 2Ыз+1 разрешающих систем линейных алгебраических уранений (СЛАУ) 2Ы^е1-го порядка относительно коэффициентов в (18) имеет вид:
К^ = О?, кз =-Ыз,...,Ыз .
(19)
Перемещения внутри области определяются формулами на основе (7):
N ы2 ы3
и(х1,х2,хз) = X X X ик1к2кз^к1(х1)^к2(х2)^кз(х3) ; (20)
к1 =-N1 к2 =-N2 кз =-N3
и _[и и и V V V ]Т
(21)
Определение напряжений внутри области производится на основе соотношений теории упругости, (20) и дифференцирования по Ланцошу.
Численная реализация ПДКМГЭ с использованием рядов Фурье, аналогичная той, что рассматривалась для НДКМГЭ и описана, например, в [4].
3.2. Численные реализации ДКМГЭ с использованием смешанной аппроксимации рядами Фурье и полиномами. Обоснованность выбора смешанной аппроксимации рядами Фурье и полиномами в рамках ДКМГЭ заключается в том, что аппроксимируемые функции могут содержать разрывы первого рода, а их приближение рядами Фурье не везде является удовлетворительным, вследствие явления Гиббса [14]. Для устранения этого фактора вводятся дополнительные полиномиальные составляющие в формулах аппроксимации, а рядами Фурье приближаются регулярные части решений, не имеющие разрывов. Авторами предложены схемы полиномиальной аппроксимации трех видов: двусторонняя, односторонняя справа и односторонняя слева (рис. 3), даны обоснованные рекомендации по выбору. На границе расширенной области w для полиномов задаются дополнительные условия гладкости [8].
Рис. 3. Схемы полиномиальной аппроксимации
Преимущества смешанной аппроксимации проявляются уже на модельной задаче об изгибе балки на упругом основании, для которой выполнено сравнение аналитического решения и решений, полученных по методу В.А. Даревского [12] и ДКМГЭ с использованием смешанной аппроксимации. Графики перемещений, углов поворота, изгибающих моментов и поперечных сил по длине балки при использовании ДКМГЭ совпадали с аналитическими, в то время как методика В.И. Даревского дала близкие результаты для функции перемещений, постепенно «ухудшаясь» по мере ее дифференцирования. Наиболее критичными ожидаемо стали графики поперечных сил (рис. 4).
ДКМГЭ с использованием смешанной аппроксимации распространен на решение задач теории упругости (рис. 4.5) [6].
Рис. 4. Сравнение аналитического решения и решений, полученных по методу В.А. Даревского и
ДКМГЭ для модельной задачи изгиба балки
Рис. 5. Расчетная схема бесконечной полосы на упругом основании
3.3. Численные реализации ДКМГЭ с использованием преобразований Фурье и часпшчной аппроксимации рядами Фурье. Априори использование интегрального преобразования Фурье по основному направлению кажется наиболее естественным. Этот выбор в сочетании с аппроксимациями рядами по «поперечным» переменным позволяет перейти от (8) к СЛАУ относительно образа Фурье
вектор-функции q(x3) [10]:
A(s)W(s) = B(s) , где W(s) = FX3 [q(x3)] . (22)
В дальнейшем для решения (8) появляется потребность в вычислении обратного преобразования Фурье от W(s) по s, что является главной проблемой при реализации. Актуальна задача рационального и численно обоснованного выбора набора относительно небольшого числа параметров s, при которых решается (22). С вычислительной точки зрения операция вычисления обратного преобразования Фурье является некорректной, некоторые методы ее решения описаны в работах В.И. Крылова, Л.Г. Кругликовой и др. [13].
Дискретное преобразование Фурье (алгоритм быстрого преобразования Фурье [10]) в ДКМГЭ может эффективно приметаться не на этапе дискретизации граничных псевдодифференциальных уравнений, а при определении перемещений и напряжений внутри области, обеспечивая выигрыш в скорости.
Вариант ДКМГЭ с использованием частичной аппроксимации рядами Фурье [10] предполагает применение рядов Фурье лишь по основному направлению и оперирование с точными аналитическими зависимостями по остальным. Преимущество - отсутствие погрешностей в граничных псевдодифференциальных уравнениях при аппроксимациях рядами дельта-функций и их производных, недостаток - вычисление коэффициентов разложения по приближенным формулам, для определенного количества коэффициентов соответствующий алгоритм связан с большим числом шагов интегрирования.
3.4. Численные реализации ДКМГЭ с использованием вейвлет-анализа. В ДКМГЭ вейвлеты привлекаются для построения базиса, с помощью которого производится многоуровневый анализ решений краевых задач, позволяющий качественно и количественно оценить степень локальности различных явлений. Под многоуровневым анализом понимается частичное разложение решения по локальному вейвлет-базису и рассмотрение компонентов решения на каждом из уровней базиса. Предложено две вейвлет-версии ДКМГЭ [9]. Первая заключается в комбинировании рядов Хаара [9] по основному направлению с аппроксимациями рядами Фурье по остальным. Разрешающая СЛАУ для граничных неизвестных имеет вид
AYH = В, (23)
где Y_ - вектор коэффициентов ряда Хаара искомой граничной вектор-функции; А, В - матрица коэффициентов и вектор правых частей СЛАУ
Вследствие более сложного алгоритма вычисления свертки в рядах Хаара порядок СЛАУ (23) равен M-Nbel-NHaar, где М - количество граничных неизвестных в точке; Nbel - количество ДКГЭ; NHaar - число учитываемых членов ряда Хаара. Матрица А является заполненной и здесь уже невозможна декомпозиция СЛАУ на набор подсистем меньшего порядка и возникает проблема решения СЛАУ большой размерности. По мере увеличения числа уравнений она становится все более критичной и приводит к временным проигрышам по сравнению с прочими варианта-
Рис. 6. Задача о расчете гравитационной плотины в трехмерной постановке.
ми ДКМГЭ. «Тяжеловесны» и формулы вычисления основных неизвестных внутри области. Метод эффективен при решении двумерных задач и трехмерных с не слишком густыми сетками ДКГЭ и не очень большим числом членов ряда Хаара.
Второй вариант ДКМГЭ, использующий вейвлет-аппроксимацию
предполагает приближение рядами Хаара по основному направлению с одновременным сохранением континуальных зависимостей по прочим переменным. Преимущества выбора: отсутствие проблем, связанных с явлением Гиббса и иными паразитическими эффектами на границе, меньший объем вычислительной работы по сравнению с первым вариантом, точное определение коэффициентов разложения ядер и, как следствие, ускорение при машинном счете, недостаток - проблема решения СЛАУ большой размерности.
4. Программные реализации
ДКМГЭ и примеры расчетов.
Описанный в статье ДКМГЭ для решения задач расчета строительных конструкций реализован в программном комплексе DCBEM 2D/3D [5]. Текущие версии вычислительных модулей составлены авторами на языке FORTRAN стандарта FORTRAN-90/95 c использованием Compaq Visual Fortran 6.6B Professional и Intel Fortran Compiler 7.0. Комплекс DCBEM 2D/3D функционирует в операционных системах Microsoft Windows 98/NT/2000/ ME/XP/ 2003. Помимо расчетных функций, программное обеспечение предоставляет возможности для различного численного анализа решаемых задач, предусмотрена диагностика и визуализация вводимых пользователем исходных данных и сервисное обеспечение.
В ходе исследовательских и хозяйственно-договорных работ решены модельные и практические задачи расчета конструкций. Для сопоставления ре-
зультатов использовались программные комплексы промышленного типа СТА-ДИО 2003 и Лп^/СЫШМ [2].
Фрагментарная галерея результатов, полученных при расчетах гравитационной плотины в трехмерной постановке показана на рис. 6.
Литература
1. Золотое А.Б., Акимов П.А. Некоторые аналитико-численные методы решения краевых задач строительной механики: Монография - М.: Издательство АСВ, 2004. - 200 стр.
2. Акимов П.А., Золотое А.Б. Численно-аналитические методы расчета строительных конструкций: перспективы развития и сопоставления. // САПР и графика, 2005, №1, с. 78-82.
3. Акимов П.А. Дискретно-континуальные методы расчета сооружений. // «НТТ - наука и техника транспорта», 2005, №1, с. 56-59.
4. Золотое А.Б., Акимов П.А. Прямой дискретно-континуальный метод граничных элементов для определения напряженно-деформированного состояния трехмерных конструкций. // «НТТ - наука и техника транспорта», 2004, №3, с. 70-77.
5. Акимов П.А. Программные комплексы, реализующие дискретно-континуальные методы расчета строительных конструкций, зданий т сооружений. // Вопросы прикладной математики и вычислительной механики: Сб. науч. тр. №8. - М.: МГСУ, 2005, с. 7-13.
6. Акимов П.А., Золотов А.Б. Дискретно-континуальный метод граничных элементов на основе смешанной аппроксимации рядами Фурье и полиномами для расчета бесконечной полосы на упругом основании. // Сб. докладов 13-го Словацко-Польско-Российского семинара «Теоретические основы строительства». - М., МГСУ, 2004, с. 103-108.
7. Акимов П.А., Золотов А.Б., Кай-туков Т.Б., МсхалаяЖ.И. Исследование
и регуляризация ядер основных псевдодифференциальных операторов дискретно-континуального метода граничных элементов. // Вопросы прикладной математики и вычислительной механики: Сб. науч. тр. №7. - М.: МГСУ, 2004, с. 15-23.
8. Золотое А.Б., Акимов П.А. Использование смешанной аппроксимации рядами Фурье и полиномами в рамках дискретно-континуального метода граничных элементов. // Вопросы прикладной математики и вычислительной механики: Сб. науч. тр. №7. - М.: МГСУ, 2004, с. 130-144.
9. Золотое А.Б., Акимов П.А. Численные реализации дискретно-конти-ну-ального метода граничных элементов с использованием вейвлет-анализа. // Вопросы прикладной математики и вычислительной механики: Сб. науч. тр. №7. - М.: МГСУ, 2004, с. 145-156.
10.Золотое А.Б., Акимов П.А., Кай-туков Т.Б., Ширинский В.И. Численные реализации дискретно-континуального
метода граничных элементов с использованием преобразований Фурье и частичной аппроксимации рядами Фурье. // Вопросы прикладной математики и вычислительной механики: Сб. науч. тр. №7. - М.: МГСУ, 2004, с. 157-167.
11. Sidorov V.N., Zolotov A.B., Akimov P.A. Discrete-continual Boundary Element Methods of Structural Analysis. // International Journal for Computational Civil and Structural Engineering. Volume 1, Number 5, Begell House Inc. Publishers & ASV, 2003, p. 84-99.
12.Дарееский B.M., Фридман А.З. Аналитический метод расчета тонкостенных элементов конструкций с применением ЭВМ. - СНИО, Николаев. 1989 г. - 77 с.
13. Крылов В.И., Крутикова Л.Г. Справочная книга по численному гармоническому анализу. - Минск: Наука и техника, 1968. - 165 с.
14.Ланцош К. Практические методы прикладного анализа. - М.: Гос. изд-во физ.-мат. лит-ры. 1961. - 524 с.
Публикуемая статья создана с использованием результатов выполнения работ на средства Гранта Президента Российской Федерации для государственной поддержки молодых российских ученых МД-1785.2006.8.