Содержание к диссертации
Введение
ГЛАВА I. Метод решения задачи обтекания крылового профиля 21
1. Введение 21
2. Метод граничных элементов для решения стационарных задач газовой динамики 23
3. Адаптивные расчетные сетки. Введение 38
4. Вариационный принцип построения расчетных сеток 42
5. Декартова адаптивная сетка 43
6. Криволинейная адаптивная сетка 48
7. Иные вариационные принципы 49
8. Способы учета вязкости 53
9. Тестовые расчеты 56
Рисунки к ГЛАВЕ I 58
ГЛАВА II. Способы представления варьируемых границ для задач аэродинамического проектирования 78
1. Введение 85
2. Параметрические полиномы четвертого порядка 86
3. Параметрические полиномы произвольного порядка 97
4. Интерполяционные свойства параметрических полиномов 99
5. Примеры представления геометрии границы 114
Рисунки к ГЛАВЕ II 116
ГЛАВА III. Методы и стратегия решения задач оптимизации и проектирования дозвуковых профилей 130
1. Общая постановка задачи оптимизации 130.
2. Алгоритмы оптимизации 133.
3. Примеры построения оптимальных профилей при различных геометрических ограничениях 136
Рисунки к ГЛАВЕ III 144
ГЛАВА IV. Оптимизация экспериментальных профилей 151
1 Оптимизация профиля NACA 642 —215 151
2. Оптимизация профиля П-2.9/31 - 17 155
3. Оптимизация профиля П- 2.3/75 - 17 157
Рисунки к ГЛАВЕ IV 160
ГЛАВА V. Многоточечная оптимизация крыловых профилей 179
1. Трехточечная ортимизация дозвукового профиля 179
2. Проектирование профиля типа "летающеелсрыло" 182
Рисунки к ГЛАВЕ V 185
Заключение 192
Литература
- Вариационный принцип построения расчетных сеток
- Параметрические полиномы четвертого порядка
- Примеры построения оптимальных профилей при различных геометрических ограничениях
- Проектирование профиля типа "летающеелсрыло"
Вариационный принцип построения расчетных сеток
Поскольку в предлагаемой работе для решения задачи обтекания используется метод граничных элементов (МГЭ), остановимся на некоторых вопросах с ним связанных. Подробно же история развития этих методов, начиная с работ Грина, Фредгольма, Михлина, Купрадзе по теории интегральных уравнений и до создания численных алгоритмов на их основе, которые стали интенсивно развиваться в 80-х годах и их приложениям к решению прикладных задач, представлена в монографиях [35,39]. Там же есть ссылки и на работы, посвященные применению МГЭ для решения задач аэрогидродинамики, относящиеся к течениям несжимаемой жидкости и линеаризованным трансзвуковым течениям газа. В отечественной литературе упоминания о разработке и применению МГЭ для расчета течений сжимаемого газа отсутствуют. А примеры использования этого метода при решении задач оптимизации отсутствуют и в литературе зарубежной.
Одним из вариантов метода граничных элементов (МГЭ) являются так называемые непрямые методы, в которых интегральные уравнения полностью выражаются через фундаментальное сингулярное решение исходных дифференциальных уравнений, распределенное с неизвестной плотностью по границам рассматриваемой области. Фундаментальное сингулярное решение может быть, например, функцией Грина для неограниченной области. Сами по себе функции плотности не имеют определенного физического смысла, но когда они найдены (численным решением интегральных уравнений), значения параметров решения во всей области могут быть получены из них простым интегрированием. Прямые же методы оперируют с неизвестными функциями, входящими в интегральное уравнение, имеющими физический смысл. Выбор метода, как правило, диктуется постановкой задачи.
Если в задаче для однородной области должны быть учтены распределенные объемные силы, или основные дифференциальные уравнения квазилинейны, то к граничным интегралам следует добавить объемный интеграл, включающий произвольное подразделение области. В этих случаях, однако, разбиение на подобласти не приводит к увеличению порядка окончательной системы алгебраических уравнений, подлежащей решению, и преимущества МГЭ сохраняются. Платой за это преимущество является заполненность матрицы, порождаемой при помощи МГЭ системы в отличие от МКЭ. И вычисление каждого элемента матриц при решении МГЭ приводит к большим арифметическим вычислениям, чем в методе конечных элементов, что компенсирует некоторое количество машинного времени, сэкономленного при решении системы. Однако, из предпринятых различными авторам исследований следует, что по мере роста размера задачи расходы для схем МГЭ, связанные с ЭВМ, растут менее резко, чем для иных методов решения. Эта разница еще больше для тех классов задач, которые особенно благоприятны для МГЭ, например, для систем, границы которых частично находятся в бесконечности. Поскольку процедура решения МГЭ автоматически удовлетворяет допустимым граничным условиям на бесконечности, разбиение этих границ не требуется.
Важным преимуществом является и то обстоятельство, что после решения интегрального уравнения, могут быть вычислены значения переменных, описывающие решение, в любой точке области. Более того, решение полностью непрерывно всюду в области. Эти особенности присущи только МГЭ и выделяют его среди возможных альтернатив.
Само по себе граничное интегральное уравнение является формулировкой поставленной задачи, ведущей к точному ее решению, и погрешности вследствие дискретизации и численных аппроксимаций возникают на границе и в области (если есть ее разбиение) из - за численного интегрирования. Но при использовании криволинейных элементов на границе и непрерывно меняющихся функций на ней эта процедура может быть весьма точной и погрешности, ей привносимые, достаточно малыми. И, конечно же, численное интегрирование всегда представляет собой более устойчивый и точный процесс, чем численное дифференцирование, и ни прямой, ни непрямой МГЭ не требуют никакого дифференцирования численных величин.
Еще один аспект, связанный с варьированием контура профиля и состоящий в выборе способа его аппроксимации и который практически не обсуждается в публикациях, состоит в том, что математическое представление геометрии тела должно исключать "паразитные осцилляции" при варьировании определяющих параметров приводящие к некорректности оптимизационной задачи. Эта проблема далеко непроста, поскольку практически все традиционные способы аналитического описания и представления в ЭВМ криволинейных границ имеют свои достоинства и недостатки. К последним можно отнести: тенденцию к осцилляциям у полиномиальных представлений, невысокую степень гладкости сопряжения для кусочно - полиномиальных представлений, необходимость задания производных в узлах для эрмитовой интерполяции. У сплайнов также могут возникать осцилляции на интервалах с большими градиентами и в точках разрыва кривизны, у параметрических сплайнов свободных от традиционных осцилляции возможно появление петель. Кривые Фергюссона - Бернштейна-Безье не могут одновременно удовлетворять требованиям высокой гладкости сопряжения и хорошего приближения к заданным узлам.
И конечно, успех в решении собственно оптимизационной задачи зависит от эффективной программы минимизации функции многих переменных при наличии функциональных ограничений в виде равенств и неравенств [45,4о]. Так в работах [148,151] использовался метод возможных направлений, в [149] одна из модификаций метода градиента, в [84] - метод проекции градиента. В работе [80] дается сравнительный анализ различных классов методов решения задач нелинейного программирования и предлагается подход, который является неградиентным методом с адаптацией и с использованием элемента случайности. Он включает в себя методы сканирующего конуса, покоординатного спуска на серии случайных испытаний, модификацию матрично - векторного спуска [59] и спуска в заданном напрвлении с использованием метода склона [78,80].
Параметрические полиномы четвертого порядка
Численные расчеты и аналитические исследования показывают, что адаптивные сетки могут существенно повышать точность и экономичность вычислительных алгоритмов. Они используются при расчетах задач в различных областях науки, но особенно значительный эффект достигается при решении задач механики. Это связано с тем что, как правило, в решении таких задач имеются зоны с большими градиентами и для получения приемлемого по точности решения необходимо иметь в этих зонах значительное число узлов расчетной сетки, что неприемлемо с точки зрения экономичности расчетов при равномерном распределении узлов сетки. Адаптацию применяют также для устранения осцилляции в решении и для уменьшения искусственной вязкости, характерной для многих численных методов решения уравнений. Интенсивные исследования в этой области ведутся фактически с начала активного использования ЭВМ при решении задач математической физики и продолжаются в настоящее время. Описание ряда алгебраических, дифференциальных и вариационных методов построения адаптивных сеток дано в монографиях, обзорах и статьях [81,82,112,133,146,147,171-174]. В данной работе анализ существующих методов построения адаптивных расчетных сеток дается с позиции их возможного применения для решения задач аэродинамического проектирования и оптимизации.
Можно сформулировать два основных требования к расчетным сеткам с учетом специфики этих задач.
Во-первых, при решении задач проектирования и оптимизации прямыми методами приходится многократно (до нескольких тысяч раз) решать задачу обтекания для получения информации о текущих значениях функционала качества и контроля выполнения ограничений. Следовательно, требование к вычислительной экономичности алгоритма генерации сетки в этих условиях приобретает особый вес.
Во-вторых, в процессе проектирования и особенно оптимизации, допускающей весьма широкий диапазон изменения параметров задающих контур обтекаемого тела, возможны варианты, не соответствующие стандартному априорному представлению о поведении решения. Следовательно, в алгоритм генерации сетки не следует закладывать информацию о конкретной качественной структуре решения, а основываться должен он на максимально общих критериях качества получаемого решения.
Естественно, что эти требования являются дополнением к общим требованиям, предъявляемым к расчетным сеткам, а именно: возможность глобальной и (или) локальной адаптации к решению, близость к ортогональности координатных линий выходящих с границ области для более точной реализации краевых условий Неймана, равномерность изменения шагов сетки [51,61,83,97,98].
Существует три основных класса сеток, применяемых для решения многомерных задач: 1) структурные, 2) неструктурные, 3) гибридные. Наиболее распространенными среди структурных сеток при решении задач математической физики являются координатные сетки. В двумерном случае их узлы определяются пересечениями линий некоторой координатной системы. Она может быть согласованной с границей или не согласованной в зависимости от того состоит ли граница из конечного числа координатных линий. Ясно, что координатными сетками, согласованными с границей в задачах обтекания тел произвольной формы являются криволинейные сетки, тогда как декартовые -не согласованы с границей. Естественно, что сетки,хорошо аппроксимирующие границу эффективны при численных расчетах задач, решения которых существенным образом зависят от точности аппроксимации краевых условий. Построение такой сетки происходит путем последовательного распространения дискретного преобразования с границы внутрь области. Это аналогично процессу интерполяции функции с границы или решению дифференциальной краевой задачи. Поэтому основными при построении этого типа сеток являются три группы методов: 1) алгебраические, основанные на разных видах интерполяции, 2) дифференциальные, основанные на решении, как правило эллиптических уравнений (хотя встречаются и гиперболические с параболическими), 3) вариационные, основанные на решении оптимизационных задач.
Алгебраические методы просты и экономичны, однако в сложных областях координатные линии, полученные существующими алгебраическими методами, могут пересекаться, ячейки вырождаться и, кроме того, они, как правило, сохраняют особенности границ, например изломы и т. д. [135,146,174].
Для построения сеток в областях с произвольными границами обычно используют дифференциальные методы, основанные на решении эллиптических и параболических уравнений [50-54,160,185]. Координатные линии в области являющиеся решениями этих уравнений представляют собой
БИБЛИОТЕК гладкие линии, что исключает распространение граничных изломов внутрь области. Однако итерационный характер алгоритмов построения этих сеток требует существенно больших затрат времени нежели методы алгебраические и весьма непрост в реализации.
Вариационные методы находят все более широкое применение при конструировании сеток, т. к. позволяют учесть одновременно несколько требований к сеткам: невырожденность, гладкость, близость к ортогональности, адаптивность и т. д. Однако, как и для дифференциальных методов, построение сеток на основе методов вариационных является сложной задачей, требующей специальных алгоритмов и, хотя известны предложения по использованию функционалов, минимизирующих погрешности решения и аппроксимаций уравнений, позволяющие получить уравнения для преобразования координат из необходимых условий минимума, реализации этого подхода в вычислительных процедурах для конкретных практических задач весьма сложны [67,68,152,164,182,184].
Для неструктурных сеток характерно нерегулярное расположение узлов, а также вложенность ячеек и их пересечение. Они удобны для локальной адаптации, осуществляемой за счет удаления или добавления узлов. Однако для многих численных алгоритмов решения краевых задач применение неструктурных сеток ведет к значительному усложнению алгоритма. Это происходит из-за ухудшения структуры матриц разностных операторов, необходимости обработки и хранения информации о нумерации узлов, их взаимном расположении. Они главным образом применяются в численных алгоритмах,основанных на методах конечных элементов[101,123,169].
И, наконец, в настоящее время широко распространяются гибридные сетки, сочетающие черты, как структурных сеток, так и неструктурных [43,53,114,157,181]. Такое сочетание осуществляется либо простым объединением структурных и неструктурных сеток, заданных в различных частях области, либо их композицией.
На основании изложенного можно сделать вывод о том, какие алгоритмы генерации сеток наиболее отвечают сформулированным выше требованиям: алгоритм построения сетки должен основываться на вариационном принципе, реализуемом прямыми методами численной оптимизации на множестве параметров, задающих сжимающие (и растягивающие) преобразования сетки. Это позволяет учитывать, если требуется, произвольное число дополнительных требований к генерируемой сетке, естественным образом включить процедуру минимизации функционала в итерационный процесс решения системы нелинейных уравнений в частных производных и в принципе исключает возможность возникновения нерасчетных ситуаций при решении оптимизационной задачи или задачи проектирования.
Были разработаны: декартова структурная адаптивная сетка не согласованная с границей сетка, криволинейная, структурная адаптивная сетка, согласованная с границей, неструктурная сетка, а также гибридные сетки, основанные на их композиции.
Примеры построения оптимальных профилей при различных геометрических ограничениях
В этих методах параметрами определяющими плоскую кривую являются координаты набора заданных точек (х0, у0), (х,, j,),..., (х v, уу), гдех, х2 - ... xN. Причем кривая, как правило, должна проходит через них, исключением являются задачи сглаживания исходных данных. При наличии iV+І точек можно записать интерполяционный полином, например, Лагранжа Р(х) степени N, который проходит через эти точки[32,36]. Однако при большом значении N возможны нежелательные осцилляции Р(х), которые могут иметь до N- 1 максимумов и минимумов, если все нули производной полинома Лагранжа Р (х) окажутся вещественными. Ясно, что для задания произвольной формы контура, или интерполяции какой-либо конкретной его формы с приемлемой для решения задачи обтекания точностью, требуется не один десяток таких узлов, а тенденция к осцилляциям увеличивается с ростом степени полинома Р(х), т. е. с увеличением числа интерполируемых точек. Кроме того, построение полинома носит глобальный характер, не в полной мере отвечающий задачам проектирования. Следовательно, для рассматриваемого класса задач этот и эквивалентные ему способы представления геометрии искомой границы не подходят. Один из способов устранения осцилляции - это построение составной кривой, в которой полиномы низкой степени (например, кубические) последовательно применяются для интерполяции групп точек. Полученная в результате кусочно - полиномиальная функция будет непрерывной, но в общем случае может иметь разрывы производных в точках соединения последовательных отрезков кривых. Для рассматриваемых приложений (проектирование дозвуковых профилей) это неприемлемо, т.к. наличие изломов на контуре профиля (за исключением задней кромки) приводит в них к нарушению дозвукового режима обтекания. Использование во всем диапазоне одного полинома невысокой степени, построенного по методу наименьших квадратов, в большинстве случаев не дает достаточно хорошего соответствия кривой заданным точкам. Другая возможность построения кривой связана с применением полиномов Эрмита, но для этого помимо значений функции у\ в точках Xj должны быть заданы и производные у, . Как правило, эта информация отсутствует. Если речь идет об интерполяции точек соответствующих конкретной модели, то это связано с трудностями точного измерения градиента. Если же речь идет о представлении кривой в задачах проектирования и оптимизации, то априорная информация контуре имеет лишь качественный характер: знаки кривизны, число перегибов и т.п. Вместе с тем, кусочная интерполяция полиномами Эрмита привлекательна тем, что обеспечивает непрерывность первых производных, не порождая неприятной проблемы осцилляции. Однако следует отметить, что вторые производные и, следовательно, кривизна будут разрывными в узлах. Таким образом, ясно, что в классическом виде эрмитова интерполяция также не решает ряд проблем, возникающих при задании варьируемого контура
И, наконец, существует еще один способ гладкой интерполяции данных с использованием полиномов невысокой степени при самых минимальных требованиях к знанию значений производных - это полиномиальные сплайны и их различные обобщения: сплайны с натяжением (экспоненциальные), рациональные сплайны, сплайны с дополнительными узлами и др.[1,38,44,48,57,58,65,66,70-73,87,106,141,161,163]. С физической точки зрения сплайн минимизирует энергию внутренних напряжений, с математической точки зрения кривая, ему соответствующая, имеет минимальную среднюю квадратичную кривизну. В декартовых координатах кривая минимальной энергии минимизирует интеграл, взятый между ее концами:
Определение соответствующей кривой у-у(х) из условия минимальности этого интеграла является вариационной задачей, не имеющей элементарного решения. Однако если считать, что \у\ 1 вдоль всей длины сплайна, получается более простая задача минимизации интеграла jy"2dx решением которой оказывается кусочно - кубическая функция, непрерывная вплоть до вторых производных, т.е. кубический сплайн в математическом смысле. Если задача ограничивалась бы только интерполяцией табличных данных, то с помощью местного переопределения системы координат можно рассматривать и кривые с большим наклоном, так как в локальной системе координат всегда удается обеспечить малость первой производной. Но увязать эту процедуру с расчетом течения в заданной системе координат весьма непросто. Можно также применять параметрический сплайн, который имеет то преимущество, что он не зависит от выбора осей координат, но необходимость пересчета локальных характеристик полученной кривой, необходимых при формулировании краевых условий в задаче обтекания, в исходной системе координат может породить нежелательные осцилляции. В практическом плане кубический интерполирующий сплайн конструируется по N+1 значению функции ) І и двум дополнительным условиям, задающим производные на концах интервала интерполяции. По существу интерполяция кубическим сплайном является приложением кусочной интерполяции Эрмита, в которой по исходному набору данных так определяются первые производные в узлах, что это обеспечивает непрерывность и вторых производных. Существуют также сплайны и более высокого порядка, чем кубический сплайн, однако повышение порядка гладкости производных в его узлах сопровождается и повышением вероятности возникновения нежелательных осцилляции. Для улучшения геометрических характеристик сплайнов, с тем чтобы их поведение согласовывалось с качественными характеристиками исходных данных, разработаны их различные обобщения, связанные с введением в структуру сплайна параметров для управления качественным поведением получаемой кривой. На этом пути достигнут заметный прогресс, но ценой значительной по объему предварительной работы по анализу исходных данных, построению промежуточной интерполяции для определения производных в узлах, их коррекции, введению дополнительных точек перегиба, и т.д. Кроме того, если не вводить в структуру сплайна как частный случай кусочно - линейную интерполяцию, то в общем случае участок кривой изогеометрической интерполяции, проходящий через три точки лежащие на одной прямой отличен от прямолинейного, что неприемлемо при интерполяции кормового участка крылового профиля. К тому же эти алгоритмы разработаны именно для параболических и кубических сплайнов, а с точки зрения возможности получения достаточно широкого набора конфигураций при использовании сплайнов не выше третьего порядка требуется значительное число узлов, что противоречит одному из требований к способу представления геометрии. Кроме того, недостатком классических сплайнов является то обстоятельство, что локальное изменение заданных значений в узлах влечет за собой вычисление заново всего сплайна. Это ограничение можно устранить используя В - сплайны. В общем в виде В - сплайн Mmi (х) порядка т (или степени т - 1 ) на данном множестве узлов везде равен нулю, кроме т последовательных отрезков х,. т х х,- и определяется этим множеством и одной величиной у. Принято исключать и эту степень свободы, фиксируя амплитуду В - сплайна некоторым стандартным образом, например, (Кокс и Де Бур, 1972): щ iMmi(x)dx=m \ Практическое значение В - сплайнов объясняется тем (Карри и Шенберг, 1966), что любой сплайн порядка т на множестве узлов хо ,х\,- ...,х# может быть выражен в виде линейной комбинации В - сплайнов, определенных на том же множестве узлов, расширенном дополнительными узлами на каждом из концов интервала: N+m-l S(x)= ] ,Mm;(x), ;=1 где S(x) - любой сплайн степени т - 1, с,- - числовые коэффициенты, определяемые из условий прохождения кривой через заданные точки. Ввиду того, что В - сплайн является только локально ненулевым изменение коэффициента с,- у одного из В - сплайнов влечет за собой изменение только на т отрезках кривой, позволяя производить локальные ее изменения без полного пересчета.
Проектирование профиля типа "летающеелсрыло"
Современные постановки задач характеризуются стремлением обеспечить максимальную степень адекватности модели отображаемому реальному процессу или объекту. Это порождает значительные трудности при решении оптимизационных задач, связанные с ростом размерности пространства поиска, сложной структурой ограничений и сложным поведением функционалов. Многообразие оптимизационных задач и отсутствие универсального, эффективного алгоритма их решения обусловило появление большого количества работ, посвященных созданию алгоритмов и программ оптимизации для различных задач, в том числе и для задач нелинейного программирования к которым относится и сформулированная выше общая задача оптимизации крылового профиля. При решении многих технических задач, требующих большого объема вычислений для расчета значения функционала, большую важность приобретает свойство быстродействия выбранного алгоритма, находящееся для таких задач в прямой зависимости от потребного количества вычислений значения функционала. По мнению ряда авторов[59, 78-80], при решении задач поиска в пространстве параметров полностью детерминированные методы не могут быть эффективными при реализации на ЭВМ. Причина, по - видимому, кроется в противоречии детерминированных предпосылок метода с ошибками реализации вычислительных процедур и точностью представления числа в ЭВМ. Примером может служить метод вращающихся координат Г. Розенброка, модифицированный введением элемента случайности в выборе направления спуска, что повысило его эффективность. С другой стороны, методы Монте-Карло, основанные на случайных выборках или методы поиска на равномерно -распределенных последовательностях, обладая преимуществами "глобальности" могут потребовать очень большого количества вычислений функционала, чтобы обеспечить достаточную уверенность в правильности полученного решения.
Представляется, что главным средством увеличения эффективности алгоритмов оптимизации является адаптация параметров алгоритмов в процессе поиска. Методы, используемые в данной работе можно охарактеризовать как неградиентные методы поиска с адаптацией и использованием элемента случайности. Они реализованы в виде специализированного комплекса программ синтезированного на основе модифицированных методов вращающихся координат, направляющего конуса и матричного спуска [59,95].
При решении оптимизационной задачи ограничения учитываются методом штрафных функций. Для этого вводится составной функционал, причем для коэффициентов штрафа используются оценки, гарантирующие существование локального экстремума составного функционала в области, где ограничения удовлетворяются с заданной точностью, что позволяет избежать "овражной " структуры минимизируемой функции [79], т.е. решается задача нелинейного программирования в следующей постановке:
Здесь є и є - задаваемые точности выполнения ограничений, К и К -коэффициенты штрафа, начальные значения которых задаются и уточняются в процессе работы программы. Отметим, что Ф0 должно быть больше единицы, например, Ф0 = (1 + Ф), где Ф - целевой функционал.
Задачи проектирования крыловых профилей, обладающих заданными свойствами, отличаются от задач оптимизации отсутствием в составном функционале целевой функции. Решение такой задачи определяется множеством допустимых значений варьируемых параметров удовлетворяющих заданному набору ограничений, не являясь, следовательно, единственным.
Ьолее того, так как и целевой функционал и ограничения нелинейным образом зависят от определяющих параметров задачи, то решения задач оптимизации также не единственны. Они зависят от выбранной начальной точки в пространстве варьируемых параметров, от требуемой точности выполнения ограничений, критериев прекращения работы программы и т. д.
В конечном счете, получение точных единственных решений экстремальных задач в аэрогазодинамике, как уже отмечалось, возможно лишь в случае простейших моделей и схем течения и ценность таких решений состоит в том, что они могут быть использованы для тестирования алгоритмов численного решения вариационных задач в предельных постановках. Примером такого тестового решения экстремальной задачи является классическая задача о максимуме подъемной силы гладкого контура, обтекаемого идеальной несжимаемой жидкостью. Точное решение - окружность. Максимальный коэффициент подъемной силы - 4/г, а передняя и задняя критические точки совпадают.
