осложненных условиях путем применения магнитных устройств // Нефтяное хозяйство. - 1997. - № 7. - С. 46-47.
5. Пивоварова Н.А., Унгер Ф.Г., Туманян Б.П. Влияние обработки постоянным магнитным полем на парамагнитную активность нефтяных систем // Химия и технология топлив и масел.
- 2002. - № 1. - С. 30-32.
6. Лесин В.И., Дюнин А.Г., Хавкин А.Я. Изменение физико-химических свойств водных растворов под влиянием электромагнитного поля // Журнал физической химии. - 1993. - Т. 67. -№ 7. - С. 1561-1562.
7. Борсуцкий З.Р., Ильясов С.Е. Исследования механизма магнитной обработки нефтей на основе результатов лабораторных и промысловых испытаний // Нефтепромысловое дело. - 2002.
- № 8. - С. 28-37.
8. Пат. 2153126 РФ. МКИ5 С01Б 111/22. Устройство для защиты трубопроводов от коррозии / В.М. Кондаков, А.Н. Качуров-
ский, А.Л. Бушковский, В.А. Кольцов, Л.В. Прасс, В.Н. Лялин. Заявлено 28.09.1998; Опубл. 20.07.2000, Бюл. № 20. - 16 с.: 7 ил.
9. Спектроскопия оптического смешения и корреляция фотонов / Под ред. Г. Камминса, Э. Пайка. - М.: Мир, 1978. - 574 с.
10. Посадов И.А., Поконова Ю.В. Структура нефтяных асфальте-нов. - Л.: ЛТИ, 1977. - 76 с.
11. Лоскутова Ю.В., Юдина Н.В. Влияние постоянного магнитного поля на реологические свойства высокопарафинистых неф-тей // Коллоидный журнал. - 2003. - Т. 65. - № 4. -С. 510-515.
12. Лоскутова Ю.В., Юдина Н.В. Реологическое поведение нефтей в магнитном поле // Инженерно-физический журнал. - 2006. - Т. 79. - № 1. - С. 102-110.
13. Лоскутова Ю.В. Влияние магнитного поля на реологические свойства нефтей: Дис. ... канд. хим. наук. - Томск, 2003. -144 с.
УДК 533.6.011
О ВЫБОРЕ ФУНКЦИОНАЛА ДЛЯ ОДНОЙ ВАРИАЦИОННОЙ ЗАДАЧИ ГАЗОВОЙ ДИНАМИКИ
В.М. Галкин
Томский политехнический университет E-mail: vlg@tpu.ru
Рассматривается численное решение прямым методом вариационной задачи о построении сверхзвукового сопла с равномерным потоком на выходе. Предложен способ выбора минимизируемого функционала. Проведено сравнение с результатом, полученным другим методом.
1. Введение
Известно, что если возможна численная (программная) реализация некоторой задачи, то эта реализация может быть сделана неединственным способом. В этом смысле не являются исключением и прямые численные методы решения вариационных задач, когда решение получается в результате минимизации определенного функционала. При этом для класса вариационных задач, к числу которых относятся и газодинамические задачи, некоторые численные реализации могут приводить к тому, что область определения будет являться нео-дносвязным множеством.
В качестве характерного примера рассмотрим задачу, решаемую прямым методом, о численном построении сверхзвукового сопла с равномерным потоком на выходе. Отметим, что к этому близка задача о сопле максимальной тяги [1]. Пусть в некотором сопле, профиль которого определяется варьируемыми переменными, рассчитывается поле течения, а по найденному полю вычисляется функционал, характеризующий неравномерность потока на выходе из сопла и имеющий минимальное значение тогда, когда поток равномерный.
В основу численной реализации указанной задачи могут быть положены два основных подхода. Первый использует тот факт, что поток остается
сверхзвуковым и возможно использование простых в реализации и быстрых в расчете маршевых схем; второй подход допускает существование дозвуковых течений и, соответственно, требуется применение гораздо более сложных и более медленных численных методов. Пусть в основу первой численной реализации положен метод характеристик [2]. Тогда очевидно, что если в некоторых контурах не будет реализовываться полностью сверхзвуковое течение, то это приведет к аварийному останову (авосту). Следствием этого будет неодносвязность области определения, поскольку функционал, использующий параметры потока на выходе из сопла, вычислить нельзя. Второй подход, учитывающий появление дозвукового течения и использующий, например, метод установления и схему Годунова [3], лишен это недостатка, однако интеллектуальные и временные затраты по сравнению с методом характеристик возрастают на порядок.
В работе [4] наряду с методом характеристик, позволяющим быстро получать решение, предложен функционал, который, в частности, использует значения параметров течения, найденные на каждой характеристике С+. Это позволяет даже при наличии авостов продолжить функционал на одно-связную область.
В данной статье, которая является развитием работы [4], предлагается вычислять функционал
более простым способом. Для оценки правильности полученных результатов использовалось решение, полученное по методике [2].
2. Постановка задачи
Дано установившееся незакрученное изоэнтро-пическое и изоэнергетическое течение идеального совершенного газа в осесимметричном сопле, рис. 1. Характеристические уравнения и условия совместности имеют вид:
^ = ^0±а),
ах
а0±-
С082(а)
а 5"1(а)5"1(0) ах = 0, (1)
(у +1)/2 - С082(а) у со8(0±а)
где у - показатель адиабаты, далее 7=1,4; х и у -продольная и поперечная координаты, отнесенные к радиусу сопла в минимальном сечении; а=агап(1/М) - угол Маха, М - число Маха; 0 - угол наклона вектора скорости (линии тока) к оси х; знак +(-) соответствует характеристике С+(С~); полагается, что во входном сечении ха=0. На стенке реализуется условие непротекания: tg(0)=/'(x), где /(х) - функция, описывающая контур сопла, а штрих обозначает производную по х. Во входном сечении ха поток однородный: 0п=0, Мп=1, где индекс «ш» соответствует входу в сопло. Далее индекс «0» соответствует выходу из сопла.
Рис. 1. Схема сопла: оа - минимальное сечение; аЬ - сопло с угловой точкой; аё - характеристика С-, принадлежащая пучку волн разрежения; ае - начальная характеристика ; к1 - характеристика С+ пт - характеристика С++; сЬ - конечная характеристика С
При вычислении параметров течения в угловой точке «а» используется формула Прандтля-Май-ера, следующая из (1):
62 +^(«2) = д1 +
у/(а) - -а + агС^
у — 1
у +1
tg(«)
(2)
где нижние индексы 1 и 2 соответствуют параметрам до и после разворота в угловой точке; 01=0, а=ап.
Необходимо найти контур сопла (рис. 1), на выходе из которого в сечении х=хь поток должен иметь однородные поперек сопла параметры 0=0 и число Маха М0>1. Считается, что давление на выходе из сопла больше давления окружающей среды.
3. Тестовая задача
При заданном показателе адиабаты рассмотрим однопараметрическое, зависящее от М0, семейство сопел с угловой точкой и равномерной характеристикой на выходе. Известно [1], что эти сопла имеют максимальную тягу, являются кратчайшими и обеспечивают нулевые потери на рассеяние в выходном сечении. Поскольку на выходе поток будет однородный и параллельный, то указанные сопла будут удовлетворять поставленным выше условиям.
Использование методики [2] при заданных М0>1 и показателе адиабаты позволяет найти единственное сопло, принадлежащее указанному семейству и удовлетворяющее условиям задачи. Контур сопла, полученный таким образом, считается эталоном, и с ним будет проводиться сравнение прямого метода, а координаты начальной точки «а» и найденной конечной точки «Ь» будут являться исходными данными для прямого метода.
4. Прямой метод
Даны координаты точек «а», «Ь» и тангенс угла наклона сопла в точке «Ь»:
/(ха ) = Уа, /(хЬ ) = УЬ, /'(хЬ ) = У^
ха = 0, Уа = 1, У' = 0. (3)
Последнее равенство следует из условия 0=0 на выходе из сопла.
Для аппроксимации искомого контура в качестве базисных функций используются степенные полиномы:
N+3
/(Х) = X 1 = (2Х - ХЪ - Ха )/(X - Ха ) ,
1=1
г е [-1,1], Х е[Ха, Хъ]. (4)
Так как должны выполняться условия (3), то коэффициенты ст, Сдт+2, ст и профиль сопла/(х) выражаются через N линейно независимых коэффициентов с1, ..., Сдт. Требуется найти профиль сопла, удовлетворяющий геометрическим условиям (3) и постановке задачи.
Прямой расчет. В рамках прямого метода принято называть прямым расчетом единичный расчет поля течения, который заканчивается вычислением функционала. Выбор функционала будет обсуждаться ниже. В прямом расчете система уравнений (1) решается по схеме, рис. 2, где цифрами «1» и «2» обозначены точки с известными параметрами, из которых выходят характеристики С + и С-. Цифрой «3» обозначена точка их пересечения. Записывая уравнения (1) в разностном виде, получим систему уравнений относительно неизвестных параметров а/, 0/, х/, у3 в точке «3»:
У3 - У 2
Уз - У1
= ^ё(023 -«2зХ
= 1Б(01з +«13),
соэ2а
(у +1)12 - соэ2 а
-(а]3 -а1) +
зта13 зт013
У13 СО§(013 +«13)
соэ а
(у +1)/ 2 - соэ2 а
(^ -*2) = 0,
*) = о,
К -а2) -
э1па23 зт023
У23 СО§(023 «23 )
здесь 7=1,2,... - номер итерации. Обозначим р=а,в,х,у, тогда Р1з=(Р1+Рзу-1)/2, Р2з=(Р2+Рз-1)/2. Полученная система решалась итерационно до выполнения условия тахр-рз^^О-8. На начальной итерации полагалось р30=(р1+р2)/2. Так как на оси 0=0, у=0, а на стенке 0=аг^(Дх)) и у=Цх), то в данной системе уравнений производились очевидные упрощения.
Рис. 2. Схема расчета: 13 - отрезок характеристики С+;
23 - отрезок характеристики С
Таким образом, от оси к стенке, по известной характеристике С+ рассчитывается следующая характеристика С ++1, где индексы I и /+1 - порядковые номера вычисляемых характеристик. Расчет производится до тех пор, пока очередная характеристика С не приходит в точку «Ь». В качестве начальных условий используются 0Ы=О и М=1,001 на начальной характеристике ае. Такие начальные условия достаточно часто применяют вместо плоской звуковой линии [5]. Отметим, что для вычисления интегралов использовалась формула трапеции.
Выбор функционала. В [4] использовался функционал, который рассчитывался вдоль характеристики (рис. 1):
I=л .
(5)
Если при вычислении характеристики С++1 происходил авост, то формула (5) заменялась выражением:
Г/ ь
I = (хь - х,) + П02й/ + |ф(х)йх,
У к а
ОI(х)|, I'(х) < о
ф( х) =
О, I'(х) > О
(6)
где первый интеграл вычислялся вдоль последней рассчитанной характеристики С+, а второй - суммировал площадь участков контура сопла с отрицательными углами наклона. Если авост произошел при расчете характеристики С+, то х=0 и 0=О вдоль начальной характеристики С2+. В результате функционал (6) приобретает вид:
ь
I = хь + |ф( х )йх.
(7)
Как видно, способ вычисления функционала (5-7), предложенный в [4], достаточно сложен. И, кроме этого, замена формулы (5) на формулу (6) приводит к тому, что функционал, сохраняя непрерывность, в общем случае может стать недиффе-ренцируемым.
В данной работе предлагается другой подход. Рассмотрим два функционала
11 =
12 = Ь1 ((а - «О )/а )2й/, Ь =| й/,
\ С+ С+
(8)
которые вычисляются вдоль произвольной характеристики С+. Выбор функционалов (8) обусловлен следующими соображениями.
Пусть сопло с известным Мо удовлетворяет условиям задачи. Тогда, как видно на рис. 3, первый функционал из (8) имеет два минимума 11=0. Первому минимуму на характеристике ае соответствует 0=О и ай=агс8т(1/Мя). Второму минимуму на характеристике сЬ соответствует 0=О и аО=аг8т(1/М0). Как следует из уравнений (1), только второй минимум является необходимым и достаточным условием решения задачи. Очевидно, что использование первого функционала не гарантирует единственности решения.
Рис. 3. Значение функционалов вдоль характеристик С+, / - номер характеристики
Второй функционал из (8), в отличие от первого, имеет единственный минимум 12=0 на характеристике сЬ (рис. 3), чему соответствует а=а0. Однако этот минимум является необходимым, но не достаточным условием решения задачи, поскольку может не выполняться условие 0=О. Комбинация I и 12 дает следующий функционал:
а
с
«о)/ а )2+02] dl.
(9)
который имеет единственный минимум /3=0 на характеристике сЬ (рис. 3). Так как этот минимум является необходимым и достаточным условием решения задачи и ему соответствует 6=0 и а=а0, то далее будет использоваться функционал (9).
При вычислении функционала (9) воспользуемся приемом, предложенным в [4], и если при вычислении характеристики С+ происходит авост, то функционал (9) вычисляется вдоль уже рассчитанной характеристики С+.
Если сопло имеет отрицательный наклон сопла Дха)<0 в точке «а», то авост происходит при вычислении характеристики С+, и функционал (9) вычисляется на начальной характеристике ае. Поскольку при вычислении параметров течения в угловой точке «а» используется формула Прандтля-Майера (2), то можно использовать эту формулу и для вычисления а2 при 62<0. В этом случае найденные значения а2 не имеют физического смысла, т. к. этому соответствует дозвуковое течение. Однако в силу своей непрерывности (рис. 4) функция а2(62) позволяет продолжить функционал (9) на односвя-зную область в случае отрицательного наклона сопла Дха)<0 и обеспечить функционалу непрерывность и дифференцируемость. Как показали расчеты, при вычислении по формуле трапеции функционала на начальной характеристике ае достаточно брать только две точки - «а» и «е».
Зависимость а2=а2(62)
Поскольку на входе в сопло граничные условия не меняются, а параметры течения определяются только профилем сопла (4), то значение функционала (9) будет неявно зависеть от этого профиля. Таким образом, задача нахождения профиля сопла У(х), доставляющего экстремум функционалу (9), сводится к поиску точки (с1,^,сЛГ), в которой функция многих переменных /3=/3(с1,^,сж) имеет экстремум. Для нахождения минимума этой функции использовался квазиньютоновский метод Бройдена из [6].
5. Численные результаты
Число N варьируемых переменных с менялось от 1 до 10, начальное значение с=0. На характери-
стике ае задавалось 50 точек. Координаты точек «а» и «Ь»: (0; 1) и (3,576; 1,299). Этим значениям соответствует найденное из [2] эталонное сопло с угловой точкой и числом Маха на выходе М0=2.
Сравнение решения для обоих функционалов с эталонным соплом показало, что при N=10 ординаты сопел отличались в четвертом знаке после запятой, а максимальная относительная погрешность по составила 0,02 %.
На рис. 5 представлены результаты минимизации в виде окончательных значений функционалов для различного числа N. Таблица демонстрирует влияние числа коэффициентов и используемого функционала как на число прямых расчетов, так и на количество авостов. Из представленных результатов видно, что предложенный функционал (9), несмотря на простоту, по эффективности не уступает функционалу (5-7), предложенному в [4].
Рис. 5. Минимизация функционалов: 1) формулы (5-7) из [4], 2) предлагаемая формула (9)
Таблица. Минимизация функционалов. В числителе - число прямых расчетов, в знаменателе - авостов
№ Формулы (5-7) [4] Предлагаемая формула (9)
1 13/2 17/2
2 36/3 39/4
3 63/4 69/5
4 98/4 99/5
5 135/5 115/7
6 150/7 166/7
7 199/7 203/6
8 273/6 304/7
9 385/8 320/8
10 444/7 492/9
Заключение
Проведенные численные исследования показывают, что предложенный функционал достаточно прост в реализации, при наличии авостов имеет продолжение на односвязную область и сохраняет дифференцируемость, а по эффективности не уступает ранее рассмотренному в работе [4] функционалу. Сравнение профиля найденного сопла с эталонным профилем показало, что при 10 варьируемых переменных максимальная относительная погрешность по ординате составляет 0,02 %.
СПИСОК ЛИТЕРАТУРЫ
1. Крайко А.Н. Вариационные задачи газовой динамики. - М.: Наука, 1979. - 448 с.
2. Кацкова О.Н. Расчет равновесных течений газа в сверхзвуковых соплах. - М.: ВЦ АН СССР, 1964. - 61 с.
3. Годунов С.К., Забродин А.В., Иванов М.Я. и др. Численное решение многомерных задач газовой динамики. - М.: Наука, 1976. - 400 с.
4. Галкин В.М., Волков Ю.С. Сравнение базисных функций в прямой задаче профилирования сверхзвуковой части сопла //
Сибирский журнал индустриальной математики. - 2004. - Т. 7. - № 4(20). - С. 48-58.
5. Крайко А.Н., Тилляева Н.И. Оптимальное профилирование контура сверхзвуковой части тарельчатого сопла // Механика жидкости и газа. - 2000. - № 6. - С. 172-184.
6. Дэннис Дж. мл., Шнабель Р. Численные методы безусловной оптимизации и решения нелинейных уравнений. - М.: Мир, 1988. - 440 с.
УДК 536.21+692.2:691.327:666.973.2:666.64-492.3
НЕСТАЦИОНАРНЫЙ ПРОСТРАНСТВЕННЫЙ ТЕПЛОПЕРЕНОС В НЕОДНОРОДНОЙ КЕРАМЗИТОБЕТОННОЙ СТЕНЕ
А.Н. Хуторной, А.Я. Кузин, Н.А. Цветков, Т.А. Мирошниченко, А.В. Колесникова
Томский государственный архитектурно-строительный университет E-mail: kaftgs@tsuab.ru
С помощью математического моделирования исследовано тепловое состояние неоднородного трехслойного керамзитобетон-ного фрагмента стены с фасадным утеплением на гибкой связи. Установлен характер распределения полей температуры в ке-рамзитобетонной конструкции. Выполнена оценка зоны влияния гибкой связи на температурное поле стены. Разработанная численная методика позволяет проводить тепловую экспресс-диагностику наружных неоднородных керамзитобетонных стен с различными теплофизическими и геометрическими характеристиками материалов системы в реальных условиях эксплуатации.
Вопросы энергосбережения относятся к приоритетным направлениям развития науки, технологий и техники Российской Федерации. Очевидно, что решение этих вопросов должно иметь комплексный характер. В работе [1] выполнена оценка энергоэффективности комплекса энергосберегающих мероприятий в нескольких зданиях общественного назначения г. Москвы. Суммарная экономия энергии оказалась весьма значительной (56...63 %), при этом доля утепления несветопро-зрачных ограждений составила 16.20 %. Там же [1] подчеркнута важность энергосберегающих мероприятий, связанных с утеплением несветопрозрач-ных ограждений, как обеспечивающих суммарное снижение энергопотребления здания не менее чем в 2 раза. В этой связи разработка способов повышения теплоэффективности наружных стен существующих зданий и разработка новых конструкций с повышенными теплозащитными свойствами является объективной необходимостью.
Повысить уровень теплозащиты наружных стен зданий можно за счет устройства в толще конструкции пустот, заполненных эффективным утеплителем [2]. Однако, как показано в работе [3], для ряда климатических зон такое решение может быть не достаточным. С целью обеспечения требованиям СНиП 11-3-79* «Строительная теплотехника» возникает необходимость в устройстве дополнительного утепления ограждающих конструкций, которое может быть реализовано с помощью фасадных систем утепления на гибких связях.
Физико-математическая постановка задачи. Рассмотрим нестационарный теплоперенос через трехслойную неоднородную конструкцию, состоящую из керамзитобетона - 1 с вертикальной пустотой, заполненной утеплителем - 4, утеплителя фасадной системы утепления - 2 и обшивки - 3 (рис. 1). Керамзитобетон с обшивкой соединены гибкой связью (коннектором) - 5. Форма керамзитобето-на, внутреннего утеплителя, коннектора, утеплителя фасадной системы утепления и обшивки - прямые параллелепипеды, поперечные и продольные сечения которых в общем случае представляют собой разносторонние прямоугольники. Известны те-плофизические характеристики (А;, р, с,, /=1,5) материалов системы, ее геометрические размеры, температуры наружной (^£) и внутренней (г1^) сред, коэффициенты теплоотдачи на наружной (о,) и внутренней (а0) поверхностях ограждения [4].
Необходимо рассчитать поля температур и плотностей тепловых потоков в сечениях неоднородной керамзитобетонной стены с фасадным утеплением.
Теплоперенос в неоднородном керамзитобетон-ном фрагменте стены описывается в декартовой системе координат трехмерными нелинейными нестационарными уравнениями теплопроводности
(рс) _ д п ) + 5 п ) + 5 (, д)
i = 1,5.
(1)