НАУКА и ОБРАЗОВАНИЕ

УДК 534.1

Россия, МГТУ им. Н.Э. Баумана

olgabatseva@yahoo.com

dim.sm2@yandex.ru

            При расчете динамики конструкций методом конечного элемента (МКЭ)  обычно рассматриваются два варианта матриц массы конечного элемента (КЭ) – диагональная матрица массы, характерная для многомассовой схемы, и согласованная матрица массы. При получении согласованной матрицы массы для аппроксимации распределения ускорений по КЭ используются  те же функции (их принято называть функциями формы), что и для аппроксимации перемещений при получении матрицы жесткости [1], поэтому матрица массы получается недиагональной. В классической монографии [1], посвященной методу конечного элемента, О. Зенкевич, обсуждая матрицы масс, пишет, «что на ранних этапах исследования динамических задач … масса элементов обычно предполагалась произвольно сконцентрированной в узлах, что приводило всегда к диагональной матрице,  даже если не существовало сосредоточенных масс».  И далее «данная процедура в действительности не нужна и приводит к плохой аппроксимации». Окончательный вывод: согласованная матрица масс «является единственно допустимой матрицей, используемой в расчете». Хотя с момента выхода монографии [1] прошло свыше сорока лет, метод сосредоточенных масс продолжает широко использоваться на практике. Во многих пакетах [2], реализующих метод конечного элемента, по умолчанию используется именно такой подход, хотя возможно и использование согласованных матриц. Одной из причин использования диагональных матриц являются, по-видимому, упрощения при реализации численных алгоритмов. Выигрыш, который это дает, может быть использован, например, для измельчения сетки и, следовательно, для компенсации дополнительных погрешностей, вносимых диагональными матрицами.

Возможно получение недиагональных и при этом несогласованных матриц массы. В этом случае для матриц массы и жесткости используются разные аппроксимирующие функции. О. Зенкевич [1] обсуждает этот вопрос в связи с учетом моментов сил инерции при изгибных колебаниях пластины, указывая, что «подобный подход применяется редко». Несогласованные матрицы масс могут использоваться и в других случаях. Так, в описании пакета программ NASTRAN[2] приводится матрица массы стержневого элемента, которая, как видно, не согласована по функциям формы с приведенной там же матрицей жесткости (как будет показано ниже, она равна полусумме диагональной и согласованной матриц). Расчеты с этой матрицей (приводятся так же и в настоящей статье) дают точность, значительно превосходящую результаты, даваемые стандартной  согласованной матрицей массы при определении собственных частот колебаний.

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

1) диагональная матрица, полученная по методу сосредоточенных масс;

2)согласованная матрица, для получения которой используются те же функции формы, что и для матрицы жесткости;

3)несогласованная матрица с функциями формы отличными от функций формы, которые использовались для  матрицы жесткости;

4)матрица массы, полученная из условия совпадения низшей ненулевой частоты колебаний с частотой колебаний свободного стержня;

5)матрица массы, являющаяся полусуммой диагональной и согласованной матриц.

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

1.Матрицы массы стержневого элемента

Матрица массы  двух узлового стержневого конечного элемента постоянного поперечного сечения может быть найдена по формуле [1], полученной с использованием принципа Даламбера - Лагранжа:

.                                              (1)

Здесь  - плотность материала стержня, - площадь поперечного сечения, -длина конечного элемента,  и  - функции формы - функции однозначно определяющие перемещения внутри и на границах конечного элемента через перемещения его узловых точек. Функции формы  принимают единичное значение в своем узле и нуль в чужом [1]. Для двух узлового стержневого элемента

, ; , .                                                         (2)

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

.                                                                                      (3)

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

1.1.Диагональная матрица массы

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

Рис.1 Многомассовая схема

Матрица имеет вид:

                                                                                          (4)

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

1.2.Согласованная матрица массы

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

 .                                                                                    (5)

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

Использование стандартных формул [1] для матрицы жесткости  и массы , полученных на основе принципа возможных перемещений, приводит к результату

,                                                                                          (6)

.                                                                                            (7)

Эта матрица массы дает так же точное значение массы элемента. Значения частот получаются строго завышенными [1]

1.3.Матрица  массы, полученная с использованием тригонометрических функций в качестве функций формы

Рассмотрим функции (в данной статье они подобраны эвристически)

; .                                                 (8)

Очевидно, что условия (2) и (3) для этих функций выполняются, и они могут рассматриваться в качестве функций формы стержневого элемента. Используем, как и ранее, стандартную формулу (1) для получения матрицы массы :

.

Найдем входящие в формулу интегралы

,

,

.

Соответственно матрица массы имеет вид

.                                                                                         (9)

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

1.4.Матрица массы, дающая точное значение собственных частот стержня

