ГОСТ 30319.3-96
МЕЖГОСУДАРСТВЕННЫЙ СТАНДАРТ
ГАЗ ПРИРОДНЫЙ
МЕТОДЫ РАСЧЕТА ФИЗИЧЕСКИХ СВОЙСТВ
ОПРЕДЕЛЕНИЕ ФИЗИЧЕСКИХ СВОЙСТВ
ПО УРАВНЕНИЮ СОСТОЯНИЯ
МЕЖГОСУДАРСТВЕННЫЙ
СОВЕТ
ПО СТАНДАРТИЗАЦИИ, МЕТРОЛОГИИ И СЕРТИФИКАЦИИ
Минск
Предисловие
1 РАЗРАБОТАН Всероссийским научно-исследовательским центром стандартизации, информации и сертификации сырья, материалов и веществ (ВНИЦ СМВ) Госстандарта России; фирмой «Газприборавтоматика» акционерного общества «Газавтоматика» РАО «Газпром»
ВНЕСЕН Госстандартом Российской Федерации
2 ПРИНЯТ Межгосударственным Советом по стандартизации, метрологии и сертификации (протокол № 9-96 от 12 апреля 1996 г.)
За принятие проголосовали:
|
Наименование государства |
Наименование национального органа по стандартизации |
|
Азербайджанская Республика |
Азгосстандарт |
|
Республика Армения |
Армгосстандарт |
|
Республика Беларусь |
Госстандарт Беларуси |
|
Республика Грузия |
Грузстандарт |
|
Республика Казахстан |
Госстандарт Республики Казахстан |
|
Киргизская Республика |
Киргизстандарт |
|
Республика Молдова |
Молдовастандарт |
|
Российская Федерация |
Госстандарт России |
|
Республика Таджикистан |
Таджикгосстандарт |
|
Туркменистан |
Главная государственная инспекция Туркменистана |
|
Украина |
Госстандарт Украины |
3 ПОСТАНОВЛЕНИЕМ Государственного комитета Российской Федерации по стандартизации, метрологии и сертификации от 30 декабря 1996 г. № 723 межгосударственный стандарт ГОСТ 30319.3-96 введен в действие непосредственно в качестве государственного стандарта Российской Федерации с 1 июля 1997 г.
4 ВВЕДЕН ВПЕРВЫЕ
5 ПЕРЕИЗДАНИЕ
СОДЕРЖАНИЕ
ГОСТ 30319.3-96
МЕЖГОСУДАРСТВЕННЫЙ СТАНДАРТ
Газ природный
МЕТОДЫ РАСЧЕТА ФИЗИЧЕСКИХ СВОЙСТВ
Определение физических свойств по уравнению состояния
Natural gas. Methods of calculation of
physical properties.
Definition of physical properties by equation of state
Дата введения 1997-07-01
1 Назначение и область применения
Настоящий стандарт предназначен для определения физических свойств природного газа. Стандарт устанавливает метод расчета плотности, показателя адиабаты, скорости звука, динамической вязкости природного газа, основанный на использовании его уравнения состояния. Метод расчета физических свойств природного газа, приведенный в настоящем стандарте, рекомендуется применять для аттестации других методов расчета.
Используемые в настоящем стандарте определения и обозначения приведены в соответствующих разделах ГОСТ 30319.0.
2 Нормативные ссылки
В настоящем стандарте использованы ссылки на следующие стандарты:
ГОСТ 30319.0-96 Газ природный. Методы расчета физических свойств. Общие положения
ГОСТ 30319.1-96 Газ природный. Методы расчета физических свойств. Определение физических свойств природного газа, его компонентов и продуктов его переработки
ГОСТ 30319.2-96 Газ природный. Методы расчета физических свойств. Определение коэффициента сжимаемости
3 Уравнение состояния природного газа
3.1 Вид уравнения состояния
Во Всероссийском научно-исследовательском центре по стандартам, информации и сертификации сырья, материалов и веществ (ВНИЦ СМВ) для расчета физических свойств природного газа разработано уравнение состояния (УС)
, (1)
где сkl - коэффициенты УС;
rп = rм/rпк - приведенная плотность;
Тп = Т/Тпк - приведенная температура;
rм - молярная плотность, кмоль/м3;
rпк и Тпк - псевдокритические параметры природного газа.
Формулы расчета коэффициентов УС и псевдокритических параметров природного газа приведены в ГОСТ 30319.2 (см. п. 3.2.5).
3.2 Пределы применения уравнения состояния и погрешности расчета свойств
Исходными данными для расчета свойств по УС (1) являются давление, температура и компонентный состав природного газа, который выражен в молярных или объемных долях компонентов.
УС (1) предназначено для работы в интервале параметров:
по давлению - до 12 МПа;
по температуре - 240-480 К;
по составу в молярных долях:
метан ³ 0,50
этан £ 0,20
пропан £ 0,05
н-бутан £ 0,03
и-бутан £ 0,03
азот £ 0,30
диоксид углерода £ 0,30
сероводород £ 0,30
остальные компоненты £ 0,01
по плотности газа при стандартных условиях - 0,66 - 1,05 кг/м3 (плотность газа при стандартных условиях рассчитывают по формуле (16) ГОСТ 30319.1);
по высшей удельной теплоте сгорания газа - 20 - 48 МДж/м3 (высшую удельную теплоту сгорания рассчитывают по 7.2 ГОСТ 30319.1, допускается рассчитывать высшую удельную теплоту сгорания по формуле (52) ГОСТ 30319.1)»;
Погрешности расчета плотности, показателя адиабаты, скорости звука по УС (1) и динамической вязкости природного газа по уравнению (15) в указанных диапазонах параметров определены в соответствии с рекомендациями работ [1-3] и с использованием данных по скорости звука [4]. Погрешности приведены в таблице 1 без учета погрешностей исходных данных.
(Измененная редакция, Изм. № 1).
4 Определение физических свойств природного газа
4.1 Определение плотности
4.1.1 Алгоритм определения плотности r из уравнения (1) при заданных давлении (р, МПа) и температуре (Т, К) приведен в ГОСТ 30319.2 (см. п. 3.2.5).
Плотность r, кг/м3, вычисляют по формуле
. (2)
Таблица 1 - Погрешности расчета свойств природного газа
|
Область параметров состояния |
Примечание |
|||
|
240 £ Т < 270К |
Т = (270 - 480) К и р < 12 МПа |
|||
|
р £ 6МПа |
6 < р £ 12 МПа |
|||
|
Плотность |
0,3 % |
0,4 % |
0,2 % |
Природный газ не содержит сероводород |
|
Показатель адиабаты |
0,9 % |
1,0 % |
0,6 % |
|
|
Скорость звука |
0,3 % |
1,0 % |
0,5 % |
|
|
Вязкость |
2,0 % |
3,0 % |
2,0 % |
|
|
Плотность |
0,6 % |
(1,0-1,5) % |
0,4 % |
Природный газ, содержащий сероводород |
|
Показатель адиабаты |
0,6 % |
1,1 % |
0,6 % |
|
|
Скорость звука |
0,3 % |
1,0 % |
0,5 % |
|
|
Вязкость |
2,0 % |
3,0 % |
2,0 % |
|
4.1.2 Если компонентный состав природного газа задан в молярных долях, молярную массу природного газа вычисляют по формуле
, (3)
где молярные массы i-го компонента природного газа (Mi) приведены в таблице 1 ГОСТ 30319.1 (см. п. 3.2.3).
4.1.3 Если компонентный состав природного газа задан в объемных долях, то молярные доли компонентов рассчитывают по формуле (12) ГОСТ 30319.1 и далее молярную массу природного газа вычисляют по 4.1.2*.
4.1.3. (Новая редакция, Изм. № 1).
4.2 Определение показателя адиабаты
Показатель адиабаты природного газа при использовании УС (1) вычисляют по формуле
, (5)
где сp и cv - изобарная и изохорная теплоемкости,
A1 - безразмерный комплекс УС (1).
Безразмерный комплекс А1 УС (1) имеет вид
А1 = 
. (6)
Изобарную и изохорную теплоемкости рассчитывают по следующим выражениям:
, (7)
, (8)
где cvom - изохорная теплоемкость природного газа в идеально газовом состоянии, а безразмерные комплексы А2 и А3 имеют вид:
, (9)
. (10)
Изохорную теплоемкость в идеально газовом состоянии вычисляют по формулам:
; (11)
. (12)
Изобарную теплоемкость (cроi) i-го компонента в идеально газовом состоянии определяют из соотношения
, (13)
где qi = T/Tni.
Температура Тni, пределы суммирования N1i и N2i, а также константы (aj)i, и (bj)i уравнения (13) для i-го компонента природного газа приведены в таблице 2.
Таблица 2 - Константы уравнения (13)
|
j |
(aj)i |
(bj)i |
|
|
Метан N1i = 10 N2i = 6 Tni = 100 К |
0 |
1,46696186 × 102 |
|
|
1 |
-6,56744186 × 101 |
-2,09233731 × 102 |
|
|
2 |
2,02698132 × 101 |
2,06925203 × 102 |
|
|
3 |
-4,20931845 × 100 |
-1,35704831 × 102 |
|
|
4 |
6,06743008 × 10-1 |
5,64368924 × 101 |
|
|
5 |
-6,12623969 × 10-2 |
-1,34496111 × 101 |
|
|
6 |
4,30969226 × 10-3 |
1,39664152 × 100 |
|
|
7 |
-2,06597572 × 10-4 |
||
|
8 |
6,42615810 × 10-6 |
||
|
9 |
-1,16805630 × 10-7 |
||
|
10 |
9,40958930 × 10-10 |
||
|
Этан N1i = 6 N2i = 5 Tni = 100 К |
0 |
6,81209760 × 101 |
|
|
1 |
-3,06340580 × 101 |
-8,74070840 × 101 |
|
|
2 |
9,52750290 × 100 |
7,84813740 × 101 |
|
|
3 |
-1,69471020 × 100 |
-4,48658590 × 101 |
|
|
4 |
1,76305850 × 10-1 |
1,46543460 × 101 |
|
|
5 |
-9,95454020 × 10-3 |
-2,05183930 × 100 |
|
|
6 |
2,35364300 × 10-4 |
||
|
Пропан N1i = 6 N2i = 4 Tni = 100 К |
0 |
-9,209726737 × 101 |
|
|
1 |
3,070930782 × 101 |
1,748671280 × 102 |
|
|
2 |
-4,924017995 × 100 |
-1,756054503 × 102 |
|
|
3 |
5,045358836 × 10-1 |
8,874920732 × 101 |
|
|
4 |
-3,140446759 × 10-2 |
-1,720610207 × 101 |
|
|
5 |
1,076680079 × 10-3 |
||
|
6 |
-1,556890669 × 10-5 |
||
|
н-Бутан N1i = 6 N2i = 5 Tni = 100 К |
0 |
-2,096096482 × 102 |
|
|
1 |
6,877783535 × 101 |
4,055272850 × 102 |
|
|
2 |
-1,228650555 × 101 |
-4,457015773 × 102 |
|
|
3 |
1,413691547 × 100 |
2,743667350 × 102 |
|
|
4 |
-1,002920638 × 10-1 |
-8,643867287 × 101 |
|
|
5 |
3,985571861 × 10-3 |
1,070428636 × 101 |
|
|
6 |
-6,786460870 × 10-5 |
||
|
и-Бутан N1i = 5 N2i = 2 Tni = 300 К |
0 |
-3,871419306 × 101 |
|
|
1 |
4,711104578 × 101 |
2,171601450 × 101 |
|
|
2 |
-1,758225423 × 101 |
-4,492603200 × 100 |
|
|
3 |
4,183494309 × 100 |
||
|
4 |
-5,520042474 × 10-1 |
||
|
5 |
3,034658409 × 10-2 |
||
|
Азот N1i = 6 N2i = 6 Tni = 100 К |
0 |
0,113129000 × 102 |
|
|
1 |
-0,215960000 × 101 |
-0,174654000 × 102 |
|
|
2 |
0,352761000 × 100 |
0,246205000 × 102 |
|
|
3 |
-0,321705000 × 10-1 |
-0,217731000 × 102 |
|
|
4 |
0,167690000 × 10-2 |
0,116418000 × 102 |
|
|
5 |
-0,467965000 × 10-4 |
-0,342122000 × 101 |
|
|
6 |
0,542603000 × 10-6 |
0,422296000 × 100 |
|
|
Диоксид углерода N1i = 6 N2i = 4 Tni = 300 К |
0 |
-9,508041394 × 10-1 |
|
|
1 |
7,008743711 × 100 |
1,087462263 × 100 |
|
|
2 |
-3,505801670 × 100 |
-7,976765747 × 10-2 |
|
|
3 |
1,096778000 × 100 |
-2,837014896 × 10-3 |
|
|
4 |
-2,016835088 × 10-1 |
1,479612229 × 10-4 |
|
|
5 |
1,971024237 × 10-2 |
||
|
6 |
-7,860765734 × 10-4 |
||
|
Сероводород N1i = 5 N2i = 5 Tni = 100 К |
0 |
3,913550000 × 100 |
|
|
1 |
-6,848510000 × 10-2 |
0,0 |
|
|
2 |
5,644240000 × 10-2 |
0,0 |
|
|
3 |
-4,837450000 × 10-3 |
1,186580000 × 100 |
|
|
4 |
1,717820000 × 10-4 |
-1,907470000 × 100 |
|
|
5 |
-2,275370000 × 10-6 |
8,285200000 × 10-1 |
(Измененная редакция, Изм. № 1).
4.3 Определение скорости звука
Скорость звука природного газа при использовании УС (1) вычисляют по формуле
, (14)
где cp, cv и А1 - соответственно изобарная, изохорная теплоемкости природного газа и безразмерный комплекс УС (1), см. (6) - (13);
М - молярная масса природного газа, см. (3) или (4).
4.4 Определение динамической вязкости
Динамическую вязкость природного газа вычисляют по формуле
(16)
, (17)
, (18)
Молярную массу природного газа (М) вычисляют по формуле (3) или (4), а формулы расчета фактора Питцера (W) и псевдокритических параметров природного газа (Тп, rп, Тпк, rпк) приведены в ГОСТ 30319.2 (см. п. 3.2.5).
5 Вычисление погрешности расчета физических свойств природного газа с учетом погрешности исходных данных
При измерении расхода и количества природного газа, транспортируемого в газопроводах, давление (р), температуру (Т) и состав (хi) измеряют с определенной погрешностью. Перечисленные параметры являются исходными данными для расчета физических свойств по УС (1) и уравнению для вязкости (15).
В соответствии с рекомендациями ИСО 5168 [5] погрешность расчета физических свойств, которая появляется в связи с погрешностью измерения исходных данных, определяют по формуле
dид = 
, (19)
где dид - погрешность расчета свойства Q, связанная с погрешностью измерения исходных данных;
dqk - погрешность измерения параметра исходных данных;
= 
; (20)
. (21)
qk - условное обозначение k-го параметра исходных данных (р, Т, хi);
`qk - среднее значение k-го параметра в определенный промежуток времени (сутки, месяц, год и т.д.);
qkмакс и qkмин - максимальное и минимальное значения k-го параметра в определенный промежуток времени;
Q - условное обозначение свойства природного газа (r, к, и, m);
Nq - количество параметров исходных данных, Nq = 2 + N (N - количество основных компонентов природного газа, которыми являются: метан, этан, пропан, бутаны, азот, диоксид углерода, сероводород).
При вычислении частных производных по формуле (20)
свойства Qqk+ и Qqk- рассчитывают при средних параметрах 
и параметрах qk+ = 
+ D
и qk- = - D,
соответственно. Рекомендуется выбирать D = 0,5×10-2 dqk.
Свойство Q (среднее значение) рассчитывают при средних параметрах `qk.
Общую погрешность расчета физических свойств определяют по формуле
, (22)
где dQ - погрешность расчета физических свойств по УС (1) и по уравнению для вязкости (15), значение которой для каждого свойства приведено в таблице 1.
(Измененная редакция, Изм. № 1).
6 Применение уравнения состояния для аттестации других методов расчета физических свойств природного газа
Приведенный в настоящем стандарте метод расчета физических свойств природного газа необходимо применять для аттестации других методов расчета. Алгоритм проведения такой аттестации состоит в следующем:
Таблица 3
|
Концентрация компонентов, мол.%, при rс, кг/м3 |
||||
|
0,67 - 0,70 |
0,70 - 0,76 |
0,76 - 0,88 |
свыше 0,88 |
|
|
Метан |
90,40 - 99,60 |
86,35 - 98,50 |
73,50 - 92,00 |
74,20 - 81,53 |
|
Этан |
0,0 - 4,10 |
0,0 - 8,40 |
1,57 - 10,91 |
6,29 - 12,19 |
|
Пропан |
0,0 - 1,16 |
0,0 - 3,35 |
0,18 - 5,00 |
3,37 - 5,00 |
|
н-Бутан |
0,0 - 0,48 |
0,0 - 1,54 |
0,12 - 1,50 |
0,51 - 1,98 |
|
н-Пентан |
0,0 - 0,32 |
0,0 - 1,00 |
0,10 - 1,00 |
0,10 - 1,00 |
|
Азот |
0,0 - 4,60 |
0,12 - 8,47 |
0,22 - 16,30 |
0,56 - 4,40 |
|
Диоксид углерода |
0,0 - 1,70 |
0,0 - 3,30 |
0,0 - 5,60 |
0,10 - 14,80 |
|
Сероводород |
0,0 |
0,0 - 6,50 |
0,0 - 5,30 |
0,0 - 24,00 |
1) используя данные, приведенные в таблице 3, подбираются 5-6 тестовых смесей природного газа таким образом, чтобы сумма молярных долей компонентов этих смесей была равна 1;
2) в заданных интервалах давления и температуры по УС (1) и уравнению для вязкости (15) насчитываются массивы физических свойств для выбранных тестовых смесей, рекомендуемое количество тестовых точек в массивах - не менее 100;
3) вычисляются систематическое и стандартное отклонения рассчитанных по аттестуемым методам физических свойств от тестовых данных, которые получены в перечислении 2) алгоритма
, (23)
, (24)
в формулах (23) и (24) N - количество тестовых точек в массивах
, (25)
где Qрасч и Qтест - условное обозначение, соответственно, расчетного по аттестуемым методам и рассчитанного в перечислении 2) алгоритма тестового значений физического свойства природного газа (r, к, и, m);
4) определяется погрешность расчета свойства Q по аттестуемым методам согласно ИСО 5168 [5]
, (26)
где dQ - погрешность расчета физических свойств по УС (1) и по уравнению для вязкости (15), значение которой для каждого свойства приведено в таблице 1.
Если для аттестуемых методов в качестве исходных данных используют плотность смеси природного газа при стандартных условиях (rc), ее значение для тестовых смесей необходимо рассчитывать по УС (1). Допускается также рассчитывать плотность rc по формуле (16) ГОСТ 30319.1 (см. п. 3.3.2).
ПРИЛОЖЕНИЕ А
Листинг программы расчета физических свойств природного газа
Расчет физических свойств природного газа по уравнению состояния (1) и по уравнению для вязкости (15) реализован на ПЭВМ, совместимых с IBM PC/AT/XT, на языке программирования ФОРТРАН-77.
C **************************************************************
C * *
С * Программа расчета физических свойств (плотности, показателя *
С * адиабаты, скорости звука и вязкости) природного газа по *
С * уравнению состояния ВНИЦ СМВ. *
С * *
C **************************************************************
IMPLICIT REAL*8(A-H,O-Z)
CHARACTER*26 AR
DIMENSION PI(100),TI(100),ROP(100,100),PAP(100,100),
*WP(100,100),ETAP(100,100)
COMMON/P/P/T/T/RON/RON/YI/YC(25)/NPR/NPR/Z/Z/TS/RO,PA,W
*/ETA/ETA/AR/AR(25)
200 WRITE(*,300)
300 FORMAT(18(/))
WRITE(*,400)
400 FORMAT(
*’ Расчет физических свойств природного газа’/
*’ по уравнению состояния’/////)
WRITE(*,1)
1 FORMAT(’ Введите исходные данные для расчета.’/)
WRITE(*,35)
35 FORMAT(’ Введите 0, если состав задан в молярных долях’/
*’ или 1, если состав задан в объемных долях ’\)
READ(*,*)NPR
IF(NPR.EQ.1)THEN
WRITE(*,’(A\)’)
*’ Плотность при 293.15 К и 101.325 кПа, в кг/куб.м’
READ(*,*)RON
WRITE(*,33)
33 FORMAT(’ Значение объемной доли, в об.%’)
ELSE
RON=0D0
WRITE(*,3)
3 FORMAT(’ Значение молярной доли, в мол.%’)
ENDIF
DO 5 I=1,25
WRITE(*,’(A\)’) AR(I)
READ(*,*)YC(I)
5 YC(I)=YC(I)/100.
WRITE(*,’(A\)’)
*’ Введите количество точек по давлению: ’
READ(*,*)NP
WRITE(*,’(A\)’)
*’ Введите количество точек по температуре: ’
READ(*,*)NT
WRITE(*,’(A\)’)
*’ Введите значения давлений в МПа: ’
READ(*,*)(PI(I),I=1,NP)
WRITE(*,’(A\)’)
*’ Введите значения температур в К: ’
READ(*,*)(TI(I),I=1,NT)
WRITE(*,’(A\)’)
*’ Ввод исходных данных завершен.’
P=.101325D0
T=293.15D0
ICALC=1
CALL EOSVNIC(ICALC)
IF(Z.EQ.0D0) THEN
CALL RANGE(NRANGE)
IF (NRANGE) 134,134,200
ENDIF
ICALC=2
NTS=0
DO 7 I=1,NP
P=PI(I)
DO 7 J=1,NT
T=TI(J)
CALL EOSVNIC(ICALC)
IF(Z.NE.0D0) THEN
NTS=NTS+1
ROP(I,J)=RO
PAP(I,J)=PA
WP(I,J)=W
ETAP(I,J)=ETA
ELSE
ROP(I,J)=0D0
PAP(I,J)=0D0
WP(I,J)=0D0
ETAP(I,J)=0D0
ENDIF
7 CONTINUE
500 WRITE(*,100)
100 FORMAT(25(/))
IF(NTS.EQ.0) THEN
CALL RANGE(NRANGE)
IF (NRANGE) 134,134,200
ELSE
1=1
9 IS=0
DO 11 J=1,NT
IF(ROP(I,J).EQ.0D0) IS=IS+1
11 CONTINUE
IF(IS.EQ.NT) THEN
IF(I.NE.NP) THEN
DO 13 J=I,NP-1
PI(J)=PI(J+1)
DO 13 K=1,NT
ROP(J,K)=ROP(J+1,K)
PAP(J,K)=PAP(J+1,K)
WP(J,K)=WP(J+1,K)
13 ETAP(J,K)=ETAP(J+1,K)
ENDIF
NP=NP-1
ELSE
I=I+1
ENDIF
IF(I.LE.NP) GO TO 9
J=l
15 JS=0
DO 17 I=1,NP
IF(ROP(I,J).EQ.0D0) JS=JS+1
17 CONTINUE
IF(JS.EQ.NP) THEN
IF(J.NE.NT) THEN
DO 19 I=J,NT-1
TI(I)=TI(I+1)
DO 19 K=1,NP
ROP(K,I)=ROP(K,I+1)
PAP(K,I)=PAP(K,I+1)
WP(K,I)=WP(K,I+1)
19 ETAP(K,I)=ETAP(K,I+1)
ENDIF
NT=NT-1
ELSE
J=J+1
ENDIF
IF(J.LE.NT) GO TO 15
CALL PROP(NPROP)
IF(NPROP.EQ.5) GO TO 134
IF(NPROP.EQ.l) CALL TABL(PI,TI,ROP,NP,NT,NPROP)
IF(NPROP.EQ.2) CALL TABL(PI,TI,PAP,NP,NT,NPROP)
IF(NPROP.EQ.3) CALL TABL(PI,TI,WP,NP,NT,NPROP)
IF(NPROP.EQ.4) CALL TABL(PI,TI,ETAP,NP,NT,NPROP)
WRITE(*,’(A\)’)
*’ Продолжить вывод рассчитанных свойств ? 0 - нет, 1 - да ’
READ(*,*)NCONT
IF(NCONT.EQ.l) GO TO 500
ENDIF
134 STOP
END
SUBROUTINE PROP(NPROP)
WRITE(*,1)
1 FORMAT(//
*10X,’¾¾¾¾ Рассчитаны следующие физические свойства ¾¾¾’/
*10Х,’ ’/
*10Х,’ 1. Плотность ’/
*10Х,’ ’/
*10Х,’ 2. Показатель адиабаты ’/
*10Х,’ ’/
*10Х,’ 3. Скорость звука ’/
*10Х,’ ’/
*10Х,’ 4. Коэффициент динамической вязкости ’/
*10Х,’ ’/
*10Х,’¾¾¾¾¾¾¾¾¾¾¾¾¾¾¾¾¾¾¾¾¾¾¾¾¾¾¾’/)
WRITE(*,5)
5 FORMAT(/,3X,
*’Введите порядковый номер свойства для вывода результатов расче’,
*’та’/
*’ или 5 для выхода в ДОС ’\)
READ(*,*)NPROP
RETURN
END
SUBROUTINE RANGE(NRANGE)
IMPLICIT REAL*8(A-H,O-Z)
COMMON/Z/Z
WRITE(*,1)
1 FORMAT(//
*’ Метод расчета при заданных параметрах “не работает”’/
*’ Продолжить работу программы ? 0 - нет, 1 - да ’\)
READ(*,*)NRANGE
RETURN
END
SUBROUTINE TABL(PI,TI,ZP,NP,NT,NPROP)
IMPLICIT REAL*8(A-H,O-Z)
CHARACTER*26 AR,FNAME
CHARACTER PROP(4)*58,A*6,LIN1(5)*9,LIN2(5)*9,LIN3(6)*9,LIN4*9,
*AT(6)*28,RAZM(4)*39
CHARACTER*70 F,FZ(11,2),FW(11,2)
DIMENSION PI(100),TI(100),ZP(100,100),ZPP(6)
COMMON/YI/YC(25)/NPR/NPR/AR/AR(25)
DATA PROP/
*’ Плотность природного газа.’,
*’ Показатель адиабаты природного газа.’,
*’ Скорость звука природного газа.’,
*’ Коэффициент динамической вязкости природного газа.’/
DATA RAZM/
*’ (в кг/куб.м)’,’ ’,
*’ (в м/с)’,
*’ (в мкПа*с)’/
DATA LIN1/5*’¾¾¾¾¾ ’/,LIN2/5*’¾¾¾¾¾’/,LIN3/6*’¾¾¾¾¾’/,
*LIN4/’¾¾¾¾¾’/,А/’ - ’/
DATA AT/
*’ Т,К’,’ Т,К’,’ Т,К’,’ Т,К’,
*’ Т,К’,’ Т,К’/
DATA FZ/
*’(3Х,F5.2,2X,6(3X,F6.2))’,’(3X,F5.2,5X,A6,5(3X,F6.2))’,
*’(3X,F5.2,2X,2(3X,A6),4(3X,F6.2))’,’(3X,F5.2,2X,3(3X,A6),
*3(3X,F6.2))’,
*’(3X,F5.2,2X,4(3X,A6),2(3X,F6.2))’,’(3X,F5.2,2X,5(3X,A6),
*3X,F6.2)’,
*’(3X,F5.2,2X,5(3X,F6.2),3X,A6)’,’(3X,F5.2,2X,4(3X,F6.2),
*2(3X,A6))’,
*’(3X,F5.2,2X,3(3X,F6.2),3(3X,A6))’,’(3X,F5.2,2X,2(3X,F6.2),
*4(3X,A6))’,
*’(3X,F5.2,5X,F6.2,5(3X,A6))’,’(3X,F9.6,1X,F6.2,5(3X,F6.2))’,
*’(3X,F9.6,1X,A6,5(3X,F6.2))’,’(3X,F9.6,1X,A6,3X,A6,4(3X,F6.2))’,
*’(3X,F9.6,1X,A6,2(3X,A6),3(3X,F6.2))’,’(3X,F9.6,1X,A6,3(3X,A6),
*2(3X,F6.2))’,
*’(3X,F9.6,1X,A6,4(3X,A6),3X,F6.2)’,’(3X,F9.6,1X,F6.2,4(3X,F6.2),
*3X,A6)’,
*’(3X,F9.6,1X,F6.2,3(3X,F6.2),2(3X,A6))’,’(3X,F9.6,1X,F6.2,
*2(3X,F6.2),3(3X,A6))’,
*’(3X,F9.6,1X,F6.2,3X,F6.2,4(3X,A6))’,’(3X,F9.6,1X,F6.2,5(3X,A6))’/
DATA FW/
’(3X,F5.2,2X,6(4X,F5.1))’,’(3X,F5.2,5X,A6,5(4X,F5.1))’,
*’(3X,F5.2,2X,2(3X,A6),4(4X,F5.1))’,’(3X,F5.2,2X,3(3X,A6),
*3(4X,F5.1))’,
*’(3X,F5.2,2X,4(3X,A6),2(4X,F5.1))’,’(3X,F5.2,2X,5(3X,A6),
*4X,F5.1)’,
*’(3X,F5.2,2X,5(4X,F5.1),3X,A6)’,’(3X,F5.2,2X,4(4X,F5.1),
*2(3Х,А6))’,
*’(3X,F5.2,2X,3(4X,F5.1),3(3X,A6))’,’(3X,F5.2,2X,2(4X,F5.1),
*4(3X,A6))’,
*’(3X,F5.2,6X,F5.1,5(3X,A6))’,’(3X,F9.6,2X,F5.1,5(4X,F5.1))’,
*’(3Х,F9.6,1X,A6,5(4X,F5.1))’,’(3X,F9.6,1X,A6,3X,A6,4(4X,F5.1))’,
*’(3X,F9.6,1X,A6,2(3X,A6),3(4X,F5.1))’,’(3X,F9.6,1X,A6,3(3X,A6),
*2(4X,F5.1))’,
*’(3X,F9.6,1X,A6,4(3X,A6),4X,F5.1)’,’(3X,F9.6,2X,F5.1,4(4X,F5.1),
*3X,A6)’,
*’(3X,F9.6,2X,F5.1,3(4X,F5.1),2(3X,A6))’,’(3X,F9.6,2X,F5.1,
*2(4X,F5.1),3(3X,A6))’,
*’(3X,F9.6,2X,F5.1,4X,F5.1,4(3X,A6))’,’(3X,F9.6,2X,F5.1,5(3X,A6))’/
22 WRITE(*,44)
44 FORMAT(//’ Устройство вывода результатов расчета ?,’)
WRITE(*,’(A\)’)
*’ 0 - дисплей, 1 - принтер, 2 - файл на диске ’
READ(*,*)NYST
IF(NYST.EQ.0) OPEN(1,FILE=’CON’)
IF(NYST.EQ.l) OPEN(1,FILE=’PRN’)
IF(NYST.EQ.2) WRITE(*,’(A\)’) ’ Введите имя файла ’
IF(NYST.EQ.2) READ(*,’(A)’)FNAME
IF(NYST.EQ.2) OPEN(1,FILE=FNAME)
IF(NYST.EQ.0) WRITE(*,100)
100 FORMAT(25(/))
IF(NYST.EQ.l) PAUSE
*’ Включите принтер, вставьте бумагу и нажмите <ВВОД> ’
WRITE(1,88)PROP(NPROP),RAZM(NPROP)
88 FORMAT(A58/A39/)
NW=3
IF(NPR.EQ.0) WRITE(1,3)
3 FORMAT(’ Содержание в мол.%’)
IF(NPR.EQ.l) WR1TE(1,33)
33 FORMAT(’ Содержание в об.%’)
NW=NW+1
I=1
9 J=I+1
13 CONTINUE
IF(YC(J).NE.0D0) THEN
WRITE(1,5)AR(I),YC(I)*100.,AR(J),YC(J)*100.
5 FORMAT(2(A26,F7.4))
NW=NW+1
DO 11 I=J+1,25
IF(YC(I).NE.0D0.AND.I.NE.25) GO TO 9
IF(YC(I).NE.0D0.AND.I.EQ.25)THEN
WRITE(1,5)AR(I),YC(I)*100.
NW=NW+1
GO TO 99
ENDIF
11 CONTINUE
ELSE
J=J+1
IF(J.LE.25) THEN
GO TO 13
ELSE
WRITE(1,5)AR(I),YC(I)*100.
NW=NW+1
ENDIF
ENDIF
99 CONTINUE
IF(NW.GT.12.AND.NYST.EQ.0) THEN
WRITE(*,7)
7 FORMAT(/)
PAUSE ’ Для продолжения вывода нажмите <ВВОД> ’
WRITE(*,100)
NW=0
ENDIF
DO 15 I=1,NT,6
IF(NW.GT.12.AND.NYST.EQ.0) THEN
WRITE(*,7)
PAUSE ’ Для продолжения вывода нажмите <ВВОД> ’
WRITE(*,100)
NW=0
ENDIF
IF(NW.GT.46 AND.NYST.NE.0) THEN
WRITE(1,7)
WRITE(*,7)
IF(NYST.EQ.l) PAUSE
*’ Для продолжения вывода вставьте бумагу и нажмите <ВВОД> ’
NW=0
ENDIF
IF(I+5.LE.NT) THEN
NL=6
ELSE
NL=NT-I+1
ENDIF
WRITE(1,7)
IF(NL.GT.1)WRITE(1,17)LIN2(1),(LIN1(K),K=1,NL-1)
IF(NL.EQ.l) WRITE(1,17)LIN2(1)
17 FORMAT(’¾¾¾’,6A9)
WRITE(1,19)AT(NL)
19 FORMAT(’ ú’,A28)
IF(NL.GT.1)WRITE(1,21)LIN4,(LIN2(K),K=1,NL-1)
IF(NL.EQ.l) WRITE(1,21)LIN4
21 FORMAT(’ p, МПа ’,6А9)
WRITE(1,23)(TI(K),K=I,I+NL-1)
23 FORMAT(10X,6(:,’ ’,F6.2))
WRITE(1,17)(LIN3(K),K=1,NL)
NW=NW+6
DO 25 J=1,NP
JP=1
IF(PI(J).EQ.0.101325D0) JP=2
NL1=0
NLN=0
DO 27 K=I,I+NL-1
NL1=NL1+1
IF(ZP(J,K).EQ.0D0) THEN
ZPP(NL1)=A
NLN=NLN+1
ELSE
ZPP(NL1)=ZP(J,K)
ENDIF
27 CONTINUE
IF(NLN.EQ.NL) GO TO 133
IF(NLN.EQ.0) THEN
IF(NPROP.NE.3) F=FZ(1,JP)
IF(NPROP.EQ.3) F=FW(1,JP)
ELSE
IF(ZP(J,I).EQ.0D0.AND.NPROP.NE.3) F=FZ(NLN+1,JP)
IF(ZP(J,I+NL-1).EQ.0D0.AND.NPROP.NE.3) F=FZ(NLN+12-NL,JP)
IF(ZP(J,I).EQ.0D0.AND.NPROP.EQ.3)F=FW(NLN+1,JP)
IF(ZP(J,I+NL-1).EQ.0D0.AND.NPROP.EQ.3) F=FW(NLN+12-NL,JP)
ENDIF
IF(NL1.EQ.1) WRITE(1,F)PI(J),ZPP(1)
IF(NL1.EQ.2) WRITE(1,F)PI(J),ZPP(1),ZPP(2)
IF(NL1.EQ.3) WRITE(1,F)PI(J),ZPP(1),ZPP(2),ZPP(3)
IF(NL1.EQ.4) WRITE(1,F)PI(J),ZPP(1),ZPP(2),ZPP(3),ZPP(4)
IF(NL1.EQ.5)
*WRITE(1,F)PI(J),ZPP(1),ZPP(2),ZPP(3),ZPP(4),ZPP(5)
IF(NL1.EQ.6)
*WRITE(1,F)PI(J),ZPP(1),ZPP(2),ZPP(3),ZPP(4),ZPP(5),ZPP(6)
NW=NW+1
133 CONTINUE
IF(NW.EQ.20.AND.NYST.EQ.0) THEN
IF(J.EQ.NP.AND.I+NL-1.EQ.NT) GO TO 29
WRITE(*,7)
PAUSE ’ Для продолжения вывода нажмите <ВВОД> ’
WRITE(*,100)
NW=0
WRITE(1,7)
IF(NL.GT.1) WRITE(1,17)LIN2(1),(LIN1(K),K=1,NL-1)
IF(NL.EQ.l) WRITE(1,17)LIN2(1)
WRITE(1,19)AT(NL)
IF(NL.GT.1) WRITE(1,21)LIN4,(LIN2(K),K=1,NL-1)
IF(NL.EQ.l) WRITE(1,21)LIN4
WRITE(1,23)(TI(K),K=I,I+NL-1)
WRITE(1,17)(LIN3(K),K=1,NL)
NW=NW+6
ENDIF
IF(NW.EQ.54.AND.NYST.NE.0) THEN
IF(J.EQ.NP.AND.I+NL-1.EQ.NT) GO TO 29
WRITE(1,7)
WRITE(*,7)
IF(NYST.EQ.l) PAUSE
*’ Для продолжения вывода вставьте бумагу и нажмите <ВВОД> ’
NW=0
IF(NL.GT.1) WRITE(1,17)LIN2(1),(LIN1(K),K=1,NL-1)
IF(NL.EQ.l) WRITE(1,17)LIN2(1)
WRITE(1,19)AT(NL)
IF(NL.GT.1) WRITE(1,21)LIN4,(LIN2(K),K=1,NL-1)
IF(NL.EQ.l) WRITE(1,21)LIN4
WRITE(1,23)(TI(K),K=I,I+NL-1)
WRITE(1,17)(LIN3(K),K=1,NL)
NW=NW+6
ENDIF
25 CONTINUE
15 CONTINUE
29 CLOSE(1)
WRITE(*,7)
PAUSE ’ Вывод завершен, для продолжения работы нажмите <ВВОД>’
WR1TE(*,66)
66 FORMAT(/’ Назначить другое устройство вывода ?’,
*’, 0 - нет, 1 - да ’\)
READ(*,*)NBOLB
IF(NBOLB.EQ.l) GO TO 22
RETURN
END
SUBROUTINE EOSVNIC(ICALC)
IMPLICIT REAL*8(A-H,O-Z)
REAL*8 LIJ(8,8)
DIMENSION VC(8),TC(8),PII(8),DIJ(8,8)
COMMON/PARCD/VCD(8),TCD(8),PIID(8)/ABIJ/AIJ(10,8),BIJ(10,8)
*/B/B(10,8)/RM/RM/Y/Y(8)/BM/BM(8)/NI/NI(8)/NC/NC/RON/RON/PIM/PIM
COMMON/CPCI/CPC1(20,5),CPC2(20,3)/IDGFD/TOID(8),MCOD(8),MCPD(8)
*/IDGF/CPC(20,8),TOI(8),MCO(8),MCP(8)
COMMON/P/P/T/T/Z/Z/TS/RO,PA,W/ETA/ETA
RM=8.31451D0
IF(ICALC.NE.1)GOTO 1
CALL COMPON
IF(Z.EQ.0D0) GO TO 133
DO 11111 J=l,8
DO 11111 I=1,20
IF(J.LE.5) CPC(I,J)=CPC1(I,J)
IF(J.GT.5) CPC(I,J)=CPC2(I,J-5)
11111 CONTINUE
CALL DDIJ(DIJ,LIJ)
DO 75 I=1,NC
TC(I)=TCD(NI(I))
VC(I)=BM(I)/VCD(NI(I))
PII(I)=PIID(NI(I))
MCO(I)=MCOD(NI(I))
MCP(I)=MCPD(NI(I))
TOI(I)=TOID(NI(I))
MP=MCO(I)+MCP(I)+1
DO 23 J=1,MP
23 CPC(J,I)=:CPC(J,NI(I))
DO 123 J=1,NC
IF(I.GE.J) GO TO 123
DIJ(I,J)=DIJ(NI(I),NI(J))
LIJ(I,J)=LIJ(NI(I),NI(J))
123 CONTINUE
75 CONTINUE
CALL PARMIX(DIJ,LIJ,TC,VC,PII)
DO 27 I=1,10
DO 27 J= 1,8
27 B(I,J)=AIJ(I,J)+BIJ(I,J)*PIM
IF(RON.NE.0D0) THEN
CALL PHASE
RON=0D0
GO TO 133
ENDIF
1 CALL PHASE
133 RETURN
END
SUBROUTINE COMPON
IMPLICIT REAL*8(A-H,O-Z)
DIMENSION BMI(25),ROI(8),GI(8),YI(25)
COMMON/Y/Y(8)/BMM/BMM/BM/BM(8)YI/YC(25)/NI/NI(8)/NC/NC/RON/RON
DATA BMI/16.043D0,30.07DO0,44.097D0,2*58.123D0,28.0135D0,
*44.01D0,34.082D0,26.038D0,28.054D0,42.081D0,3*72.15D0,
*86.177D0,78.114D0,100.204D0,92.141D0,114.231D0,128.259D0,
*142.286D0,4.0026D0,2.0159D0,28.01D0,31.9988D0/
DATA ROI/0.6682D0,1.2601D0,1.8641D0,2.4956D0,2.488D0,
*1.1649D0,1.8393D0,1.4311D0/
DO 100 I=1,25
100 YI(I)=YC(I)
IF(RON.NE.0D0) GO TO 333
BMM=0D0
DO 3333 I=1,25
3333 BMM=BMM+YI(1)*BMI(I)
333 YS=0D0
DO 55 I=9,25
55 YS=YS+YI(I)
YS1=0D0
DO 67 I=12,21
67 YS1=YS1+YI(I)
YS2=0DO0
DO 69 I=22,25
69 YS2=YS2+YI(I)
YI(2)=YI(2)+YI(9)+YI(10)
YI(3)=YI(3)+YI(11)
YI(4)=YI(4)+YS1
YS3=YI(4)+YI(5)
IF(RON.NE.0D0.AND.YI(5).LT.0.01D0.AND.YS3.LT.0.03D0) THEN
YI(4)=YS3
YI(5)=0D0
ENDIF
IF(RON.EQ.0D0.AND.YI(5).LT.0.01D0.AND.YS3.LE.0.03D0) THEN
YI(4)=YS3
YI(5)=0D0
ENDIF
YI(6)=YI(6)+YS2
IF(RON.EQ.0D0) GO TO 555
ROM=0D0
DO 7 I=1,8
7 ROM=ROM+YI(I)*ROI(I)
DO 9 I=1,8
9 GI(I)=YI(I)*ROI(I)/ROM
SUM=0D0
DO 11 I=1,8
11 SUM=SUM+GI(I)/BMI(I)
SUM=1./SUM
DO 13 I=1,8
13 YI(I)=GI(I)*SUM/BMI(I)
555 NC=0
YSUM=0D0
DO 155 I=1,8
IF(YI(I).EQ.0D0) GO TO 155
NC=NC+1
NI(NC)=I
Y(NC)=YI(I)
YSUM=YSUM+Y(NC)
BM(NC)=BMI(I)
155 CONTINUE
CALL MOLDOL(YI.YS)
DO 551 I=1,NC
551 Y(I)=Y(I)/YSUM
RETURN
END
SUBROUTINE MOLDOL(YI,YS)
IMPLICIT REAL*8(A-H,O-Z)
DIMENSION YI(25)
COMMON/Z/Z
Z=-1D0
IF(YI(1).LT.0.5D0.OR.YI(2).GT.0.2D0.OR.YI(3).GT.0.05D0.OR.
*YI(4).GT.0.03D0.OR.YI(5).GT,0.03D0.OR.YS.GT.0.01D0) Z=0D0
IF(YI(6).GT.0.3D0.OR.YI(7).GT.0.3D0.OR.YI(8).GT.0.3D0) Z=0D0
RETURN
END
SUBROUTINE DDIJ(DIJ,LIJ)
IMPLICIT REAL*8(A-H,O-Z)
REAL*8 LIJ(8,8)
DIMENSION DU(8,8)
DO 1 I=1,8
DO 1 J=l,8
LIJ(I,J)=0.D0
1 DIJ(I,J)=0.D0
DIJ(1,2)=0.036D0
DIJ(1,3)=0.076D0
DIJ(1,4)=0.121D0
DIJ(1,5)=0.129D0
DIJ(1,6)=0.06D0
DIJ(1,7)=0.074D0
DU(2,6)=0.106D0
DIJ(2,7)=0.093D0
DU(6,7)=0.022D0
DIJ(1,8)=0.089D0
DIJ(2,8)=0.079D0
DIJ(6,8)=0.211D0
DIJ(7,8)=0.089D0
LIJ(1,2)=-0.074D0
LIJ(1,3)=-0.146D0
LIJ(1,4)=-0.258D0
LIJ(1,5)=-0.222D0
LIJ(1,6)=-0.023D0
LIJ(1,7)=-0.086D0
LIJ(6,7)=-0.064D0
LIJ(7,8)=-0.062D0
RETURN
END
SUBROUTINE PARMIX(DIJ,LIJ,TC,VC,PII)
IMPLICIT REAL*8(A-H,O-Z)
REAL*8 L1J(8,8)
DIMENSION Y(8),DIJ(8,8),VCIJ(8,8),TCIJ(8,8),V13(8),TC(8),VC(8),
*PII(8),PIIJ(8,8)
COMMON/PARCM/TCM,VCM/Y/Y/NC/NC/PCM/PCM/PIM/PIM
DO 1 I=1,NC
1 V13(I)=VC(I)**(1.D0/3.D0)
DO 3 I=1,NC
VCIJ(I,I)=VC(I)
PIIJ(I,I)=PII(I)
TCIJ(I,I)=TC(I)
DO 3 J=1,NC
IF(I.GE.J) GO TO 3
VCIJ(I,J)=(1.D0-LIJ(I,J))*((V13(I)+V13(J))/2.)**3
PIIJ(I,J)=(VC(I)*PII(I)+VC(J)*PII(J))/(VC(I)+VC(J))
TCIJ(I,J)=(1.D0-DIJ(I,J))*(TC(I)*TC(J))**0.5
VCIJ(J,I)=VCIJ(I,J)
PIIJ(J,I)=PIIJ(I,J)
TCIJ(J,I)=TCIJ(I,J)
3 CONTINUE
VCM=0.D0
PIM=0.D0
TCM=0.D0
DO 5 I=1,NC
DO 5 J=1,NC
VCM=VCM+Y(I)*Y(J)*VCIJ(I,J)
PIM=PIM+Y(I)*Y(J)*VCIJ(I,J)*PIIJ(I,J)
5 TCM=TCM+Y(I)*Y(J)*VCIJ(I,J)*TCIJ(I,J)**2
PIM=PIM/VCM
TCM=(TCM/VCM)**0.5
PCM=8.31451D-3*(0.28707D0-0.05559*PIM)*TCM/VCM
RETURN
END
SUBROUTINE PHASE
IMPLICIT REAL*8(A-H,O-Z)
COMMON/Z/Z/RM/RM/T/T/P/P/PCM/PCM/RON/RON/BMM/BMM
*/AI/AO,A1,A2,A3
IF(T.LT.240D0.OR.T.GT.480D0.OR.P.LE.0D0.OR.P.GT.12D0) THEN
Z=0D0
GO TO 134
ENDIF
PR=P/PCM
RO=9D3*P/(RM*T*(1.1*PR+0.7D0))
CALL FUN(RO)
CALL OMTAU(RO,T)
IF(Z.EQ.0D0) GO TO 134
Z=1.D0+AO
IF(RON.NE.0D0) THEN
BMM=1D-3*Z*RON*RM*T/P
GO TO 134
ENDIF
NPRIZ=2
CALL COMPL(RO,T,NPRIZ)
CALL TP(RO)
CALL ETAS(RO)
134 RETURN
END
С Подпрограмма, реализующая итерационный процесс определения
С плотности из уравнения состояния (метод Ньютона)
SUBROUTINE FUN(X)
IMPLICIT REAL*8(A-H,O-Z)
COMMON/P/P/RM/RM/T/T/AI/AO,A1,A2,A3
ITER=1
1 CONTINUE
NPRIZ=0
IF(ITER.NE.l) NPRIZ=1
CALL COMPL(X,T,NPRIZ)
Z=1.D0+AO
FX=1.D6*(P-(1.D-3*RM*T*Z*X))
F=1.D3*RM*T*(1.D0+A1)
DR=FX/F
X=X+DR
IF(ITER.GT.10) GO TO 4
ITER=ITER+1
IF(DABS(DR/X).GT.1.D-6) GO TO 1
4 CALL COMPL(X,T,NPRIZ)
RETURN
END
SUBROUTINE OMTAU(RO,T)
IMPLICIT REAL*8(A-H,O-Z)
COMMON/PARCM/TCM,VCM/Z/Z
Z=-1D0
TR=T/TCM
ROR=RO*VCM
IF(TR.LT.1.05D0) Z=0D0
IF(ROR.LT.0.D0.OR.ROR.GT.3.D0) Z=0D0
RETURN
END
С Подпрограмма определения безразмерных комплексов АО,А1,А2 и A3
SUBROUTINE COMPL(RO,T,NPRIZ)
IMPLICIT REAL*8(A-H,O-Z)
DIMENSION B(10,8),BK(10)
COMMON/PARCM/TCM,VCM/B/B/AI/AO,A1,A2,A3
IF(NPRIZ.NE.0) GO TO 7
TR=T/TCM
DO 1 I=1,10
BK(I)=0
DO 1 J=l,8
1 BK(I)=BK(I)+B(I,J)/TR**(J-1)
7 ROR=RO*VCM
AO=0.D0
A1=0.D0
IF(NPRIZ.EQ.l) GO TO 5
A2=0.D0
A3=0.D0
5 DO 33 I=1,10
D=BK(I)*ROR**I
AO=AO+D
A1=A1+(I+1)*D
IF(NPRIZ.EQ.1) GO TO 33
DO 3 J=1,8
D1=B(I,J)*ROR**I/TR**(J-1)
A2=A2+(2-J)*D1
3 A3=A3+(J-1)*(2-J)*D1/I
33 CONTINUE
RETURN
END
С Подпрограмма расчета плотности, показателя адиабаты, скорости
С звука
SUBROUTINE TP(ROM)
IMPLICIT REAL*8(A-H,O-Z)
COMMON/BMM/BMM/AI/AO,A1,A2,A3/RM/RM/T/T/TS/RO,PA,W/Z/Z
CALL IDGFU(T,CVOS)
RO=BMM*ROM
R=RM/BMM
A11=1.D0+A1
A21=1.D0+A2
CV=R*(A3+CVOS)
CP=CV+R*A21**2/A11
W=DSQRT(DABS(1.DЗ*R*T*CP/CV))*DSQRT(DABS(A11))
PA=CP/CV*A11/Z
RETURN
END
С Подпрограмма расчета изохорной теплоемкости в идеально газовом
С состоянии
SUBROUTINE IDGFU(T,CVOS)
IMPLICIT REAL*8(A-H,O-Z)
DIMENSION CPO(8),CVO(8)
COMMON/IDGF/CPC(20,8),TOI(8),MCO(8),MCP(8)/Y/Y(8)/NC/NC
CVOS=0.D0
DO 21 I=1,NC
M=MCP(I)
N=MCO(I)
TAU=T/TOI(I)
S1=0.D0
S2=0.D0
S3=0.D0
S1=CPC(1,I)
IF(M.EQ.0) GO TO 7
DO 9 J=1,M
9 S2=S2+СРС(J+1,I)*ТАU**J
7 IF(N.EQ.0) GO TO 11
DO 13 J=1,N
13 S3=S3+CPC(M+J+1,I)/TAU**J
11 CPO(I)=S1+S2+S3
CVO(I)=CPO(I)-1.D0
21 CVOS=CVOS+Y(I)*CVO(I)
RETURN
END
С Подпрограмма расчета вязкости
SUBROUTINE ETAS(ROM)
IMPLICIT REAL*8(A-H,O-Z)
COMMON/ETA/ETA/PARCM/TCM,VCM/BMM/BMM/T/T/PIM/PIM/PCM/PCM
DKSI=TCM*(1D0/6D0)/BMM**.5/PCM**(2D0/3D0)
ROR=VCM*ROM
TR=T/TCM
ETA=78.037D0+3.85612*PIM-29.0053*PIM**2-156.728/TR+145.519/TR**2
*-51.1082/TR**3+6.57895*ROR+(11.7452D0-95.7215*PIM**2/TR)*ROR**2+
*17.1027*ROR**3*PIM+.519623/TR**2*ROR**5
ETA=ETA/DKSI/10.
RETURN
END
BLOCK DATA BDVNIC
IMPLICIT REAL*8(A-H,O-Z)
CHARACTER*26 AR
COMMON/PARCD/VCD(8),TCD(8),PIID(8)/ABIJ/AIJ(10,8),BIJ(10,8)
COMMON/CPCI/CPC1(20,5),CPC2(20,3)/IDGFD/TOID(8),MCOD(8),MCPD(8)
*/AR/AR(25)
DATA TCD/190.67D0,305.57D0,369.96D0,425.4D0,407.96D0,
*125.65D0,304.11D0,373.18D0/
DATA VCD/163.03D0,205.53D0,218.54D0,226.69D0,225.64D0,
*315.36D0,466.74D0,349.37D0/
DATA PIID/0.0006467D0,0.1103D0,0.1764D0,0.2213D0,0.2162D0,
*0.04185D0,0.2203D0,0.042686D0/
DATA AIJ/.6087766D0,-.4596885D0,1.14934D0,-.607501D0,
*-.894094D0,1.144404D0,-.34579D0,-.1235682D0,.1098875D0,
*-.219306D-1,-1.832916D0,4.175759D0,-9.404549D0,10.62713D0,
*-3.080591D0,-2.122525D0,1.781466D0,-.4303578D0,-.4963321D-1,
*.347496D-1,1.317145D0,-10.73657D0,23.95808D0,-31.47929D0,
*18.42846D0,-4.092685D0,-.1906595D0,.4015072D0,-.1016264D0,
*-.9129047D-2,-2.837908D0,15.34274D0,-27.71885D0,35.11413D0,
*-23.485D0,7.767802D0,-1.677977D0,.3157961D0,.4008579D-2,0.D0,
*2.606878D0,-11.0б722D0,12.79987D0,-12.11554D0,7.580666D0,
*-1.894086D0,4*0.D0,
*-1.15575D0,3.601316D0,-.7326041D0,-1.151685D0,.5403439D0,
*5*0.D0,.9060572D-1,-.5151915D0,.7622076D-1,7*0.D0,
*.4507142D-1,9*0.D0/
DATA BIJ/-.7187864D0,10.67179D0,-25.7687D0,17.13395D0,
*16.17303D0,-24.38953D0,7.156029D0,3.350294D0,-2.806204D0,
*.5728541D0,6.057018D0,-79.47685D0,216.7887D0,-244.732D0,
*78.04753D0,48.70601D0,-41.92715D0,10.00706D0,1.237872D0,
*-.8610273D0,-12.95347D0,220.839D0,-586.4596D0,744.4021D0,
*-447.0704D0,99.6537D0,5.136013D0,-9.57б9D0,2.41965D0,
*.2275036D0,15.71955D0,-302.0599D0,684.5968D0,-828.1484D0,
*560.0892D0,-185.9581D0,39.91057D0,-7.567516D0,-.1062596D0,
*0.D0,-13.75957D0,205.541D0,-325.2751D0,284.6518D0,
*-180.8168D0,46.05637D0,4*0.D0,
*6.466081D0,-57.3922D0,36.94793D0,20.77675D0,-12.56783D0,
*5*0.D0,-.9775244D0,2.612338D0,-.4059629D0,7*0.D0,
*-.2298833D0,9*0.D0/
DATA CPC1/1.46696186D+02,-6.56744186D+01,2.02698132D+01,
*-4.20931845D0,6.06743008D-01,-6.12623969D-02,4.30969226D-03,
*-2.06597572D-04,6.4261581D-06,-1.1680563D-07,9.4095893D-10,
*-2.09233731D+02,2.06925203D+02,-1.35704831D+02,5.64368924D+01,
*-1.34496111D+01,1.39664152D0,3*0.D0,
*6.8120976D+01,-3.0634058D+01,9.5275029D0,-1.6947102D0,
*1.7630585D-01,-9.9545402D-3,2.353643D-4,-8.7407084D+1,
*7.8481374D+1,-4.4865859D+1,1.4654346D+1,-2.0518393D0,8*0.D0,
*-9.209726737D+1,3.070930782D+1,-4.924017995D0,5.045358836D-1,
*-3.140446759D-2,1.076680079D-3,-1.556890669D-5,1.74867128D+2,
*-1.756054503D+2,8.874920732D+1,-1.720610207D+ 1,9*0.D0,
*-2.096096482D+2,6.877783535D+1,-1.228650555D+1,1.413691547D0,
*-1.002920638D-1,3.985571861D-3,-6.78646087D-5,4.05527285D+2,
*-4.457015773D+2,2.74366735D+2,-8.643867287D+1,1.070428636D+1,
*8*0.D0,
*-3.871419306D+1,4.711104578D+1,-1.758225423D+1,4.183494309D0,
*-5.520042474D-1,3.034658409D-2,2.17160145D+1,-4.4926032D0,
* 12*0.D0/
DATA CPC2/0.113129D+2,-0.21596D+1,0.352761D0,-0.321705D-1,
*0.16769D-2,-0.467965D-4,0.542603D-6,-0.174654D+2,0.246205D+2,
*-0.217731D+2,0.116418D+2,-0.342122D+1,0.422296D0,7*0.D0,
*-9.508041394D-1,7.008743711D0,-3.50580167D0,1.096778D0,
*-2.016835088D-1,1.971024237D-2,-7.860765734D-4,1.087462263D0,
*-7.976765747D-2,-2.837014896D-3,1.479612229D-4,9*0.D0,
*3.91355D0,-6.84851D-2,5.64424D-2,-4.83745D-3,1.71782D-4,
*-2.27537D-6,2*0.D0,1.18658D0,-1,90747D0,8.2852D-1,9*0.D0/
DATA MCOD/6,5,4,5,2,6,4,5/
DATA MCPD/10,6,6,6,5,6,6,5/
DATA TOID/4*100D0,300D0,100D0,300D0,100D0/
DATA AR/’ метана (СН4)’,’ этана (С2Н6)’,’ пропана (С3Н8)’,
*’ н-бутана (н-С4Н10)’,’ и-бутана (и-С4Н10)’,’ азота (N2)’,
*’ диоксида углерода (СO2)’,’ сероводорода (H2S)’,
*’ ацетилена (С2Н2)’,’ этилена (С2Н4)’,’ пропилена (С3Н6)’,
*’ н-пентана (н-С5Н12)’,’ и-пентана (и-С5Н12)’,
*’ нео-пентана (нео-С5Н12)’,’ н-гексана (н-С6Н14)’,
*’ бензола (С6Н6)’,’ н-гептана (н-С7Н16)’,’ толуола (С7Н8)’,
*’ н-октана (н-С8Н18)’,’ н-нонана (Н-С9Н20)’,
*’ н-декана (н-С10Н22)’,’ гелия (Не)’,’ водорода (Н2)’,
*’ моноксида углерода (СО)’,’ кислорода (O2)’/
END
ПРИЛОЖЕНИЕ Б
Пример расчета физических свойств природного газа
Состав природного газа в молярных процентах:
метан ................................................................................................... 89,27
этан ...................................................................................................... 2,26
пропан ................................................................................................. 1,06
и-бутан ................................................................................................ 0,01
азот ...................................................................................................... 0,04
диоксид углерода ............................................................................... 4,30
сероводород ........................................................................................ 3,05
пропилен ............................................................................................ 0,01
Давление ............................................................................................. 1,081 МПа
Температура ....................................................................................... 323,15 К
Плотность ........................................................................................... 7,54 кг/м3
Показатель адиабаты ......................................................................... 1,29
Скорость звука ................................................................................... 429,8 м/с
Динамическая вязкость ..................................................................... 12,36 мкПа × с
Давление ............................................................................................. 9,950 МПа
Температура ....................................................................................... 323,15 К
Плотность ........................................................................................... 78,51 кг/м3
Показатель адиабаты ......................................................................... 1,44
Скорость звука ................................................................................... 427,7 м/с
Динамическая вязкость ..................................................................... 14,75 мкПа × с
ПРИЛОЖЕНИЕ В
Библиография
[2] Козлов А.Д., Кузнецов В.М., Мамонов Ю.В. Анализ современных методов расчета рекомендуемых справочных данных о коэффициентах вязкости и теплопроводности газов и жидкостей. - М.: ИВТАН СССР, 1989, № 3, с. 3-80
[3] МР 67-89. Расчет плотности, изобарной и изохорной теплоемкости, энтальпии, энтропии, скорости звука жидких и газообразных веществ, применяемых в криогенном машиностроении в интервале температур до 500 К и давлений до 50 МПа на основе уравнения Старлинга-Хана. - Методика ГСССД, Деп. ВНИИКИ, № 609, 1990
Ключевые слова: природный газ, методы расчета физических свойств, давление, температура, компонентный состав, молярные и объемные доли, плотность, показатель адиабаты, скорость звука, динамическая вязкость, погрешность, уравнение состояния, листинг программы.