1323

ЭЛЕМЕНТЫ ВЫСОКОГО ПОРЯДКА.. МАШИННАЯ РЕАЛИЗАЦИЯ

В гл. 13— 15 обсуждались общие свойства элементов высокого порядка. Был рассмотрен только один пример использования этих элементов, а именно в задаче о переносе тепла в стержне был применен одномерный квадратичный элемент. Вопрос о том, как выполнить надлежащие расчеты с помощью ЭВМ, не обсуждал­ ся. Настоящая глава завершает рассмотрение элементов высокого порядка. Здесь будет описана машинная реализация указанных элементов, приведены три конкретных примера, а также будет показано, как определить координаты узлов, расположенных на криволинейных границах.

16.1. Машинная реализация

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

Отправной точкой в расчетах является выбор точек интегри­ рования и весовых коэффициентов для численного интегрирова­ ния. Число точек интегрирования зависит от порядка интерполя­ ционного полинома, который в свою очередь определяется тем, какой элемент используется при построении дискретной модели. Информация о типе и порядке (линейный, квадратичный, кубич­ ный и т. д.) элемента должна быть введена в ЭВМ до того, как начнется вычисление матриц элемента. Эта информация обычно вводится вместе с номерами узлов элемента. Порядок элемента должен быть определен при задании геометрии элемента и при интерполировании искомой величины по ее узловым значениям.

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

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

Составление матрицы Якоби, обратной к ней матрицы и вы­ числение ее определителя составляют первый этап работы цикла, в котором вычисляются коэффициенты матриц элементов.' Необ­ ходимые для составления матрицы Якоби частные производные вычисляются в подпрограмме, которая определяет как функции формы Ri, так и их производные dRi/dt, и dRi/dx\. Выбор соответ­ ствующего множества функций формы осуществляется с помощью условных операторов IF, устанавливающих тип и порядок эле­ мента.

После того как вычисление матрицы Якоби завершено, нужно^ решить, какими функциями формы будем пользоваться для интер­ полирования величины ф. Будут ли это те же самые функции формы, которые использовались для получения матрицы Якоби, или они будут другими? Если это те же самые функции формы, т. е. элемент изопараметрический, то Ni, dNi/dg, dNJdx\ и т. д. совпадают с Ru dRi/dl и т. д., использованными для получения матрицы Якоби, и можно приступить непосредственно к вычисле­ нию частных производных dNi/dx, dNi/dy и dNJdz. Если интерпо­ ляционный полином для ф отличается от полинома, который при­ менялся для задания формы элемента (субпараметрический или суперпараметрический элемент), то необходимо вычислить Ni, dNifdc,, dNi/dx] и т. д., прежде чем приступить к следующему шагу.

Частные производные по координатам х, у и z определяются в подпрограмме, которая выполняет умножение вектор-столбца, содержащего dNijd| и т. д., на матрицу Якоби. Зная Частные про­