Получение этой матрицы основано на предложенном в статье [3] варианте метода синтеза форм колебаний [4-6], Метод синтеза форм колебаний (метод подконструкций в динамике) позволяет уменьшить число степеней свободы конечно-элементной модели конструкции благодаря использованию форм колебаний ее составных частей (подконструкций). В зависимости от используемых форм колебаний различают два основных варианта метода – метод жестких границ [4] и метод свободных границ [5, 6]. В пределе, без использования форм колебаний, метод жестких границ [4] переходит в метод Гайна-Айронса [7], в котом используются, так называемые «статические формы». Последние можно рассматривать как дискретный аналог функций формы. В работе [3] для исключения части степеней свободы предлагается использовать формы колебаний подконструкции, незакрепленной в граничных узлах. Причем в отличие от классического метода свободных границ [5] формы колебаний используются только для редукции матрицы массы. Для матрицы жесткости используются «статические формы», как и в классическом методе Гайна-Айронса [7]. При этом в [3] добиваются совпадения, как низших частот колебаний сокращенной динамической модели, так и форм колебаний в граничных узлах.  Сокращенная модель в отличие от [4-6] не содержит координат форм, а только перемещения в граничных узлах как в [7]. Применяя формулы статьи [3] и используя две низших частоты (первая равна нулю) и две формы колебаний стержня свободного на правом и левом концах приходим к следующей матрице КЭ стержня:

.                                                       (10)

Эта матрица также дает точное значение массы конечного элемента.

1.5.Матрица массы, используемая в пакете NASTRAN

По умолчанию в пакете MSCNASTRAN используется диагональная матрица, при задании дополнительной опции для стержня применяется матрица

 .                                                                 (11)

Нетрудно видеть, что эта матрица является полусуммой диагональной и согласованной матриц, приведенных выше: . Отметим, что в руководстве пользователя [2] матрица (11) называется согласованной, хотя согласно классическому определению [1], согласованная матрица массы – это матрица, для получения которой используются те же функции формы, что и для матрицы жесткости. Матрица (11), так же как и предыдущие (4), (7), (9), (10) дает точное значение массы конечного элемента.

2.Решение обобщенной проблемы собственных значений для конечно-элементной модели закрепленного на обоих концах стержня в дискретных тригонометрических рядах

            Сравнение пяти различных матриц массы проведем на примере стержня, закрепленного на правом и левом концах. Во всех случаях будем использовать матрицу жесткости (6). Дадим аналитическое решение этой задачи с помощью тригонометрического ряда, определенного на множестве узловых точек конечно-элементной модели

2.1.Диагональная матрица массы

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

.

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

Если учесть условия закрепления на правом и левом концах путем вычеркивания соответствующих строк и столбцов, получим матричное уравнение с матрицами меньшего размера (в нашем примере ):

.

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

 .                  (12)

Решение (12) будем искать в виде:

,                                                           (13)

где – номер тона колебаний, , (номера узлов соответствуют конечно-элементной модели без учета условий закрепления) - длина стержня, - длина конечного элемента. Очевидно, что . Поэтому, искомое решение можно переписать в виде

.                                                                    (14)

Подстановка (14) в любую строку матричного (12) уравнения дает

. Преобразуем последнее выражение, используя формулу для синуса суммы

.

Отсюда имеем

.                                                    (15)

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

2.2.Согласованная матрица массы

Рассмотрим конечно-элементную модель стержня при использовании согласованной матрицы массы(7) и матрицы жесткости (6). Будем искать решение в виде (14) и используем ту же процедуру, что и в подразделе 2.1. Опуская промежуточные выкладки, запишем

.                                                               (16)

2.3.Матрица массы, полученная с использованием тригонометрических функций формы

Конечно-элементная модель с использованием матрицы массы (9) и матрицы жесткости (6) дает следующие значения квадратов угловых собственных частот:

.                                                         (17)

2.4.Матрица массы на основе точного решения для стержня постоянного поперечного сечения

Для матрицы массы (10), полученной с применением метода синтеза форм [3] и форм колебаний незакрепленного стержня, применение процедуры подраздела 2,1 приводит к решению

.                                                            (18)

2.5.Матрица массы из пакета NASTRAN

Решение с использованием матрицы масс (11) и матрицы жесткости (6) для отдельного конечного элемента дает значение квадрата угловой частоты

.                                                                               (19)

3.Сравнительных анализ численных результатов

Результаты расчета частот с различными матрицами массы, полученными по формулам (1519), будем сравнивать с точным решением для стержня, закрепленного на правом и левом концах [8]:

.                                                                      (20)

3.1.Зависимость погрешности от числа разбиений

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

Рис. 2. Относительная погрешность при использовании диагональной матрицы массы

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

