ГОСТ 30319.3-96

МЕЖГОСУДАРСТВЕННЫЙ СТАНДАРТ

ГАЗ ПРИРОДНЫЙ

МЕТОДЫ РАСЧЕТА ФИЗИЧЕСКИХ СВОЙСТВ

ОПРЕДЕЛЕНИЕ ФИЗИЧЕСКИХ СВОЙСТВ
ПО УРАВНЕНИЮ СОСТОЯНИЯ

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

Минск

Предисловие

1 РАЗРАБОТАН Всероссийским научно-исследовательским центром стандартизации, информации и сертификации сырья, материалов и веществ (ВНИЦ СМВ) Госстандарта России; фирмой «Газприборавтоматика» акционерного общества «Газавтоматика» РАО «Газпром»

ВНЕСЕН Госстандартом Российской Федерации

2 ПРИНЯТ Межгосударственным Советом по стандартизации, метрологии и сертификации (протокол № 9-96 от 12 апреля 1996 г.)

За принятие проголосовали:

Наименование государства	Наименование национального органа по стандартизации
Азербайджанская Республика	Азгосстандарт
Республика Армения	Армгосстандарт
Республика Беларусь	Госстандарт Беларуси
Республика Грузия	Грузстандарт
Республика Казахстан	Госстандарт Республики Казахстан
Киргизская Республика	Киргизстандарт
Республика Молдова	Молдовастандарт
Российская Федерация	Госстандарт России
Республика Таджикистан	Таджикгосстандарт
Туркменистан	Главная государственная инспекция Туркменистана
Украина	Госстандарт Украины

3 ПОСТАНОВЛЕНИЕМ Государственного комитета Российской Федерации по стандартизации, метрологии и сертификации от 30 декабря 1996 г. № 723 межгосударственный стандарт ГОСТ 30319.3-96 введен в действие непосредственно в качестве государственного стандарта Российской Федерации с 1 июля 1997 г.

4 ВВЕДЕН ВПЕРВЫЕ

5 ПЕРЕИЗДАНИЕ

СОДЕРЖАНИЕ

1 Назначение и область применения. 2

2 Нормативные ссылки. 2

3 Уравнение состояния природного газа. 2

3.1 Вид уравнения состояния. 2

3.2 Пределы применения уравнения состояния и погрешности расчета свойств. 3

4 Определение физических свойств природного газа. 3

4.1 Определение плотности. 3

4.2 Определение показателя адиабаты.. 4

4.3 Определение скорости звука. 6

4.4 Определение динамической вязкости. 6

5 Вычисление погрешности расчета физических свойств природного газа с учетом погрешности исходных данных. 6

6 Применение уравнения состояния для аттестации других методов расчета физических свойств природного газа. 7

Приложение А Листинг программы расчета физических свойств природного газа. 8

Приложение Б Пример расчета физических свойств природного газа. 24

Приложение В Библиография. 24

ГОСТ 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_м - молярная плотность, кмоль/м³;

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 кг/м³ (плотность газа при стандартных условиях рассчитывают по формуле (16) ГОСТ 30319.1);

