CC BY
57
УДК 532.517 В. Б. Филиппов
Вестник СПбГУ. Сер. 1, 2004, вып. 1 (№ 1)
МОДЕЛЬ МОРСКОГО ГАЗОПРОВОДА С УЧЕТОМ ОЛЕДЕНЕНИЯ
1. Введение.
В работах [1-5] исследовалась модель неизотермического установившегося турбулентного течения газа по газопроводу постоянного кругового сечения, проложенному по дну моря, без учета его возможного оледенения. В модели указанных процессов принималось, что на внешней боковой поверхности газопровода выполняется граничное условие первого рода
1(Еп,г) = Т *(г),
Ь — температура в многослойной боковой поверхности, зависящая от радиуса г и координаты г вдоль оси газопровода, (г, в, г) — цилиндрическая система координат; Яп — внешний радиус газопровода, поверхность которого состоит из п слоев; Т* (г) — температура окружающей воды.
В северных морях температура воды бывает близка к температуре замерзания Т*. При низких температурах газа в потоке близость температуры окружающей воды к температуре замерзания может привести к оледенению поверхности газопровода. Представляет интерес модель транспортировки газа, позволяющая учесть это явление. В настоящей работе приведены постановка и алгоритм решение этой задачи. В стационарных условиях предлагается искать методом установления толщину слоя льда 6л, соответствующую установившемуся режиму течения газа с низкими отрицательными температурами. Граничное условие на поверхности газопровода с учетом оледенения содержит неизвестную величину ¿л — толщину слоя льда. На каждой итерации решается задача о течении газа при толщине льда, рассчитанной на предыдущей итерации. Затем по найденной температуре газа пересчитывается толщина льда в новом приближении и осуществляется переход к следующей итерации. Как показали расчеты, этот процесс сходится и его пределом является решение стационарной задачи о течении газа по морскому газопроводу с учетом оледенения. Предварительно необходимо решить задачу расчета толщины льда по заданной неизменной температуре газа.
2. Модель 1. Расчет толщины слоя льда на поверхности газопровода
при независящих от времени температурах газа Т и морской воды Т*
В стационарном режиме толщина слоя льда 6л определяется из условия равенства тепловых потоков 3л (от льда к границе) и зв (от воды к границе) на границе лед—вода.
Тепловой поток 3л определяется соотношеним
А гМ
л 3,г
, (1)
Ал — коэффициент теплопроводности льда; Яп — внешний радиус газопровода,
п
Яп = Я + ^^ ¿к,
к=1
© В. Б. Филиппов, 2004
3
л
Д — внутренний радиус газопровода, 5к —толщина к-го слоя боковой поверхности, п — число слоев.
Тепловой поток со стороны воды к границе лед—вода определяется соотношением
Л ¿г 3в = -Ав ~т
ат
, (2)
Ав — коэффициент теплопроводности соленой морской воды.
Постановка задачи расчета теплового потока ]л состоит из стационарного уравнения теплопроводности (3), дополненного граничными условиями первого рода (4) и (6) на внутренней и внешней (с учетом льда) поверхностях газопровода и условиями
к
равенства тепловых потоков (5) на стыке слоев при т = Дк = Д к = 1,...,п.
Пусть 5л задано.
а ( аг\
(3)
г(Д) = т, (4)
аг
Лк~Г
ат
аг
- Ад;+1 —
, к =1,...,п (5)
Пк
г(Дп + 5Л) = Т*, (6)
Ак — коэффициент теплопроводности к-го слоя (Ап+1 = Ал); Т — температура газа; Т* — температура фазового перехода вода—лед.
Граничное условие на поверхности льда (6) отражает тот факт, что при наличии льда температура на границе вода—лед должна быть равна температуре фазового перехода Т*, независимо от значений температуры газа Т и температуры окружающей воды Т*, которые влияют лишь на величины тепловых потоков ]л и ]в.
Интегрирование уравнения теплопроводности (3) в области т € [Д, Дп + 5л] при
условиях (4)—(6) приводит к следующему выражению для ]л (1):
=___ (7)
М (Дп + 5Л)(А + Ы(1 + 5Л/Дп)/\Л)'
п
А = ^1п(1 + 4/Дк-г) /Ак, (8)
к=1
До = Д — внутренний радиус газопровода.
Постановка задачи расчета теплового потока ]в
В постановке этой задачи необходимо учесть, что, во-первых, плотность воды при остывании в окрестности температуры фазового перехода Т* ведет себя немонотонно. Кроме того, при замерзании изменяется соленость воды, что также приводит к изменению ее плотности. Все это приводит к возникновению вблизи поверхности градиента плотности, приводящего к появлению конвекции. Во-вторых, условия контакта газопровода с грунтом не осесимметричны и условия теплообмена на части поверхности
газопровода, контактирующей с водой, существенно отличаются от условий теплообмена на части поверхности, контактирующей с грунтом.
Перечисленные факты приводят к тому, что задача расчета теплового потока ]в в реальных условиях должна учитывать 1) наличие на обнаженных участках поверхности газопровода дополнительных (к донным течениям) конвективных потоков и 2) отсутствие аксиальной симметрии процессов за счет контакта с грунтом.
В рамках осесимметричной модели возможен следующий подход: расчет теплового потока ]в проводится по модельной задаче осесимметричного обтекания газопровода, но в качестве толщины теплового погранслоя воды выбирается эффективный параметр 5,. Проинтегрированное по углу в решение модельной задачи с эффективным параметром 5, должно совпадать с решением реальной задачи расчета ]в, учитывающей перечисленные выше факторы. Вся сложность при этом переносится на отыскание величины эффективного коэффициента 5,. Достоверная информация о его величине может быть получена из решения так называемой обратной .задачи, постановка которой базируется на решении прямой задачи расчета всех параметров газового потока и на экспериментальных данных, полученных в условиях максимально приближенных к реальным.
Модельная .задача расчета состоит из уравнения теплопроводности (3), граничного условия (6) на внешней (с учетом льда) поверхности газопровода и граничного условия на границе эффективного слоя воды
1(Дп + бЛ + 5,) = Т *. (9)
Интегрирование уравнения теплопроводности (3) с граничными условиями (6), (9) приводит к следующему выражению для ]в (2):
. =__(Т*-П)_
Зв (Дп + 6л) 1п(1 + (5*/(Д„ + ¿л)) / Ав ' 1 1
В установившемся режиме тепловые потоки ]л и ]в должны быть равны, это, совместно с найденными для них выражениями (7), (10), позволяет записать уравнение, определяющее толщ,ину слоя льда 5л на внешней поверхности газопровода в установившемся режиме при осесимметричной постановке для неизменных температур газа Т и окружающей воды Т*:
Т* -Т« = %-Т
1п(1 + <5*/(Д„ + <5л))/Ав А + 1п(1 + 6Л/Д„)/Ал ' 1 ;
величина А определена равенством (8).
Уравнение (11) в неявной форме задает зависимость толщины слоя льда 5л от температуры газа Т, температуры воды Т*, температуры фазового перехода Т, и других параметров задачи Ав, Ал, , 5 ^, п, Д, 5, .В модели транспортировки газа по морским газопроводам все перечисленные параметры, кроме температуры газа Т, считаются известными величинами, возможно зависящими от координаты г. Таким образом, (11) служит уравнением для расчета толщины 5л льда в зависимости от температуры газа Т в установившемся режиме.
Оценка величины параметра 5,
Исследования динамических и тепловых погранслоев при обтекании тел различной формы обобщены в виде степенных зависимостей числа Нуссельта N4 от чисел
Рейнольдса Re, Прандтля Pr, Грасгофа Gr и Эккерта Ec, например, [6]. Приведем характерные значения чисел Re, Pr, Gr, Ec в рассматриваемой задаче обтекания газопровода, рассчитанные для воды при нуле градусов Цельсия, при характерной длине l = 0.6 м и при характерной скорости обтекания 0.01 м/с:
Re = 6 • 103, Pr =16, Gr = 3 • 107, Ec = 2 • 10-8.
Эти значения чисел Re, Pr, Gr, Ec свидетельствуют о том, что реализуется турбулентный режим обтекания и имеет место смешанная тепловая конвекция (вынужденная за счет придонных течений и естественная за счет изменений плотности воды). При этих значениях чисел Re, Pr, Gr, Ec средняя по всей боковой поверхности величина числа Нуссельта Nuc при обтекании цилиндра может меняться в диапазоне [6]
80 < Nuc < 110.
Выразим тепловой поток jB (2) через коэффициент теплопередачи в
Зв = в(Т* - T*)
и учтем связь числа Нуссельта с коэффициентом теплопередачи
Nuc = fîlK1.
Сравнение выражения (10) для теплового потока jB, найденного в модельной задаче, с выражением jB через козффициент в позволяет выразить число Нуссельта Nuc через S** — толщину теплового погранслоя при обтекании цилиндра, не соприкосающегося с грунтом.
Nuc = (ln(1 + S^R-1))-1,
R* = Rn + — внешний радиус газопровода с учетом слоя льда. При найденных значениях Re = 6 • 103 и Pr = 16 характерная величина S** при обтекании цилиндра не превышает 0.02 м [6], это позволяет упростить зависимость числа Nuc от параметра S**. Разложим логарифм в ряд по малому параметру S**R-1 и ограничемся первым членом ряда, в результате получим
Nuc = R*S-1
Пусть газопровод частично погружен в грунт и обтекается только часть S* боковой поверхности S (S* = yS, y € [0,1] — доля обтекаемой поверхности).
В осесимметричной модели принимается, что обтекается вся боковая поверхность, но с эффективной толщиной S* теплового погранслоя, связанной с S** равенством
S* = S** y-1.
Найдем величину S* , положив число Нуссельта равным 100, и считая, что газопровод на половину погружен в грунт:
S** = R*Nu-1 = 6 • 10-3м, S* = S**(0.5)-1 = 1.2 • 10-2м.
Конечно, приведенная оценка величины S* носит приближенный характер и достоверная информация о ней, как отмечалось выше, может быть получена только из решения обратной задачи по соответствующим экспериментальным данным.
В качестве примера в таблице 1 приведены результаты расчета зависимости 6л(Т) для разных значений параметра ¿*. Остальные параметры в системе единиц СИ приняты равными
Ал = 2, 3; Ав = 0.6; Т* = 272.15; Т* = 271.15;
Еп = 0.5795; А = 0.0886631. ( )
При выборе значений Т*, Ав учитывалась соленость морской воды. В равенствах (12) значения температуры даны в градусах Кельвина.
Таблица 1
St = 0.05 M St = 0.01м <5* = 0.005 м
т° с ¿л M т° с S л М т° с S л М
-7 0.53 -10 0.14 -20 0.17
-6 0.42 -9 0.11 -15 0.10
-5 0.31 -8 0.09 -12 0.06
-4 0.19 -7 0.06 -10 0.03
-3 0.05 -6 0.03 -9 0.01
-2.7 0.01 -5.2 0.003 -8.5 0.005
3. Модель 2. Модель неизотермического установившегося турбулентного течения сжимаемого газа по морским газопроводам с учетом оледенения
Сформулированная задача в осесимметричной постановке без учета оледенения исследовалась в работе [5], ее одномерный вариант с другими термодинамическими замыкающими уравнениями рассмотрен в книге [7]. Для учета оледенения выпишем, следуя работе [5], замкнутую модель процессов в осредненной по сечению форме и проследим, как влияет ситуация на внешней поверхности газопровода. Постановка указанной задачи имеет вид
Ki dp dp К2
~п ~г = ~г Н---рд cos az,
p2 dz dz p
Qpc,
dT Qk2
dz
p
+
PQ dP , „2 ( 2 ~Tz+P7<R \-Rqw
cp2
+
pQ 3c dp
2(1 + Sp)T1/2 dz '
hpT
l-Sp ~ (1 + 5p)T1/2
(13)
(14)
(15)
(16)
p(z), p(z), T(z) — плотность, давление и температура газа в z-ом сечении, (13) — осред-ненное по сечению газопровода уравнение движения в проекции на ось z. В уравнении (13) слагаемое (^r^f) отражает влияние инерционных сил, слагаемое — сил давления, слагаемое — сил турбулентных касательных напряжений, слагаемое (-pgcos az) — сил тяжести. Константы ki, К2 определяются профилем функции p(r) = p(z)u(r, z), где u(r, z) — скорость потока (расчет функции может быть отделен от решения общей задачи [3]):
p|z =
z=0
po,
TI
z=0
T0
ki
2
Д2
^2(r)r dr,
(17)
R
= (18) 4тг2Д5 ^ '
Коэффициент сопротивления Л, входящий в (18), равен
А = 8^, Риа
иа — средняя по сечению скорость потока, тт — касательное напряжение по стенке. Коэффициент Л определяет величину константы «2 (18) и может быть найден либо из экспериментального закона сопротивления Л = Л( Ив ,к), либо рассчитан по профилю функции у: Ив — число Рейнольдса, определенное по средней скорости иа
2Q р Ив
Q — весовой расход газа, считающийся в рассматриваемых задачах постоянной величиной, р — вязкость газа, к = кэ/К — безразмерный коэффициент относительной шероховатости, кэ — коэффициент эквивалентной равномерно-зернистой шероховатости. В уравнение движения (13) входят также д — ускорение силы тяжести и аг — угол между направлением силы тяжести и осью газопровода.
В уравнении баланса внутренней энергии (14) слагаемое С}рсу ^ отражает изменение внутренней энергии газа за счет конвективного переноса, слагаемое Qк2/p — за счет диссипации, слагаемое ^ —за счет работы сил давления, слагаеме /Э7гД2( — qw) —
за счет теплообмена с внешней средой, слагаемое 2{1+бр)т1/2 Лг — за счет неидеальности газа.
К обсуждению величины теплового потока , отражающего условия теплообмена с внешней средой, вернемся далее. Вид последнего слагаемого правой части уравнения (14) обусловлен уравнением состояния Редлиха—Квонга (15), в котором константы Н, с, 6 определяются химическим составом газа. Кроме перечисленных величин в уравнение (14) входит коэффициент удельной теплоемкости су газа.
При известной величине теплового потока система уравнений (13)—(15) с граничным условием на входе в газопровод (16) представляет собой замкнутую постановку задачи расчета плотности р(г), температуры Т(г) и давления р(г) сжимаемого неидеального газа в установившемся неизотермическом турбулентном течении по газопроводу постоянного кругового сечения.
В модели (13)—(16) константы к\, к считаются известными. Их расчет, зависящий от принятой реологической модели, рассмотрен, например, в нашей работе [3].
Остановимся подробнее на величине характеризующей теплообмен газа с окружающей средой. Допустим, что один из коэффициентов теплопроводности теплоизолирующих слоев поверхности газопровода меньше коэффициента Лт турбулентной теплопроводности газа. В этой ситуации лимитирующей стадией теплообмена газа с внешней средой является теплопроводность через многослойную боковую поверхность. Это позволяет учесть теплообмен с внешней средой интегрально, в форме объемного источника (стока) внутренней энергии ш(г), выраженной через — радиальную составляющую вектора потока тепла на внутренней поверхности газопровода
аг
Чгю — —Лт~Г
ат
к
по формуле
Постановка задачи расчета (19) в установившемся режиме практически совпадает с постановкой (3)—(6) задачи расчета теплового потока (1) с той лишь разницей, что здесь искомой величиной является поток тепла не на внешней границе боковой поверхности, а на внутренней при т = К, причем из условия равенства тепловых потоков следует
сИ Хт-г
ат
_ А
я 1<1г
, (21)
к
Лх — коэффициент теплопроводности первого слоя.
Решение задачи расчета для трех слоев без учета возможного оледенения при граничном условии первого рода на внешней поверхности газопровода получено в работе [1].
В настощей работе найдем выражения для величины в двух ситуациях.
I. Боковая поверхность газопровода состоит из п слоев с известными геометрическими и теплофизическими свойствами и из слоя льда заданной толщины §л; на внешней
поверхности слоя льда, граничащей с водой, выполняется условие первого рода
= Т*, (22)
Т* — температура фазового перехода вода—лед.
II. Боковая поверхность газопровода имеет ту же конструкцию, что и в I случае, но слой льда отсутствует, и на внешней поверхности газопровода выполняется граничное условие третьего рода с заданным коэффициентом теплопередачи в:
= -в[Аи - Т*) , (23)
аг
ат
Т* — температура окружающей воды.
Значение коэффициента теплопередачи в оценивалось по данным работ [6], [8]. В I случае решение уравнения теплопроводности (3) в области т € [К, Кп + £л] относительно потока тепла при условиях (4)—(6) приводит к выражению
, -__сТ*~Т)__(24)
Величина А определена равенством (8).
Решение той же задачи во II случае приводит к выражению
_ РЕп(Т*-Т) <Ь~ Щ1 + 1ЗПпА)' [ '
В этих равенствах температура газа Т является функцией г, температура окружающей воды Т* может зависеть от г. Выражение для (24) при наличии льда явно не содержит температуру Т* окружающей среды, она в этом случае входит через зависимость от нее толщины льда 6л. Выражение для , полученное в работе [1], является частным случаем (25) при п = 3, в
При отсутствии льда постановка задачи расчета установившегося неизотермического турбулентного течения неидеального сжимаемого газа по морскому газопроводу (13)—(16) замкнута, так как величина входящая в уравнение баланса внутрененй энергии (14), задается выражением (25) и является известной функцией температуры газа Т(г) при заданных в, Я, Яп, А, Т*(г).
Иначе обстоит дело при появлении ледяного покрова. Величина в той области газопровода, где есть слой льда, определяется выражением (24), в которое входит неизвестная величина 6л. Ее нельзя рассчитать непосредственно по модели 1 (уравнение (11)), так как в рассматриваемой задаче температура газа Т сама зависит от 6л. Предлагается решать задачу методом установления с помощью следующей итерационной процедуры.
4. Итерационная процедура решения задачи (13)—(16)
с учетом оледенения
При заданной толщине ледяного покрова задача (13)—(16) может быть решена численно, например, методом Рунге—Кутты. Уравнение (11) модели 1 позволяет рассчитать толщину льда 6л(г) по заданной температуре газа Т(г) в стационарных условиях. Наличие моделей 1 и 2 и алгоритмов расчета по ним позволяет предложить следующую итерационную процедуру решения задачи методом установления.
Нулевое приближение. Постулируем отсутствие льда, ставим граничное условие (23) на внешней поверхности газопровода и задачу (13)—(16) замыкаем выражением для (25). В результате решения этой задачи находим в нулевом приближении распределение температуры Т(0)(г), плотности р(0)(г) и давления р(0)(г) газа во всем газопроводе г € [0,Ь], Ь — длина газопровода. Уравнение (11) и найденное распределение температуры Т(0)(г) газа позволяют определить координату г*, начиная с которой на поверхности газопровода должен быть лед, кроме того, они позволяют рассчитать толщину льда г) в нулевом приближении для г € [г*,Ь].
Первое приближение. В газопроводе выделяются две области: [0,г*] и [г*,Ь]. В области [0, г*] решение задачи во всех приближениях совпадает с ее решением в нулевом приближении:
г € [0, г*] : р(к)(г) = р(0)(г), р(к)(г) = р(0)(г), Т(к)(г) = Т(0)(г), (26)
к — номер приближения.
Во второй области [г*, Ь] задача (13)—(16) замыкается для уравнением (24), учитывающем наличие льда, при этом в качестве 6л в уравнении (24) берется найденное в нулевом приближении распределение 6(0)(г). В результате решения задачи (13)—(15) в области [г*, Ь] с граничным условием
Р(1) (г*) = р(0)(г*), Т(1)(г*) = Т(0)(г*) (27)
рассчитываются распределения плотности р(1)(г), температуры Т(1)(г) и давления р(1) (г) газа в первом приближении. По найденному распределению температуры Т(1)(г) газа из уравнения (11) находится распределение толщины слоя льда ¿л^г) в первом приближении.
Второе приближение. Как и в первом приближении, решение задачи уточняется только в области [г*, Ь]. В уравнение (24) подставляется распределение толщины льда
5Л1}( г), найденное в первом приближении. После этого задача (13)—(15), (24) решается во второй области с граничным условием типа (27), в результате чего находятся распределения плотности р(2)(г), температуры Т(2)(г) и давления р(2)(г) газа во втором приближении.
По распределению температуры Т(2)(г) из уравнения (11) находится распределение толщины слоя льда 5Л2)( г) во втором приближении.
Структура последующих приближений аналогична рассмотренным.
Условием окончания итерационного процесса служит выполнение неравенств
max
z* <z<L
max
z* <z<L
T(k)(z) - T(k-1) (z) S(k)(z) - Sik-1) (z)
< £
(28)
< £9
в которых k — номер приближения; £i, £9 — малые заданные величины.
5. Выводы
Предложена модель неизотермического установившегося турбулентного течения газа по морскому трубопроводу при его частичном оледенении. Получено решение изотермической задачи расчета толщины ледяного покрова морского газопровода. Рассмотрено решение методом установления неизотермической задачи о транспортировке газа с учетом оледенения. Приведен алгоритм расчета параметров течения и толщины льда на каждой итерации.
Summary
Philippov V. B. Sea gas pipe line model with surface freezing effect consideration.
Low gas tempreture in a sea pipeline can be a reason for surface freezing. A model that describes this process is given and solving procedure that uses the establishing method is described.
Литература
1. Дерцакян А. К., Курбатова Г. И., Неизвестное Я. В., Филиппов Б. В. Некоторые научно-технические проблемы освоения шельфа арктических морей России // Труды XIII сессии междунар. школы по моделям механики сплошных сред. СПб., 1996. С. 99-109.
2. Скробач А. В., Филиппов Б. В., Шевцов В. Д. Турбулентное неизотермическое движение неидеального газа в трубопроводе с песчано-зернистой шероховатостью // Физическая механика №7. СПб., 1998. С. 35-49.
3. Филиппов Б. В., Филиппов В. Б. Профиль скорости турбулентного течения сжимаемого газа в шероховатых трубах // Нелинейные проблемы механики и физики деформируемого твердого тела. СПб., 2001. С. 73-84.
4. Курбатова Г. И., Макаров М.В., Филиппов В. Б. Анализ тепловых режимов течения газа в донных трубопроводах // Вестн. С.-Петерб. ун-та. Сер. 1, 2002, вып. 1 (№1). С. 61-67.
5. Курбатова Г. И., Филиппов Б. В., Филиппов В. Б. Неизотермическое турбулентное течение сжимаемого газа // РАН Математическое моделирование. 2003. Т. 15. №3. С. 92-108.
6. Шлихтинг Г., Теория пограничного слоя. М., 1969.
7. Васильев О. Ф., Бондарев Э. А., Воеводин А. Ф., Каниболотский М. А. Неизотермическое течение газа в трубах. Новосибирск СО, 1978.
8. Крейт Ф., Блэк У. Основы теплопередачи. М., 1983.
Статья поступила в редакцию 10 июня 2003 г.
