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

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

№ 3

УДК 629,735.33.015.4.025.1

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

ЭЛЕРОНОВ

С. В. Селюгин,

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

Задача расчета рациональных толщин силового набора конструкции крыла с заданными геометрическими параметрами при учете различных условий решалась рядом автором (см., например, обзор работ в [1]). В [1] рассмотрены также работы Г. В. Украинцева и В. М. Фролова, где предлагается варьировать и относительной толщиной профиля. В настоящей работе рассматривается задача расчета рациональных толщин силового набора и относительной толщины профиля конструкции стреловидного крыла с целью получения минимума ее веса при удовлетворении условий статической прочности и заданной эффективности элеронов. Традиционный метод решения задач оптимизации конструкций заключается в сведении этих задач к задачам нелинейного программирования. Размерность последних задач при этом является достаточно высокой (десятки, сотни варьируемых параметров). Следует отметить, что решение возникающих при оптимизации конструкций задач нелинейного программирования большой размерности „в лоб“ известными методами практически весьма затруднительно вследствие достаточно большого времени счета

на ЭВМ. Описываемая в работе декомпозиционная методика оптимизации ориентирована на решение соответствующих задач большой размерности.

1. Постановка задачи. Основные предположения и исходные соотношения работы заимствованы из [1]. Рассматривается стреловидное крыло достаточно большого удлинения. Приняты следующие допущения:

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

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

— эффективность элеронов определяется в соответствии с линейной аэродинамической теорией, при этом крыло заменяется срединной поверхностью и поточные сечения крыла считаются недеформируемыми;

— идеализированное сечение кессона плоскостью, нормальной к оси жесткости, представляет собой прямоугольник, причем силовой материал распределен равномерно вдоль его сторон и симметрично относительно осей инерции [1].