умножением соответствующих значений подынтегральных выра­ жений на весовые коэффициенты и сложением полученных вели­ чин с результатами аналогичных операций, уже выполненных для других точек интегрирования. Поскольку [&(е>] и {/(е)} получаются путем суммирования, элементы этих матриц должны быть при­ равнены нулю перед началом работы внешнего цикла.

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

16.2. Примеры применения

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

16.2.1. Кручение стержня квадратного сечения

В гл. 7 были найдены сдвиговые напряжения в стержне квадрат­ ного сечения со стороной 1 см, подверженном действию крутящего момента величиной 196 Н-см. Из соображений симметрии рассмат­

симплекс-элементов (892 Н/см2) до применения теории согласо­ ванных результантов элементов. В случае разбиения области на четыре кубичных треугольных элемента максимальное значение хгу равно 929 Н/см2. Эта величина только на 1,7% меньше теорети­ ческого максимума 945 Н/см2.

Фиг. 16.3. Значения напряжения тZy в отдельных точках области (напряжение т2„ вычислялось в узловых точках). G=8-106 Н/см2, о=1° на 100 см, 2G0 = =2791 Н-рад/см3.

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

Главное отличие, которое получается при вычислении резуль­ тантов элемента с использованием симплекс-элементов и элемен­

тов высокого порядка, иллюстрируется на фиг. 16.4, где показано изменение величины хгу вдоль оси х. При использовании симплексэлементов напряжение получается постоянным в каждом элемен­ те, что соответствует ступенчатому изменению напряжения при пе­ реходе от одного элемента к другому. В случае применения тео­ рии согласованных результантов распределение хуг получается не-

Фиг. 16.4. Распределение по оси х значений х2у, вычисленных тремя различными методами.

------------ си м плекс-элем енты ;---------- — согласованные напряжения; — • — четыре кубич­ ных треугольных элемента.

прерывным, но оно отклоняется от истинных значений в центре и на границе области. Значения, полученные с помощью модели из четырех кубичных треугольников, тоже непрерывны, но они рав­ ны нулю в центре квадрата и правильно отражают распределение напряжений вблизи границы области. Значение в точке х = 0,25 см получено усреднением по трем значениям xyz.

Применение элементов высокого порядка уменьшает не только количество требуемых данных, но и размер результирующей си­ стемы уравнений. Десять уравнений были^ решены в случае одно­ элементной модели, для четырехэлементной модели решалась систе­ ма из 28 уравнений. В обоих случаях число уравнении меньше 4о, т. е. меньше числа уравнений, полученных при использовании симплекс-элементов. Кроме того, при использовании элементов высокого порядка отпадает необходимость в применении теории согласованных результантов элемента. Таким образом,^ исключает­ ся из рассмотрения еще одна система из 45 уравнении. В резуль­ тате сокращения числа уравнений уменьшаются время решения их на ЭВМ и объем требуемой машинной памяти.

16.2.2. Распределение температуры вблизи кабеля

Задача, в которой рассматривался кабель, помещенный в теп­ лопроводящую среду, была решена в гл. 8 (фиг. 8.6), причем ис­ ходная область разбивалась на 96 элементов с 65 узлами. Эта «область была рассмотрена с применением разбиения на четыре

yt

«Фиг. 16.5. Четырехугольные элементы, используемые при исследовании переноса тепла вокруг кабеля.

квадратичных четырехугольных элемента, а также разбиения на два кубичных четырехугольных элемента. Эти элементы были рас­ положены (фиг. 16.5) так, чтобы кабель находился в узле.

Сравнение результатов, полученных для трех различных моде­ лей, проведено на фиг. 16.6, где показано распределение темпера­ туры вдоль линии у = 4 см. Модель, использующая симплекс-эле­ менты, и модель из четырех четырехугольных элементов дают по существу совпадающие результаты, тогда как двухэлементной мо­ дели соответствуют меньшие значения температуры в точке рас­

модели, использующей симплекс-элементы, четырехэлементной и двухэлементной моделей. Как видно, различие в приведенных значениях температуры кабеля несущественно, а значения темпе­ ратуры в точках верхней границы области, вычисленные по трем моделям, отличаются не более чем на 0,5 °С.

Модель из четырех четырехугольных элементов, видимо, точнее двухэлементной, так как на большой части первого элемента в двухэлементной модели температура получается постоянной

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

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

Фиг. 16.6. Распределение температуры в горизонтальной плоскости, содержащей точку расположения кабеля.

------------ симплекс-элем енты ;---------- четыре четырехугольных элемента; — — два четырех­ угольных элемента.

четырехэлементной модели необходимо решить систему из 23 урав­ нений, для двухэлементной модели это число равно 19 вместо 65 уравнений для модели, сконструированной с использованием симплекс-элементов.

16.23. Концентрация напряжений в зоне выточки

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

в зоне выточки получаются сомнительные значения а. Во-первых,

значения ахх и ауу в вершине прямого угла не равны нулю, как это должно быть на самом деле, а составляют величины порядка 2000 Н/см2. Узловому значению верхнего из узлов, общих для третьего и четвертого элементов, соответствуют два числа, отли­ чающихся между собой примерно на 15 000 Н/см2. Все это пока-

Фиг. 16.7. Некоторые значения СГ| для модели из четырех элементов в виде четырехугольников.

зывает, что данное разбиение области на элементы еще не доста­ точно мелкое. Величины напряжений вдоль криволинейной грани­ цы выточки очень сильно изменяются, и для моделирования этих значений двух уравнений явно недостаточно. Наибольшее значе­ ние CTI получается порядка 67 000 Н/см2. Это соответствует коэф­ фициенту концентрации напряжений 1,52, что отличается от при­ нятого значения 1,42 на 7%.

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

16.3. Криволинейные границы

Размещение узлов при задании формы элемента и определение узловых значений {/(е>} не составляют большого труда, когда гра­ ницы элемента прямолинейные. Однако наличие напряжений или конвективного теплообмена на криволинейных границах резко