Рис. 3. Относительная погрешность при использовании согласованной матрицы массы

            Для согласованной матрицы массы значения частот строго завышенные. В области реально применяемых разбиений (от пяти точек на одну полуволну) нарастание погрешности с ростом номера тона и уменьшением числа разбиений – монотонное. Значения погрешностей в области реально применяемых разбиений примерно такие же, как и для диагональной матрицы.

Рис. 4. Относительная погрешность при использовании матрицы массы, полученной с использованием тригонометрических функций формы

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

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

            Значения частот получаются завышенными. В области реально используемых разбиений погрешность нарастает.

Рис. 6. Относительная погрешность, при использовании матрицы массы из пакета NASTRAN

            Частоты получаются заниженными, в области реально используемых разбиений по точности превосходит все другие варианты. Для очень грубых разбиений точность может быть даже ниже чем у двух предыдущих случаев (рис. 4, 5)

3.2.Сравнение результатов, полученных с использованием пяти матриц массы

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

Рис. 7. Сравнение относительных погрешностей при 10 узловых точках

Рис. 8. Сравнение относительных погрешностей при 15 узловых точках

            Из рис. 7, 8 видно, что при достаточном разбиении матрица массы из пакета NASTRAN дает наилучшую точность.

Выводы

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

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

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

- Наличие некоторой тестовой задачи, имеющей точное аналитическое решение

            - Наличие нескольких различных матриц массы (способы получения их рассмотрены в настоящей статье)

- Между погрешностями, которые дают эти матрицы (нужно хотя бы две разных матрицы) должно наблюдаться  устойчивое соотношение при различных размерах конечных элементов. Как показывают расчеты по приведенным в статье формулам, для нашего примера это имеет место. При числе разбиений от 3 до 14 узловых точек на одну полуволну отношение погрешностей для различных матриц к погрешности, которая получается при использовании диагональной матрицы (4) таковы: для матрицы (9) – 0,49, с максимальным среднеквадратичным отклонением около 9%; для матрицы (10) -  0,78 с максимальным среднеквадратичным отклонением не более 2%; для согласованной матрицы (7) - 1,0 с максимальным среднеквадратичным отклонением около 6%

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

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

Список литературы

1.               Зенкевич О. Метод конечных элементов в технике: пер. с англ. М.:Мир, 1975. 542 с.

2.               MSC. Nastran. Basic Dynamic Analysis User's Guide. MSC. Software Corporation, 2004.

3.               Дмитриев С.Н., Калугин И.Ю. Построение сокращенной динамической модели в стандартной форме метода конечного элемента // Аэрокосмические технологии, 2004-2007 : сб. тр. Всероссийских и Междунар. научно-техн. конф. (Реутов - Москва). М.: Изд-во МГТУ им. Н.Э Баумана, 2008. С. 117-121.

4.               Голдмен Р. Исследование колебаний методом расчленения // Ракетная техника и космонавтика (русск. перевод журнала American Institute of Aeronautics and Astronautics). 1969. Т. 7, № 6. С. 191-193. [Goldman R.L. Vibration analysis by dynamic partitioning // AIAA Journal. 1969. Vol. 7, no. 6. P. 1152-1154. DOI: 10.2514/3.5290   ].

5.               Крэйг мл. Р.Р., Бэмптон М. Сочленение подконструкций при динамическом расчёте конструкций // Ракетная техника и космонавтика (русск. перевод журнала American Institute of Aeronautics and Astronautics). 1968. Т. 6, № 7. С. 113-121. [Craig, Jr. R.R., Bampton V.C.C. Coupling of substructures for dynamic analyses // AIAA Journal. 1968. Vol. 6, no. 7. P. 1313-1319. DOI: 10.2514/3.4741  ].

6.               Крейг Р., Чжан К. Методы свободных границ для связывания субконструкций при исследовании динамики // Ракетная техника и космонавтика (русск. перевод журнала American Institute of Aeronautics and Astronautics). 1978.Т. 14, № 11. С. 154-156. [Craig R., Chang C.-J. Free-interface methods of substructure coupling for dynamic analysis // AIAA Journal. 1976. Vol. 14, no. 11. P. 1633-1635. DOI: 10.2514/3.7264 ].

7.               Гайан Р. Приведение матриц жёсткости и массы // Ракетная техника и космонавтика (русск. перевод журнала American Institute of Aeronautics and Astronautics). 1965. Т. 3, № 2. С. 287. [Guyan R.J.  Distributed mass matrix for plate element bending // AIAA Journal. 1965. Vol. 3, no. 3. P. 567-568.].

8.               Вибрации в технике. Справочник. В 6 т. Т. 1 / под. ред. В.Н. Челомея. М.: Машиностроение, 1978. 351 с.