Оптимизация формы элементов авиационных конструкций с концентраторами напряжений Текст научной статьи по специальности «Физика»

УЧЕНЫЕ ЗАПИСКИ ЦАГ И Том XVI 1985

№ 6

УДК 539.4.013.3 : 624.07

ОПТИМИЗАЦИЯ ФОРМЫ ЭЛЕМЕНТОВ АВИАЦИОННЫХ КОНСТРУКЦИЙ С КОНЦЕНТРАТОРАМИ НАПРЯЖЕНИЙ

В. И. Гришин, Ф. В. Рыбаков

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

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

В процедуре оптимизации конструкции, приведенной в работе [3], предполагается, что оптимальной формой контура концентратора является такая форма, при которой распределение тангенциальных напряжений о"( по контуру близко к равномерному. Величина напряжений определяется с помощью МКЭ. Градиент целевой функции в предложенном алгоритме не используется.

В процедуре оптимизации, описанной в работе [4] и использующей метод Давидона—Флетчера—Пауэлла, градиент целевой функции вычисляется из уравнений, полученных путем дифференцирования основных соотношений МКЭ. Как показали результаты расчета, указанная методика работает надежно, но затраты времени на вычисления в процессе получения итерационного решения значительны.

В настоящей работе развит метод оптимизации, изложенный в [4], на основе эффективного применения метода подконструкций. Метод реализован в специализированном комплексе программ СКП— ФИТИНГ [5]. Приводятся алгоритм выбора оптимального шага на каждой итерации и решения задач об оптимальной форме контура для пластины с отверстием, галтельного перехода и проушины.

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

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

При численной реализации метода оптимизации напряжения рассматриваются в отдельных контрольных точках г=1, 2,...,/, граница определяется вектором конструктивных параметров t={t ь ?2,

Тогда целевая функция имеет вид

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

Математическую формулировку задачи кратко можно записать следующим образом:

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

Рассмотрим новую целевую функцию

где С[£т (/)] — штрафная функция для ограничения gm (і), г3 — пара-

Ф* (В) = гпах

/?

£т(*)<0, т = 1,2, . . . , М.

(1)

т

метр штрафа на я-й итерации по штрафу. При переходе от 5-й итера-ции по штрафу к я+ І-й параметр г8 изменяется по формуле г^+і = — .

В качестве штрафной функции б возьмем комбинированную функцию, предложенную Кавли [8]:

1 Для £т(*)< —

О [ёт (*)] =

£т V)

2вд + 5”(*)- для &,,(*)>■

•?

где параметр г8 = г8/б определяет переход от одного типа штрафной функции к другому.

При изменении вида штрафной функции и на границе допустимой области штраф равен соответственно б и 26. Такая функция непрерывна и имеет непрерывную производную в точке перехода.

Для минимизации функции г{() без ограничений воспользуемся методом наискорейшего спуска [7], согласно которому на к-й итерации значение оптимизируемых параметров вычисляется по следующей формуле

tk+\ = 2* + X*

где /гА = — — антиградиент функции г (^) в точке tk, параметр

шага ХА определяется из условия

г(** + М*) = п11п2(*, + ХЛ1к). (3)

Х>0

Для определения Хк приближенно представим функции

г(** + М*) и —~^-ХЙй) = IVг (<» + >•**), Ьк\

тейлоровскими разложениями по степеням X с использованием соответственно четырех и трех членов разложения:

Н-ХЛА) як г(/й)-{-ХА +-~5+-^-С, (4)

(V* (** + М*), А*)« А + \В + С, (5)

где

А-(у*(**), Нк). (6)

Введем обозначения

6 (X, *А) = 2 (*, + Мк) — 2 (**),

М*. ^) = (уг(^ + ХА4), Л*).

Из системы уравнений (4), (5) получаем

В = -~(3в (X, **) - Хбх (X, 4) - 2ХА), (7)

Л3

С = А(Х0х(Х, *А)-29(Х, г,) + ХА). (8)

Фиксируя некоторое значение Я=Х, по формулам (6—8) определяем А, В, С. Подставляя А, В, С в (4), находим значение параметра

Лк, приближенно удовлетворяющее условию (3).

Если при этом 2 (£* + X* ЛА) > 2 (**), то принимаем:

0<р<1, (9)