Вдоль оси жесткости вводится ось г, направленная по размаху крыла; 2 — 0 — соответствует борту фюзеляжа. Варьируемыми функциями, зависящими от z) являются приведенная толщина обшивки 8пр (г), толщина обшивки _&обш (г), толщина стенки ост (г), относительная толщина профиля с (г). В работе из соображений вычислительного удобства вместо функции ост (г) используется площадь стенки 5СТ (г). Жесткость крыла на изгиб Е1 (г) (Я— модуль упругости, /—момент инерции, соответствующий изгибу) определяется по инженерной теории изгиба, жесткость крыла на кручение С/Кр (г) ((?— модуль сдвига, /кр~ момент инерции, соответствующий кручению) — по формуле Бредта [1]. Вес силовой конструкции крыла Окр определяется через интеграл, в который 5пр (г) и 5сх(г) входят линейно. Нормальные напряжения в верхних и нижних панелях а (г) также определяются по известной формуле (см. [1], стр. 104) через изгибающий момент, Е1(х) и высоту кессона.

Обусловленное упругостью конструкции крыла распределение углов атаки вдоль оси жесткости аупр (г) при отклонении элеронов на единичный угол удовлетворяет интегральному уравнению равновесия [1]

Ф (2) ^ 1аУпР ^ ~

г

*кр гкр

\ \ Л2 Кпр

г г

'упр

где Ф (г)—некоторая функция, зависящая от геометрических параметров крыла; /кр — длина консоли крыла вдоль оси жесткости; Л15 Л2 — линейные по первому аргументу операторы, определяющие соответственно погонные аэродинамические крутящий момент и нагрузку; х— угол стреловидности оси жесткости,

Считается, что эффективность элеронов 5 есть отношение производных моментов крена по углу отклонения рулевых поверхностей на упругом и жестком крыле, величина I определяется (с точностью до коэффициента, зависящего от геометрии крыла) через интеграл, в который полученная из (1.1) функция аупр (г) входит линейно. Значение I должно удовлетворять условию

где £0 —заданная величина эффективности элеронов.

На варьируемые функции накладываются некоторые локальные ограничения, например

где опро, оОбш0, §ст0 — соответственно минимальные значения приведенной толщины обшивки, толщины обшивки и толщины стенки, а последнее условие из (1.4) описывает диапазон изменения величины ообщ (г)/8пр (г), принятый в практике проектирования из соображений устойчивости панелей;

где одоп (2)-максимально допустимый уровень нормальных напряжений, вообще говоря, зависящий от толщины силового материала.

Должны быть удовлетворены также и некоторые интегральные ограничения, например, следующего вида [1]:

где п — натуральное число, Ь (z)— хорда крыла.

Рассматриваемая в работе задача заключается в минимизации веса силовой конструкции крыла GKp по опр (z), 80бш (z), Scr(z), с (г) при ограничениях вида (1.3)—(1.6). Таким образом, в постановочном плане в настоящей работе осуществляется дальнейшее развитие работ Г. В. Украинцева и В. М. Фролова, где предложено варьировать также и относительной толщиной профиля; более точно учитывается условие на управляемость самолета по крену.

Gf

Функция аупр (z) удовлетворяет граничным условиям «упр (°) = Vp (*кр) =

1кр

о

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

2. Сведение задачи оптимизации конструкции крыла к задаче нелинейного программирования. Следуя общепринятому подходу в теории оптимизации, сведем рассматриваемую задачу к задаче нелинейного программирования. Вдоль оси 2 вводятся £ + 1 упорядоченные в порядке возрастания точки 2(-, причем 2*.ц =

= /кр. Варьируемые функции заменяются набором их значений в этих точках. В точке 2*+1 указанные значения вследствие краевых условий (1.2) приняты равными нулю. Интегралы, входящие в выражение для веса силовой конструкции крыла, уравнение (1.1), условия (1.3) и (1.6), вычисляются по методу трапеций [2] и, таким образом, превращаются в заданные функции от значений варьируемых функций в точках гь / —1, ..., к и функции аупр (г). В предлагаемой методике функция аупр(2), наравне с 8пр С2), &обш (2), £ст (г), с (2), описывается алгебраическим полиномом от 2 с неизвестными коэффициентами а1г I = 1, . . р (р — степень полинома), удовлетворяющим граничным условиям (1.2) и полученному в результате преобразований условию (1.3). Ограничения (1.4), (1.5) и преобра* зованное в соответствии со сказанным выше уравнение равновесия (1:1) естественным образом рассматриваются лишь в точках г1У 1= 1, . . & и в дальнейшем называются локальными ограниче-

ниями.

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

к

ПИП 2 а**,.;

а <=1

к

2^Д*/Х°> 5=1, , /;

г=1

Ьц

й 1}

где л; — столбец, составленный из столбцов х1 (параметров конструкции в сечении 2г), г=1, ..., £; а — столбец, составленный из %, 1,. . ., р\ аь г — 1, . . ., £ — некоторые столбцы, е15, 5=1,

• • - , К Ьф у= 1, . . ., пи У= 1, . . ., тц— некоторые функции, г = 1,..., к, естественно возникающие при дискретизации задачи оптимизации конструкции крыла; т1у пь г= 1, . . . , / — некото-

рые натуральные числа, „т“ — знак транспонирования. Условия (2.2) соответствуют интегральным ограничениям вида (1.6), условия (2.3),

(2.4)-—локальным ограничениям.

3. Декомпозиция задачи. Описание алгоритма оптимизации.

Как видно из (2.1), (2,3), (2.4), при заданном а минимизация