по высшей удельной теплоте сгорания газа - 20 - 48 МДж/м³ (высшую удельную теплоту сгорания рассчитывают по 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, кг/м³, вычисляют по формуле

. (2)

Таблица 1 - Погрешности расчета свойств природного газа

Свойство	Область параметров состояния			Примечание
	240 £ Т < 270К		Т = (270 - 480) К и р < 12 МПа
	р £ 6МПа	6 < р £ 12 МПа	Т = (270 - 480) К и р < 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-го компонента природного газа (M_i) приведены в таблице 1 ГОСТ 30319.1 (см. п. 3.2.3).

4.1.3 Если компонентный состав природного газа задан в объемных долях, то молярные доли компонентов рассчитывают по формуле (12) ГОСТ 30319.1 и далее молярную массу природного газа вычисляют по 4.1.2*.

4.1.3. (Новая редакция, Изм. № 1).

4.2 Определение показателя адиабаты

Показатель адиабаты природного газа при использовании УС (1) вычисляют по формуле

, (5)

где с_p и c_v - изобарная и изохорная теплоемкости,

A₁ - безразмерный комплекс УС (1).

Безразмерный комплекс А₁ УС (1) имеет вид

А₁ = . (6)

Изобарную и изохорную теплоемкости рассчитывают по следующим выражениям:

, (7)

, (8)

где c_vom - изохорная теплоемкость природного газа в идеально газовом состоянии, а безразмерные комплексы А₂ и А₃ имеют вид:

, (9)

. (10)

Изохорную теплоемкость в идеально газовом состоянии вычисляют по формулам:

; (11)

. (12)

Изобарную теплоемкость (c_р_оi) i-го компонента в идеально газовом состоянии определяют из соотношения

, (13)

где q_i = T/T_ni.

Температура Т_ni, пределы суммирования N₁_i и N_2i, а также константы (a_j)_i, и (b_j)_i уравнения (13) для i-го компонента природного газа приведены в таблице 2.

Таблица 2 - Константы уравнения (13)

Компонент (i)	j	(a_j)_i	(b_j)_i
Метан N_1i= 10 N_2i = 6 T_ni = 100 К	0	1,46696186 × 10²
	1	-6,56744186 × 10¹	-2,09233731 × 10²
	2	2,02698132 × 10¹	2,06925203 × 10²
	3	-4,20931845 × 10⁰	-1,35704831 × 10²
	4	6,06743008 × 10^-1	5,64368924 × 10¹
	5	-6,12623969 × 10^-2	-1,34496111 × 10¹
	6	4,30969226 × 10^-3	1,39664152 × 10⁰
	7	-2,06597572 × 10^-4
	8	6,42615810 × 10^-6
	9	-1,16805630 × 10^-7
	10	9,40958930 × 10^-10
Этан N₁_i= 6 N_2i = 5 T_ni = 100 К	0	6,81209760 × 10¹
	1	-3,06340580 × 10¹	-8,74070840 × 10¹
	2	9,52750290 × 10⁰	7,84813740 × 10¹
	3	-1,69471020 × 10⁰	-4,48658590 × 10¹
	4	1,76305850 × 10^-1	1,46543460 × 10¹
	5	-9,95454020 × 10^-3	-2,05183930 × 10⁰
	6	2,35364300 × 10^-4
Пропан N_1i= 6 N_2i = 4 T_ni = 100 К	0	-9,209726737 × 10¹
	1	3,070930782 × 10¹	1,748671280 × 10²
	2	-4,924017995 × 10⁰	-1,756054503 × 10²
	3	5,045358836 × 10^-1	8,874920732 × 10¹
	4	-3,140446759 × 10^-2	-1,720610207 × 10¹
	5	1,076680079 × 10^-3
	6	-1,556890669 × 10^-5
н-Бутан N_1i= 6 N_2i = 5 T_ni = 100 К	0	-2,096096482 × 10²
	1	6,877783535 × 10¹	4,055272850 × 10²
	2	-1,228650555 × 10¹	-4,457015773 × 10²
	3	1,413691547 × 10⁰	2,743667350 × 10²
	4	-1,002920638 × 10^-1	-8,643867287 × 10¹
	5	3,985571861 × 10^-3	1,070428636 × 10¹
	6	-6,786460870 × 10^-5
и-Бутан N₁_i= 5 N_2i = 2 T_ni = 300 К	0	-3,871419306 × 10¹
	1	4,711104578 × 10¹	2,171601450 × 10¹
	2	-1,758225423 × 10¹	-4,492603200 × 10⁰
	3	4,183494309 × 10⁰
	4	-5,520042474 × 10^-1
	5	3,034658409 × 10^-2
Азот N_1i= 6 N_2i = 6 T_ni = 100 К	0	0,113129000 × 10²
	1	-0,215960000 × 10¹	-0,174654000 × 10²
	2	0,352761000 × 10⁰	0,246205000 × 10²
	3	-0,321705000 × 10^-1	-0,217731000 × 10²
	4	0,167690000 × 10^-2	0,116418000 × 10²
	5	-0,467965000 × 10^-4	-0,342122000 × 10¹
	6	0,542603000 × 10^-6	0,422296000 × 10⁰
Диоксид углерода N_1i= 6 N_2i = 4 T_ni = 300 К	0	-9,508041394 × 10^-1
	1	7,008743711 × 10⁰	1,087462263 × 10⁰
	2	-3,505801670 × 10⁰	-7,976765747 × 10^-2
	3	1,096778000 × 10⁰	-2,837014896 × 10^-3
	4	-2,016835088 × 10^-1	1,479612229 × 10^-4
	5	1,971024237 × 10^-2
	6	-7,860765734 × 10^-4
Сероводород N₁_i= 5 N₂_i = 5 T_ni = 100 К	0	3,913550000 × 10⁰
	1	-6,848510000 × 10^-2	0,0
	2	5,644240000 × 10^-2	0,0
	3	-4,837450000 × 10^-3	1,186580000 × 10⁰
	4	1,717820000 × 10^-4	-1,907470000 × 10⁰
	5	-2,275370000 × 10^-6	8,285200000 × 10^-1

(Измененная редакция, Изм. № 1).

4.3 Определение скорости звука

Скорость звука природного газа при использовании УС (1) вычисляют по формуле

, (14)

где c_p, c_v и А₁ - соответственно изобарная, изохорная теплоемкости природного газа и безразмерный комплекс УС (1), см. (6) - (13);

М - молярная масса природного газа, см. (3) или (4).

4.4 Определение динамической вязкости

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

, (15)

(16)

, (17)

, (18)

Молярную массу природного газа (М) вычисляют по формуле (3) или (4), а формулы расчета фактора Питцера (W) и псевдокритических параметров природного газа (Т_п, r_п, Т_пк, r_пк) приведены в ГОСТ 30319.2 (см. п. 3.2.5).

5 Вычисление погрешности расчета физических свойств природного газа с учетом погрешности исходных данных

При измерении расхода и количества природного газа, транспортируемого в газопроводах, давление (р), температуру (Т) и состав (х_i) измеряют с определенной погрешностью. Перечисленные параметры являются исходными данными для расчета физических свойств по УС (1) и уравнению для вязкости (15).

В соответствии с рекомендациями ИСО 5168 [5] погрешность расчета физических свойств, которая появляется в связи с погрешностью измерения исходных данных, определяют по формуле

d_ид = , (19)

где d_ид - погрешность расчета свойства Q, связанная с погрешностью измерения исходных данных;

d_qk - погрешность измерения параметра исходных данных;

= ; (20)

. (21)

В формулах (19) - (21):

q_k - условное обозначение k-го параметра исходных данных (р, Т, х_i);

`q_k - среднее значение k-го параметра в определенный промежуток времени (сутки, месяц, год и т.д.);

q_k^макс и q_k^мин - максимальное и минимальное значения k-го параметра в определенный промежуток времени;

Q - условное обозначение свойства природного газа (r, к, и, m);

N_q - количество параметров исходных данных, N_q = 2 + N (N - количество основных компонентов природного газа, которыми являются: метан, этан, пропан, бутаны, азот, диоксид углерода, сероводород).

При вычислении частных производных по формуле (20) свойства Q_qk₊ и Q_qk_- рассчитывают при средних параметрах и параметрах q_k₊ = + D и q_k_- = - D, соответственно. Рекомендуется выбирать D = 0,5×10^-2 d_qk.

Свойство Q (среднее значение) рассчитывают при средних параметрах `q_k.

Общую погрешность расчета физических свойств определяют по формуле

, (22)

где d_Q - погрешность расчета физических свойств по УС (1) и по уравнению для вязкости (15), значение которой для каждого свойства приведено в таблице 1.

(Измененная редакция, Изм. № 1).

6 Применение уравнения состояния для аттестации других методов расчета физических свойств природного газа

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

Таблица 3

Компонент	Концентрация компонентов, мол.%, при r_с, кг/м³
Компонент	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)