и повторяем вычисление параметра шага, удовлетворяющего условию (3). Проведенные расчеты показали, Что выбор нового значения X по формуле (9), а не по формуле 1(3, как это предлагают Флетчер и Пауэлл [7], может существенно ускорить процесс поиска оптимального значения А,к.

Описанная выше процедура поиска оптимального Як не срабатывает, если исходная точка в пространстве параметров находится слишком далеко от положения минимума. В этом случае можно аппроксимировать целевую функцию кубической параболой на отрезке А,е[Я, 2Я]. Чтобы определить значение Як, лежащее в интервале [Я, 2Д можно использовать формулы (4) — (8), в которых допускается изменение параметра X в промежутке [О, Я]. Для этого в них нужно заменить А, на I, X — на I, tl^ — па tк+Xhк и, вычислив оптимальное значение /, положить Хк = Х+1. При повторной неудаче нужно аппроксимировать функцию на отрезке А,е[2Я, ЗХ] и т. д.

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

вектор однозначно определяет производную от целевой функции г

по параметру '1,. Для одного конечного элемента вектор -др- определяем, дифференцируя известные формулы МКЭ [9]:

дяе п де,

д(у д^

где И — матрица упругости, ее — вектор деформаций.

Далее, дифференцируя соотношение, связывающее величину узловых перемещений с деформациями в произвольной точке элемента = ВЪе, получим

Л* + в ля

дЬ] ді]- дЬ]

где В — матрица деформаций, Ье — вектор узловых перемещений элемента.

дВ

Матрицу определяем так же, как в работе [4]. Для опре-

; „ дъ

деления производных от вектора узловых перемещении запишем основную систему МКЭ для подконструкции [10]

Ки К„

Кщ кппл

(10)

где К— матрица жесткости, 3 — вектор узловых перемещений, Я вектор внешних сил.

Индекс п указывает на принадлежность перемещений и нагрузок к внешним узлам, а индекс г — к внутренним.

Исключая из соотношения (10) перемещения внутренних узлов (прямой ход метода Гаусса для подконструкдий), получим систему уравнений равновесия для «внешних» узлов подконструкции

Имея выражение (11) для каждой подконструкции, объединяем их суммированием и получаем общую систему уравнений равновесия для всей конструкции

Решая систему (12), определяем 3, а зная эти перемещения и используя верхнюю строку соотношений (10), получаем значения перемещений всех внутренних узлов подконструкции

что соответствует обратному ходу метода Гаусса для системы (13). Дифференцируя систему (10), получаем

для подконструкции, которой принадлежит варьируемый контур, и

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

Таким образом, для нахождения вектора производных от узловых перемещений нужно по соотношениям (15) вычислить правую часть системы^ (14) для подконструкции с концентратором и приравнять нулю правые части остальных систем. После этого нужно повторить решение систем уравнений МКЭ, преобразуя только полученные правые части, матрицы же систем остаются треугольными после предыдущего решения.

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

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

(П)

/С8 = /?.

(12)

/СпЗг = /?г-/С1п8„

(13)

(14)

где

(15)

форма отверстия,

3 п

-к- . В качестве конструктивного параметра і

Рис. Г

выбрана ордината точки отверстия с абсциссой х=1,5 см. Форма отверстия аппроксимировалась кривой второго порядка, уравнение которой имеет вид

ау2 + Ьху + сх2 4- <Лх 4- еу + /= О

с фиксированными точками А и В (рис. 1) на концах контура, с горизонтальной касательной в точке А и вертикальной касательной в точке В.

Найденное в результате расчета оптимальное значение параметра ^ОПТ = 1» 73 см соответствует эллипсу, что согласуется с известным аналитическим решением [И]. Распределение максимального главного напряжения на контуре АВ приблизительно постоянно, коэффициент концентрации Кэп =1,72 (по отношению к номинальному напряжению по сечению брутто) составляет 43% от исходного коэффициента концентрации напряжений /СРомб = 3,96.

Значение коэффициента концентрации напряжений на контуре эллиптического отверстия с полуосями а, в (рис. 1,6) в бесконечной плоскости, определяемое по выражению [11]

Ь 2 ау 2 при — = , —— = -у-, равно 1,67 и отличается от результата

расчета на 3%. Форма отверстия и распределение коэффициента концентрации напряжений вдоль контура отверстия по итерациям (О, 1 . . .) приводятся на рис. 1, б.

В данном примере конечно-элементная сеть состояла из трех под-конструкций и имела 106 узлов, причем 62 из них принадлежали перестраиваемой подконструкции, что составляет 58% общего количества узлов. За счет того, что при повторном расчете перестраивалась и пересчитывалась лишь одна подконструкция, затраты машинного времени