в (2.1) по х при ограничениях (2.3), (2.4) [без рассмотрения условия (2.2)] эквивалентна минимизации а]х{ по х1 при этих же ограничениях для того же г, 1,..., k, т. е. минимизация веса в этом случае может проводиться независимо по сечениям ziy г == 1,..., к. Таким образом, параметры а и условия (2.2) в определенном смысле характеризуют влияние сечений друг на друга. В обзорной работе [3] подобные задачи предложено называть задачами со „связывающими14 (coupling) ограничениями вида (2.2) и „связывающими" переменными а. В [4] описан алгоритм для решения задач со „связывающими41 ограничениями. Алгоритмы решения задач со „связывающими" переменными в случае линейных ограничений специального вида также известны [5]. В настоящей работе предлагается методика, основанная на указанных подходах.

Выпишем функцию Лагранжа задачи

k ik

L {xt х, X, {х, 6) = 2 ajxt + 51 К 51 eis (xi) +

t = l 5-1 f=l

k mi k ni

tSI n d4 (*< a>+ 2 Z e<7 bD (xf a>' i=l / = 1 1-1 /=1

где Xs, 5 , /, plJt J= 1,..., mlt &t7, j = 3,. . . , nt, i~3,. . .,

k — множители Лагранжа, X, p., 8 — столбцы, составленные соответственно из этих множителей.

Вынося знак суммирования по i за скобку, имеем

к I

L(x, У., X, tb в) — £ №,Т+ 2 К #is (xi) -f

/=1 5=1

т.

+ 2^ dl1 4' 2 Ьа «)]•

/=1 >=1

Обозначим через L(i)(xiy а, X, р., б) выражение, заключенное в квадратные скобки в (3.1), Тогда (3.1) запишем в виде

* >

Цх, а, к, ц, 6) = £/. «!(.*., X, ц, «).

1 = 1

Согласно теореме Куна — Таккера [5J в искомой точке оптимума х*, а* выполнены условия (вообще говоря, необходимые, так

как выпуклость функций еи, 5=1,...,/, £/,, 7= 1, • • -, щ, йф у = 1,..., ягг, 1=1,..., А, не предполагается)

1х

1я

Хв>0, 5=1, . . ., /;

^2^

/=1

н-// >0, 7 — 1, • • •, т1У /= 1,. . ., А;

Н-<; ац

где нижние индексы х и а означают соответственно частное дифференцирование по хил. Искомые значения х*, а*, а также X, [1, б должны удовлетворять (3.2)—(3.7) и условию (2.3).

Как нетрудно видеть из структуры функции Лагранжа (3.1), указанная функция при фиксированных а, X, |х, 6 является суммой по г некоторых функций (х{, а, X, [а, 6), зависящих лишь от х1У г = 1,..., А. Следовательно, при фиксированных Х>0 и а удовлетворение условий (3.2), (3.6), (3.7), (2.3) эквивалентно решению

следующих подзадач для сечений ги г = 1, . . ., к:

I

т1п [а!л,+ 2Х5^ м] ’

*1 5=1

и (*,-, а) — 0, у=1, . . . , П{;

йц (*„ а)<0, / = тг

Действительно, необходимыми условиями оптимальности в подза-

дачах (3.8) будут условия, г=1, . . ., к:

Щ{*1> а, X, {х, 6) = О,

7 = 1,..., тг,

Рийу(Хц «) = 0, / = 1,. . ., ГП:,

совпадающие соответственно с (3.2), (3.6), (3.7), и условия (2.3) (нижний индекс х1 означает частное дифференцирование по х1з 1=1,..., £). Таким образом, при удовлетворении условий (3.2),

6—„Ученые записки ЦАГИ“ № 3

81

(3.6), (3.7), (2.3) имеет место декомпозиция (расчленение) решаемой задачи нелинейного программирования (2.1)—(2.4) на подзадачи (3.8) существенно меньшей размерности для сечений. Если рассматривать решение подзадач (3.8) как способ удовлетворения условий (3.2), (3.6), (3.7), (2.3) при заданных а и Х>0, то определение а и X должно производиться из условий (3.3)—(3.5), что и делается в работе. Далее описывается итерационная методика расчета оптимальных а, X и, следовательно, соответствующих им из (3.8) х, причем для определения X служат условия (3.4), (3.5), а для определения а — условие (3.3).

Считая а заданным, для вычисления X воспользуемся алгоритмом, аналогичным описанному в [4] (стр. 167), который основывается на теоремах Эверетта (там же, стр. 168—169).

Основной смысл этих теорем в принятых в рассматриваемой задаче (2.1)—(2.4) обозначениях заключается в следующем: если для некоторых заданных значений Х>0 из подзадач (3.8) получены их решения I,. . ., й, то значение левой части

ограничений (2.2) при х^хА можно уменьшить, увеличив соответствующие значения X и, обратно, увеличить, уменьшив X. В соответствии с этим для каждого значения 5, 5 —1, . . ., I можно посредством варьирования Х5 либо удовлетворить соответствующему условию (2.2) как равенству, либо сделать Х5 = 0 при (2.2), выполненном как строгое неравенство и, таким образом, удовлетворить условиям (3.4), (3.5). В работе считается, что в случае удовлетворения (естественно, с некоторой точностью) условий (3.5) для х~хА значения Х>0, для которых выполнены по построению алгоритма, и (3.2), (3.6), (3.7), (2.3) определены. Если же условия

(3.5) для некоторых 5 не удовлетворены, то для этих 5 осуществляется изменение Х_у по правилу:

^=>'5-Н

£

где (£*)_у = 2 еи О**)» ~ — некоторый положительный малый пара-

/=-1

метр, выбираемый с учетом условия X :^0. Описанный алгоритм (3.9) нахождения X при фиксированных а, отличных от оптимальных, рассматривается в работе как способ удовлетворения условий (3.4),

(3.5) и (3.2), (3.6), (3.7), (2.3). Подчеркнем, что для каждого значения X в алгоритме (3.9) необходимо решать подзадачи (3.8).

Последними условиями, которым в соответствии с теоремой Куна — Таккера необходимо удовлетворить для решения задачи

(2.1)—(2.4), являются условия (3.3). В предлагаемой методике удовлетворение этих условий осуществляется посредством варьирования параметров а. Отметим следующее полезное свойство: если X вычислены в соответствии с алгоритмом (3.9), а из подзадач (3.8) для этих значений X определены хА, ц, 6, тогда при условии справедливости обобщенной теоремы о неявной функции [6] для выполненных в этом случае условий (3.2), (3.4)—(3.7), (2.3), величины хА, X, |л, 0 будут функциями от а, и, таким образом, функция Лагранжа (3.1) также будет функцией лишь от а вида £[;сА(а), а, X (а), ^ (о), 0 (о)]. Изменяя а в направлении антиградиента функции Лагранжа I [хд (а), а, X (а), ^ (а), б (а)] по а, фактически уменьшаем значение целевой функции, так как все остальные члены в (3.1) равны нулю вследствие удовлетворения условий

(3.4)—(3.7), (2.3). Для вычисления указанного антиградиента в явном виде оказывается полезным следующее свойство функции Лагранжа (3.1).

Пусть выполнены условия (3.2), (3.4)—(3.7), (2.3) и предположения обобщенной теоремы о неявной функции [6]. Тогда величины хА, X, {х, 6 будут функциями от а и имеет место соотношение

I [х4 и-.), а, X (а), (1 (а), 6 (а)] =

= -±.Цх*

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

~ £ [хд

а а

где нижний индекс, как и прежде, означает частное дифференцирование по соответствующим величинам, ЬХУ £х, £а, столбцы,

-Ха, Матрицы. Из (3.2), (3.4) —(3.7), (2.3) следует, что первые

четыре члена в правой части (3.11) равны нулю, т. е. верно (3.10).

Предлагаемый итерационный алгоритм определения параметров а* записывается следующим образом:

а = а — и-1а,

где и — некоторый положительный параметр (регулирующий длину шага), Ь* вычисляется на каждом шаге по а по полученным в соответствии с (3.9) и из (3.8) значениям X и хА, ц, 6. Изменение а в соответствии с (3.12) будет продолжаться до тех пор, пока не будет удовлетворено—естественно, с некоторой точностью—условие (3.3) [условия (3.2), (3.4) — (3,7), (2.3) удовлетворены по построению алгоритма оптимизации], в результате чего задача нелинейного программирования (2.1)—(2.4) будет решена.

Итак, описанный алгоритм решения задачи (2.1)—(2.4) является декомпозиционным трехуровневым, т. е. варьирование параметрами и удовлетворение необходимых условий оптимальности (3.2)—(3.7), (2.3) осуществляется на трех уровнях оптимизации, имеющих определенный физический смысл. На первом, самом нижнем уровне посредством варьирования конструктивных параметров х(> 1—1, . . к, решаются подзадачи (3.8), в которых целевой функцией фактически является погонный вес силовой конструкции крыла плюс некоторые члены с весовыми множителями а, характеризующие вклад рассматриваемого /-го сечения в „связывающие** ограничения

(2.2), и удовлетворяются условия (3.2), (3.6), (3.7), (2.3). На втором уровне при варьировании параметров X для х, получаемых на каждой итерации по X (3.9) с первого уровня, удовлетворяются „связывающие" ограничения и условия (3.4), (3.5). На третьем, самом верхнем уровне посредством варьирования параметров а в соответствии с (3.12) и для х, X, р., 9, получаемых с предыдущих уровней при каждом значении а, удовлетворяются условия (3.3), т. е. фактически определяются оптимальные упругие деформации крыла при отклонении элеронов.

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

струкции крыла был запрограммирован на языке ФОРТРАН для ЭВМ БЭСМ-6. Расчет аэродинамических нагрузок осуществлялся по методу дискретных вихрей С. М. Белоцерковского [7] (программа расчета аэродинамической матрицы влияния была любезно предоставлена автору В. И. Емельяновым), в качестве алгоритма решения подзадач (3.8) для сечений использовалась методика, описанная в [5] (гл. 8, ч. I), безусловная минимизация в которой осуществлялась в соответствии с алгоритмом Нелдера и Мида ([5], гл. 7); интегральным ограничением вида (1.6) являлось ограничение на площадь миделя крыла [«==1 в формуле (1.6)].

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

(1.5); сдоп (г) принималось равным 37 кг/мм3 по всему размаху крыла, о„о —0,1 см, 8Обш0“0,2 см, 5про = 0,3 см, в формуле (1.6) 1^ = 20 м2, конструкция считалась изготовленной из сплава Д16Т.

Проиллюстрируем работу предложенной методики оптимизации на примере одного конкретного расчета. Функция аупр (г) описывалась полиномом 4-й степени, к=\ 1.

Сходимость алгоритма на втором уровне оптимизации при некотором фиксированном значении а

Таблица 2

Таблица I

Я ИТ А. Д - Фкр т

0 1,422 -з.з-ю-2 1,60 5

1 1,259 3,7-Ю-з 1,52 5

2 1,278 -7,9.10-* 1,54 2.5

3 1,276 —3,1-ю-* 1,53 —

М.Т X Фкр и

0 1,422 1,69 0,5

1 1,276 1,53 0,5

2 1,197 1,44 0,5

3 1,132 1,38 0,5

10 0,984 1,06 0,8

17 0,966 1,00 1

18 0,971 1,00 1

19 0,974 1 —

[условие (3.5) удовлетворялось с точностью до 5*10~4] показана в табл. 1, где пт — номер итерации по X, д — значение левой части „связывающего“ ограничения на площадь миделя крыла вида (2.2), Окр —безразмерный вес силовой конструкции крыла (за единицу принято оптимальное значение величины Окр). Сходимость алгоритма на третьем уровне оптимизации [условия (3.3) удовлетворялись по норме с точностью до 5*10~3, начальное значение [!£<*(! равнялось 0,5] проиллюстрирована в табл. 2, где А^ит— номер итерации по а, Х = Х (а)—итерированное значение, полученное для каждого а со второго уровня.

На рис. 1 показано, как по мере роста А^1Х сходится к оптимальному решению во всей области определения безразмерная изгибная жесткость Е1 (г), где г = 2]1нр [здесь и ниже безразмерные функции Е! (г) получаются путем деления их на оптимальное значение в корне крыла]. Из табл. 2 также видно, что по мере выполнения итераций (3.12) по а вес силовой конструкции крыла уменьшается. Расчет этого варианта для 44-х конструктивных параметров при 66 локальных и одном интегральном ограничении потребовал 17 минут времени счета центрального процессора БЭСМ-6, что

Рис. 1

говорит о достаточно высокой вычислительной эффективности разработанной методики.

Было проведено также исследование точности получаемого решения в зависимости от степени полинома р, аппроксимирующего функцию аупр (г), и от числа сечений вдоль размаха 6. Обозначим через в*/ вес силовой конструкции крыла, полученный в результате оптимизации при заданных величинах р и Для одного из расчетных случаев при 6=11 было получено (7”'3 = 0,62 0,ф2; Ой4 - 0,59 О^'2, о”-5 = 0,57 О!,1/ = 0,57 0,£г. _____

На рис. 2 на примере безразмерной изгибной жесткости Е1 {г) показано, как ведут себя искомые функции распределения параметров (при /7 = 5 и 6 эти кривые почти совпадают). Из приведенных результатов следует, что, во-первых, практическая сходимость по р имеется, и, во-вторых, при р—Ъ качественно правильно определяется закон изменения по 2 изгибной жесткости и при больших р происходит лишь его количественное уточнение. Приемлемый результат по весу конструкции получается уже при р — 4. По мере роста & имеет место также сближение соответствующих различным

Рис. 2

& искомых функций и веса конструкции. Так, для А =14, /7 = 8 (при этом значении р была достигнута сходимость по степени полинома) были получены изгибная жесткость, приведенная на рис. 2, и С™ = 0,56 01'р2.

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

На рис. 3 и 4 соответственно представлены распределения безразмерных изгибных жесткостей Е1 (г) и относительных толщин

профиля с {г) для следующих расчетных случаев: кривые /—по условиям статической прочности [1], кривые 2 — по условиям статической прочности и эффективности элеронов £0 = 0,2 при скоростном напоре <7 = 2,1 т/м2, числе М = 0,6 (положение элерона соответствовало 0,8—0,96 по координате г). Из приведенных графиков видно, что для обеспечения потребной величины эффективности элеронов необходимо увеличивать Е1 и с вблизи элерона за счет снижения жесткости у борта, вес продольного силового

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

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

Автор выражает благодарность В. М. Фролову за постоянное внимание и помощь в работе.

1. Б и р ю к В. И., Л и п и н Е. К., Ф р о л о в В. М. Методы проектирования конструкций самолетов. М., „Машиностроение*, 1977.

2. Бахвалов Н. С. Численные методы. М., „Наука", J975.

3. Geoffrion A. Elements of largescale mathematical programming, „Management Science", 1970, 16.

4. Ф и а к к о А., Мак-Кормак Г. Нелинейное программирование. Методы последовательной безусловной минимизации. М., 1972.

5. Гилл Ф., Мюррей У. Численные методы условной минимизации. М., „Мир", 1977.

6. Колмогоров А. Н., Фомин С. В. Элементы теории функций и функционального анализа. М., „Наука*, 1976.

7. Белоцерковский С. М., С к р и п а ч Б. К., Табачников В. Г. Крыло в нестационарном потоке газа. М., „Наука“, 1971.

Рукопись поступила 18jVJI 1980 Переработанный вариант поступил 151X11 1981