Если для аттестуемых методов в качестве исходных данных используют плотность смеси природного газа при стандартных условиях (r_c), ее значение для тестовых смесей необходимо рассчитывать по УС (1). Допускается также рассчитывать плотность r_c по формуле (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

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 кг/м³

Показатель адиабаты ......................................................................... 1,29

Скорость звука ................................................................................... 429,8 м/с

Динамическая вязкость ..................................................................... 12,36 мкПа × с

Давление ............................................................................................. 9,950 МПа

Температура ....................................................................................... 323,15 К

Плотность ........................................................................................... 78,51 кг/м³

Показатель адиабаты ......................................................................... 1,44

Скорость звука ................................................................................... 427,7 м/с

Динамическая вязкость ..................................................................... 14,75 мкПа × с

ПРИЛОЖЕНИЕ В

(справочное)

Библиография

[1] Козлов А.Д., Кузнецов В.М., Мамонов Ю.В. Построение уравнений теплофизических свойств индивидуальных веществ и материалов. - Теплофизические свойства веществ и материалов, 1988, вып. 24, с. 150-164

[2] Козлов А.Д., Кузнецов В.М., Мамонов Ю.В. Анализ современных методов расчета рекомендуемых справочных данных о коэффициентах вязкости и теплопроводности газов и жидкостей. - М.: ИВТАН СССР, 1989, № 3, с. 3-80

[3] МР 67-89. Расчет плотности, изобарной и изохорной теплоемкости, энтальпии, энтропии, скорости звука жидких и газообразных веществ, применяемых в криогенном машиностроении в интервале температур до 500 К и давлений до 50 МПа на основе уравнения Старлинга-Хана. - Методика ГСССД, Деп. ВНИИКИ, № 609, 1990

[4] B.A.Younglove, N.V.Frederick, R.D.McCarty Speed of Sound Data and Related Models for Mixtures of Natural Gas Constituents - Natl. Inst. Stand. Technol., Mono. 178, 97 p. (Washington, 1993)

[5] ISO 5168:1978 International Standard. Measurement of fluid flow - Estimation of uncertainty of a flow-rate measurement

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

Информационная система«Ёшкин Кот»

4.2 Определение показателя адиабаты

4.3 Определение скорости звука

4.4 Определение динамической вязкости

Информационная система
«Ёшкин Кот»