на вычисление оптимального значения целевой функции в данном примере сокращен приблизительно на 45% по сравнению со временем, которое потребовалось бы при расчете без использования метода под-конструкций. Общий расход времени сокращен приблизительно в 1,5 раза. С уменьшением доли степеней свободы, приходящихся на перестраиваемую подконструкцию от общего количества степенной свободы, эффект применения метода подконструкций в задаче оптимизации будет увеличиваться.

Аналогичные расчеты произведены для той же самой пластины,

но при разном соотношении нагрузок Полученные оптимальные

ау

формы отверстия, соответствующие им значения параметра 10т/Ь, а также графики коэффициента концентрации на оптимальном Копт и эллиптическом Кэл отверстиях В зависимости ОТ отношения Ох/Оу изображены на рис. 2.

Рис. 2

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

Лт]г + 52 + 2В-г) = 0, где ^ и £ —локальные координаты, показанные на рис. 3,а. Координа-

Рис. 3

ты точки стыка и и ^ полностью определяют форму перехода и являются параметрами оптимизации.

Коэффициенты уравнения кривой А я В выражаются через конструктивные параметры следующим образом:

л_ Т-и ( Т-и о и \

4 — I 4-*а *2

я - <7-_ « (.А-- ^=А.),

на параметры и и /2 накладываются ограничения ^>0, причем

ограничение t^>■0 активно.

В расчете конечно-элементная сеть состояла из 33 квадратичных изопараметрических элементов, что соответствует 256 степеням свободы. Для оптимизации потребовалось пять итераций и 20 минут машинного времени на ЭВМ БЭСМ-6. Результаты расчета для различных длин перехода приведены на рис. 3, б, где 01 — максимальное главное напряжение.

Оптимизация формы внешнего контура проушины. Геометрия рассматриваемой проушины показана на рис. 4, а. Линейные размеры указаны в метрах. Задача решалась для двух случаев нагружения; при равномерном давлении а на половину внутренней поверхности (рис. А, а) и при давлении, заданном по закону косинуса (рис. 4,6).

Форма варьируемого контура задается полиномом пятой степени

у = а0 + ах х + а2 х2 + а3 х3 + Д4 х* + аь х5,

Рис. 4

где коэффициенты а0, аь а2, а3, а4, а5 определяются через ординаты ^ и и двух внутренних точек контура. Эти ординаты и являются оптимизационными параметрами. Оптимизируемый контур заключен между фиксированными точками Л и В, касательные к контуру в которых горизонтальны. На параметры ^ и 4 накладываются ограничения ^<15 см, £г<15 см.

Исходные формы проушины и формы, полученные в процессе оптимизации, показаны на рис. 4, а и 4,6. Там же дано распределение коэффициента концентрации напряжений вдоль границы для исходной /Сисх и оптимальной /Сопт форм.

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

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

ЛИТЕРАТУРА

1. D и г а 11 i A. J., Rajainh К. Lighter and Stronger. — Exp. Mech., 1980, vol. 20, N11,

2. Steinchen W. P. Experimental and computeraided investigations and optimization of stress relieving notches. — J. Strain -Analysis, 1978, vol. 13, N 3.

3. S с h n а с k E. An optimization procedure for Stress concentrations by the finite element technique. — Int. J. for Numerical Methods in Engineering, 1979, vol. 14, N 1.

4. Fr а п с a v i 11 a A., R a fti а к г i s h п а п С. V., Z i е п к i е-w i с z О. С. Optimization of shape to minimize stress concentration, — J. Strain Analysis, 1975, vol. 10, N 2.

5. Гришин В. И., Донченко В. Ю„ Барышников В. И., Тихонов Ю. В. Применение метода конечных элементов к исследованию местной прочности элементов авиационных конструкций. — Ученые записки ЦАГИ, 1983, т. XIV, № 1.

6. Евтушенко Ю. Г. Методы решения экстремальных задач и их применение в системах оптимизации.—М.: Наука, 1982.

7. По лак Э. Численные методы оптимизации.—М.: Мир, 1974.

8. К a v 1 i е D. Optimum design of statically indeterminate structures. Thesis submitted to the University of California, Berkelley, for the degree of Ph. D„ 1971.

9. Зенкевич О. Метод конечных элементов в технике. — М.: Мир,

1975.

10. Гал’ки'йа Н. С,.; Гр й iti'H № ,В. I Иц> Д;С1 н Ч й йк о В. Ю,. Иссле-

дованиб здпряженно-деформированного состояния элементов авиационных, конструкций и их соединений. Труды, ЦАГИ. в^п. 2(312, 1979.

11. Петерсон Р. Е.‘ Коэффициенты' Концентрации напряжений.— М.: Мир, 1977.

Рукопись поступила 3/VI 1984