________УЧЕНЫЕ ЗАПИСКИ Ц А Г И
Том XVII 1986
М 2
УДК 629.7.015.4
РАСЧЕТ РАСПРЕДЕЛЕННЫХ НАГРУЗОК НА ОСНОВЕ АНАЛИЗА СЛУЧАЙНОГО СТАТИКО-ДИНАМИЧЕСКОГО НАПРЯЖЕННОГО СОСТОЯНИЯ КОНСТРУКЦИЙ.
Ч. 1 —РАСЧЕТ ЭКСТРЕМАЛЬНЫХ НАПРЯЖЕННЫХ
СОСТОЯНИЙ
В. Д. Ильичев, В. А. Клименко
Предложен вероятностный метод определения экстремального по интенсивности сложного напряженного состояния при случайном статикодинамическом нагружении летательных аппаратов. Метод позволяет определять экстремальные напряженные состояния и расчетные по прочности зоны пространственных конструкций. Эти данные в дальнейшем могут использоваться при определении общих и местных распределенных нагрузок для повторно-статических испытаний на прочность и для проектировочных расчетов.
В авиационной практике случайный характер нагружения конструкций, когда параметры внешнего воздействия на нее в условиях эксплуатации суть случайные величины или (и) случайные процессы, является не исключением, а правилом. При некоторых реализациях эксплуатационных условий в различных зонах конструкции возникают наиболее опасные для прочности конструкции напряженные состояния (НС), мгновенные при статико-динамическом случайном воздействии на летательный аппарат (ЛА). Определение таких зон, состояний и обобщенных характеристик внешних воздействий составляет первую часть рассматриваемой здесь задачи.
Распределенную статическую нагрузку, вызывающую опасное напряженное состояние в какой-либо зоне конструкции, назовем, с некоторой долей условности, расчетной нагрузкой, если это состояние является предельно допустимым по интенсивности. Определение расчетных общих и местных нагрузок составляет вторую часть задачи, имеющей очевидное практическое значение для проектирования и испытаний на выносливость и прочность ЛА, его агрегатов и конструкционных элементов. Будем искать решение комплексной статико-дина-мической задачи на основе совместного использования метода конечных элементов (МКЭ) и метода заданных форм возможных перемещений (МЗФ) (см. [1—3]).
В качестве обобщенных параметров внешнего статического и динамического- воздействия и НС- возьмем обобщенные перемещения МЗФ—Z (£). Тогда вектор узловых перемещений конечно-элементной схемы ЛА или его подконструкций С? выражается как
<Ш = ч^(0. (1.1)
где Ч» — матрица узловых значений заданных форм МЗФ.
С учетом (1.1) тензор статико-динамического НС в любом конечном элементе примет вид:
9, V) = н} <?, ц) = н} (0 = н, z (() = 0;. z (О,
где (?/ (О Я У) — подвектор узловых перемещений у-го элемента, Н] — матрица напряжений элемента, Ч1^ — матрица заданных форм перемещений в узлах элемента, 0у = Н^} — матрица модальных напряжений элемента.
Будем считать вектор 2 (/) нормальной случайной стационарной векторной функцией, полагая, что таковы же внешние нагрузки на Л А, который рассматривается как линейная система. Вектор Z (£) будет складываться из случайной статической и случайной динамической независимых составляющих, т. е.
Z(o=Zo +ДZ(0,
где ДZ (0 считаем центрированным случайным вектором.
Тогда его математическое ожидание, выражаемое оператором е[ ], равно нулю
•дг=е[дг(<)] = о
и соответственно
вг = е[г(*)]-*^0] = Ч.
В силу независимости Zo и ДZ(^) ковариационная матрица вектора Z (() равна (см. [4]) сумме ковариационных матриц слагаемых
= 8 I (Я (*) — ег) (Z{£) - г2у] = + КА2,
где .
^о = е[(^-его)(Л-г2У],
/Сдг = 8 [дг(0 ДZт(Л] = const.
В качестве меры интенсивности напряженного состояния /-го элемента будем использовать или линейный критерий общего вида
Зэк Ву.(О = ^7^(О = ^0;2(О, (1.2)
где £;- — вектор коэффициентов линейной функции, или нелинейный критерий Мизеса
«ЭКВ;(0 = К'(*)М<7,(0Р=1 ^(О6;Л10,^(ОР , (1.3)
где Ж — постоянная матрица (см. [5]):
м
2 -1 -1 0 0 0
— 1 2 —1 0 0 0
-1 -1 2 0 0 0
0 0 0 6 0 0
0 0 0 0 6 0
0 0 0 0 0 6
\
Заметим, что если
z* 6} лм, z0 ы щ 6} м% дz(o,
то линейный критерий вида (1.2) можно получить линеаризацией нелинейного (1.3). Тогда
•I вт м «.к.Д*) =------—-----т
(ет2 вТ M0j 2
Для линейных критериев вида (1.2) математическое ожидание ^экв/ и дисперсия йэкв} процесса оэкв/(£) связаны с и Кг простыми соотношениями
^экв у == ^^экв ] == в/ Ь}.
Для нелинейного критерия (1.3) более удобны следующие оценки этих вероятностных характеристик (см. [5]):
^экв ~ ®/ МЬ] ^экв У> 1 ^
3*в У = Яр [АЪ, Кг »/ лт] > с/экву, )
где матрица А определяется из соотношения АТА=М, а ] — след матрицы [ ].
Можно показать, что в этом случае
2 , . ~2 ~1 ________________
<?экв / . ^*экв у ^экв / “1“ экв / ' ^экв ут
где -/Иэкв з — второй момент 0ЭКВ у (О (см. [5]) .
Эквивалентные напряжения сГэквН^), индуцированные случайными обобщенными перемещениями Z(^}, так же, как и последние, являются случайными процессами.
Легко себе представить неимоверную трудоемкость прямого определения экстремальных значений сГэквЛ^) во всех конечных элементах расчетной схемы, например, методом статистических испытаний. Значительно сократить эту трудоемкость позволяет вероятностно-экстремальная постановка данной задачи.
Не учитывая несущественную в данном случае хронологическую последовательность мгновенных воздействий и состояний, можно рассмотреть все множество реализаций случайного нормального «-мерного вектора Z, вероятность появления которых в области его определения равна заданной величине Ро, и определить экстремальные значения эквивалентного напряжения для этих реализаций. Полагая будем считать, что появление маловероятных реализаций Z вне его области определения невозможно, т. е. примем для Z усеченный нормальный закон распределения вероятности.
Так как гиперповерхностями постоянных значений плотности распределения вероятности в пространстве Z являются гиперэллипсоиды [6], то естественна следующая постановка задачи.
Найти максимумы функции аЭКвз, когда Z принадлежит области, ограниченной гиперэллипсоидом рассеивания
(Z— *z)T Kzl(Z-Sz)<?2 (1.5)
при \Kz'\^= 0.
Чтобы определить величину параметра р, воспользуемся выражением многомерного нормального распределения вероятности для Z
Р= (2*f Т || Kz II ^ J • • • J exp [- -i- (Z— ezY K7l (Z-гг)] dZx... dZnr
где V — объем гиперэллипсоида (1.5).
После преобразования координат — поворота, сжатия и перехода к обобщенным сферическим координатам — можно получить для вероятности нахождения вектора Z внутри гиперэллипсоида (1.5) следующее выражение:
Р=(2тс) 2 jexp^-/-2)?(r)dT/-, (1.6)
о
где q(r)dr — сферически симметричный элемент объема в п-мерном пространстве.
Соотношение (1.6) позволяет для каждого заданного значения вероятности Р0 при известной размерности п вектора Z вычислять значение параметра р и наоборот.
Таким образом, задача определения экстремальных напряженных состояний сводится к следующей задаче математического программирования:
шах аэкв AZ), |
1 0-7> (Z-*zYKZ4Z-*z)<?(Pa,n). J
Так как максимумы выпуклых линейной (1.2) и квадратичной (1.3) функций, определенных на замкнутом ограниченном множестве, достигаются на границе этого множества (см. [7]), то поиск максимума достаточно проводить только на поверхности гиперэллипсоида. Поэтому можно, заменяя неравенство на равенство, перейти от задачи математического программирования (1.7) к задаче определения условного максимума:
шах оэкв у (Z), I
(1.8)
(Z-BzyKz'(Z-tz) = P2(Po, я). I
Для линейного критерия интенсивности напряженного состояния
(1.2) составим функцию Лагранжа
Z=/?/Z+X[(Z — ггУ К? (Z- •*) - Р2 (Р0, я)], где R] — L]bj,
и запишем необходимое условие экстремума
grad L = R] + 2Х (Z - ez)T KJ1 = 0, откуда определяется стационарное значение Z:
(1.9)
Подставляя (1.9) в равенство (1.8), получим уравнение относительно X, решениями которого будут
Эти значения % дадут согласно (1.9) два значения Zexty, соответствующие максимальному и минимальному значениям линейной формы
(1.2) в области (1.5)
Zexti = Sz±?(P0,n)-----^Ц-. (1.10)
\R)KzRj\ 2
Умножая (1.10) слева на /?}, согласно (1.2), получим выражение экстремальной интенсивности НС в виде
®ext j R. j Zext j == ^эвв j dt P^3tKB j •
Для квадратичного критерия (1.3) задачу на условный экстремум с учетом того, что положение тахаЭквз в пространстве параметров Z совпадает с положением шах а|кв 3, сформулируем в виде:
max За кв / (Z) =шах [Z10} Mbj Z\, (Z — ez)т Kz\Z — sz) = р2 (Р0, п).
Функция Лагранжа и необходимое условие экстремума примут соответственно вид:
L = ZT 0} Mbj Z + X [(Z- ezy Кг1 (Z - *z) - P2 (P„, n)\
и
grad L = 2 [0} Af0;Z + X/CJ1 (Z- e2)] = 0,
откуда
Zext/=[^-/Cz0jjWey + £T1eir. (1.11)
Подставляя (1.11) в (1.8), получим уравнение относительно X {Kz 0} MQj -Ь - £} /CJ1 {[у - Kz 0/T Лв, + f]- £} =
= P2(*V«)- (1.12)
Чтобы избежать численного решения этого довольно громоздкого относительно X уравнения, перейдем от параметров Z— обобщенного воздействия на конструкцию в целом — к параметрам Я,— локального воздействия на j-й конечный элемент, / = 1, % ..., 7V. является вектором статистически независимых компонент НС.
Прежде всего заметим, что для большинства типов элементов из библиотеки МКЭ тензоры напряженного состояния имеют нулевые компоненты. Поэтому для записи эквивалентных напряжений по Мизесу вместо полной матрицы М достаточно использовать ее подматрицу М, состоящую из строк и столбцов, соответствующих ненулевым компонентам тензора напряженного состояния элементов данного типа.
Так как матрица М имеет дефект, равный 1, т. е. ее ранг на 1 меньше размерности, то матрица М может иметь дефект, равный 1 или 0. Представим М в виде, аналогичном (1.4)
М = АГА,
sT
где А — прямоугольная матрица, число строк которой т равно рангу, а число столбцов--размерности матрицы М.
Ранг матриц А и М меняется от 1 (для стержней) до 5 (для пространственных элементов). Теперь выражение для квадрата интенсивности напряженного состояния представим в виде
зэкв / - zI 6} ЛТ а в, z=гт в) в} г, (1.13)
где через В} обозначено произведение АЬ}.
Если ввести вектор X; размерности т}, связанный с Z соотношением
х) = у]в]г> (1.14)
где У} — ортогональная квадратная матрица, то выражение (1.13) примет следующий вид:
ОэКВ/ = ZтJB}Я;Z=Z^ff;т V5Vr^Л^Z= Х)Х}. (1.15)
Вектор параметров локального воздействия на у'-й конечный элемент X) так же, как Z, является случайным нормальным вектором (см. [4]), и согласно (1.14) его математическое ожидание гх} и ковариационная матрица Кх; соответственно равны
гх} = V, В} Кх, = V, В} Кг Щ V],
или
КХ]=У]С]У), (1-16)
где
С} = В}-Кг-В).
Возьмем в качестве строк ортогональной матрицы V) собственные векторы матрицы С;-, удовлетворяющие уравнению
С}У1=У}81, (1.17)
тогда — диагональная матрица собственных чисел С} — равна согласно (1.16) и (1.17) Кх/, так как’5^= У} Су V] = КХ]-.
Теперь выражение гиперэллипсоида рассеивания вектора Л} примет вид:
(X; - е,/ 5.-> (ХГ- *х,) - £ (Р0, щ}).
Переход от '2 — вектора обобщенных параметров внешнего воздействия на всю конструкцию — к Х^~ вектору локального воздействия на каждый конечный элемент — позволяет сформулировать задачу поиска экстремальных воздействий следующим образом. Найти максимум функции /{Х}), где
/(*Д=о2экв (1.18)
при условии
' ? (X}) = {X} - е,;)т (Я) - %х}) - Р* {Ро, т,) = 0; (1.19)
здесь гп; — число компонент вектора Х}.
Необходимое условие экстремума функции Лагранжа примет в этом случае вид:
'ЗС,+\ХЩ-*ХУ8Т1~0, где — соответствующий множитель Лагранжа.
Отсюда можно получить два эквивалентных выражения для определения вектора Xextj, соответствующего экстремальному значению оэы]:
= гх} ~(Е + KSr'T' *xJ, (1.20а)
X'«j = K(Sj + KET1*xj- (1-206)
Подставляя (1.20а) в условие (1.19), получим
e^{Sj + lxErTSjiSj + KErlKj^(Po, rnj) ] или, с учетом диагональности матриц Sj и (Sy4-)s>£)j т>
mj s- с.
V- 'V--------?2ЛР0, mj) = 0, - (1.21)
где гх1 — элементы вектора гх}, a S( — элементы диагональной матрицы Sr
Подстановка в выражение (1.18) значения Xextj из (1.206) дает для f(Xj) следующее выражение: (
f(Xj) = X* ш'х1 (Sj + X, Е)-2 г,,. (1.22)
ИЛИ' ■
* -2 2
Уравнение (1.21) дает не более 2/га,- действительных значений Х^. — множителей, соответствующих локальным, экстремумам функции f{Xj) на гиперповерхности (1.19). Можно показать, чтд вс§ действительные корни (1.21), кроме, быть может, одного, отрицательные, причем ,
шах | ХГ(?! > шах 5,, г=1, .. ., т},
q = 1, . . ., 2
Докажем, что искомому глобальному максимуму соответствует максимальный по абсолютной величине отрицательный корень —~кхе. Для этого вычислим Д/—разность между значениями функции (1.22) при \хе и У.хд, где \^хе\^\Кя\> Ч— 1. • • •> 2т j, с учетом того, что Х^ и кхд удовлетворяют уравнению (1.2,1). В результате получим
д/=/(^г) - /(*„) = - (Ке - Кд)3 2 № + x^'(S. + x^)2 •
Так как выражение под знаком суммы всегда положительно, то Д/>0, когда | Х„ | > | | и Ххе < 0. Поэтому для определения
шахоэкв^ требуется найти лишь один максимальный по абсолютной величине корень уравнения (1.21).
Алгебраическое уравнение (1.19), порядок которого может меняться от 2 до 10, решается численно методом Ньютона. Поиск максимального по модулю корня ведется на интервале (—со, —.5>Шах); ЗДеСЬ *Smax = max Sit i — 1, . . ., ttij.
В качестве начального значения удобно взять такое Х^0), при котором заведомо ? (Х(*0)) >-0. Например, можно взять
Xj: = — -Smax I 1 + у
?x Si )
или
х(0) __ тіп X і
Так как функция на интервале (—с», — 5тах) монотонная, то итерационная процедура всегда будет сходиться, причем достаточно быстро (см. [8]).
Обратим внимание на то, что число корней уравнения (1.12) относительно X и уравнения (1.21) относительно Х„ вообще говоря, различно. Можно показать, что в силу соотношения (1.14) между ^и и равенства (1.15) абсолютные максимумы на границах гиперэллипсоидов в пространствах параметров Z и Я) совпадают и достигают при одинаковых значениях X и *.х, если
Это следует из того, что, умножая (1.11) слева на У)В) с учетом (1.14) и (1.16), получим (1.20). Однако прочие условные экстремумы °эКВ/С£) могут и не являться условными экстремумами 3экв/(^0) в пространстве параметров Я}, потому что лишь часть точек границы гиперэллипсоида (1.8) образует множество точек границы гиперэллипсоида (1.19), а остальные являются внутренними точками последнего.
Подставляя полученное значение максимального Хх в выражения (1.11), (1.20) И (1.22), получим соответственно Zmax> Яшах И таХовкв/ для каждого элемента. Зная наибольшее значение эквивалентного напряжения каждого элемента и вызывающие его обобщенные перемещения Д можно перейти ко второй части задачи-— к определению распределенных максимальных нагрузок на ЛА и его подконструкции. Этот вопрос будет рассмотрен во второй части публикуемой статьи.
1. Ильичев В. Д. Матричные методы синтеза динамических и упругих характеристик линейных неконсервативных конструкций. — Ученые записки ЦАГИ, 1975, т. 6, № 2.
2. С е р г е е в а А. А. Определение упругих перемещений летательного аппарата и оценка его напряженно-деформированного состояния при многоуровневом комплексном прочностном расчете (внешняя задача).— Труды ЦАГИ, 1982, вып. 2155.
3. Клименко В. А. Определение напряженно-деформированного состояния агрегата при многоуровневом комплексном прочностном расчете летательного аппарата (внутренняя задача). — Труды ЦАГИ, 1982, вып. 2155.
4. Себер Дж. Линейный регрессионный анализ. — М.: Мир, 1980.
5. И л ь и ч е в В. Д. Применение спектрального суммирования усталостных повреждений при сложно-напряженном состоянии конструкции. — Ученые записки ЦАГИ, 1979, т. 10, № 2.
6. Б р а й с о н А., X о ю - Ш и. Прикладная теория оптимального управления. — М.: Мир, 1972.
7. Рокафеллар Р. Выпуклый анализ. — М.: Мир, 1973.
8. Форсайт Дж., Малькольм М., М о у л е р К. Машинные методы математических вычислений.—М.: Мир, 1980.
р (Р, п) = рх(Р, яг;) = соші.
(1.23)
ЛИТЕРАТУРА
Рукопись поступила 30;У1 1983 г.