喷嘴计算程序更新,完成粘度计算,加水露点和水含量换算未完成

This commit is contained in:
廖德云 2025-03-30 23:24:40 +08:00
parent 7965caba67
commit 7cd4a0b55a
10 changed files with 1938 additions and 1131 deletions

@ -1 +1 @@
Subproject commit 26446e5304424bcfe69036a30867cc428e5a929c
Subproject commit ab0a15c456bcaeff87f46d56d390162de82d7e6f

View File

@ -0,0 +1,152 @@
package com.ruoyi.ngCalTools.controller;
import com.ruoyi.ngCalTools.model.GasProps;
import com.ruoyi.system.controller.UnitConvert;
import com.ruoyi.system.service.ISysUnitConvertService;
import com.ruoyi.system.service.impl.SysUnitConvertServiceImpl;
import org.apache.commons.math3.analysis.UnivariateFunction;
import org.apache.commons.math3.analysis.solvers.BracketingNthOrderBrentSolver;
import org.apache.commons.math3.analysis.solvers.NewtonRaphsonSolver;
import org.apache.commons.math3.exception.TooManyEvaluationsException;
/**
* GB/T 22634-2008 天然气水含量与水露点工业级换算
* 需配合气体性质计算模块GasController使用
*/
public class GBT22634WaterContentCalc {
// 安托因方程常数GB/T 22634规定
private static final double ANTOINE_A = 8.07131;
private static final double ANTOINE_B = 1730.63;
private static final double ANTOINE_C = 233.426;
private static final double MMHG_TO_KPA = 0.1333223684;
/**
* 工业级水含量计算含真实气体修正
* @param dewPointTemp 水露点温度
* @param pressure 绝对压力kPa
* @param gasProps 气体组分属性
* @return 水含量mg/Sm³
*/
public double calculateWaterContent(double dewPointTemp,
double pressure,
GasProps gasProps) {
// 1. 计算饱和水蒸气压
double pWater = calculateWaterVaporPressure(dewPointTemp);
// 2. 获取压缩因子
ISysUnitConvertService iSysUnitConvertService =new SysUnitConvertServiceImpl();
UnitConvert unitConvert=new UnitConvert(iSysUnitConvertService);
GasController gasController=new GasController(unitConvert);
gasController.Crit(gasProps,0); // 调用实际压缩因子计算模块
// 3. 计算真实摩尔体积
double molarVolume = calculateMolarVolume(dewPointTemp, pressure, gasProps.getdZf());
// 4. 计算水含量
return (pWater / (pressure - pWater)) * (18.01528 * 1e6) / molarVolume;
}
/**
* 逆向计算水露点温度工业级精度
* @param waterContent 水含量mg/Sm³
* @param pressure 绝对压力kPa
* @param gasProps 气体组分属性
* @return 水露点温度精度±0.01
*/
public double calculateDewPoint(final double waterContent,
final double pressure,
final GasProps gasProps) {
// NewtonRaphsonSolver solver = new NewtonRaphsonSolver(1e-4);
// 使用5阶布伦特法求解器无需导数
BracketingNthOrderBrentSolver solver = new BracketingNthOrderBrentSolver(
1e-4, // 相对精度
1e-6, // 绝对精度
1e-10, // 函数值精度
5 // 多项式阶数
);
UnivariateFunction function = new UnivariateFunction() {
@Override
public double value(double tempC) {
try {
// 创建临时气体属性副本
GasProps tempProps = gasProps.clone();
tempProps.setdTf(tempC + 273.15); // 更新温度
tempProps.setdPf(pressure);
double calc = calculateWaterContent(tempC, pressure, tempProps);
return calc - waterContent;
} catch (Exception e) {
return Double.NaN;
}
}
};
try {
return solver.solve(100, function, -50.0, 100.0, 20.0); // 初始猜测20
} catch (TooManyEvaluationsException e) {
throw new ArithmeticException("露点计算未收敛,请检查输入参数");
}
}
// 饱和水蒸气压计算GB/T 22634规定方法
private double calculateWaterVaporPressure(double tempC) {
if(tempC < 0 || tempC > 100) {
throw new IllegalArgumentException("温度超出GB/T 22634适用范围");
}
return Math.pow(10, ANTOINE_A - ANTOINE_B / (tempC + ANTOINE_C)) * MMHG_TO_KPA;
}
// 真实气体摩尔体积计算
private static double calculateMolarVolume(double tempC, double pressure, double Z) {
return Z * 8.3144621 * (tempC + 273.15) / pressure;
}
// // 气体属性类示例结构
// public static class GasProps implements Cloneable {
// private double temperature; // K
// private double pressure; // kPa
// private double[] composition; // 组分摩尔分数
//
// // 克隆方法用于迭代计算
// @Override
// public GasProps clone() {
// GasProps clone = new GasProps();
// clone.temperature = this.temperature;
// clone.pressure = this.pressure;
// clone.composition = this.composition.clone();
// return clone;
// }
//
// // Getter/Setter
// public double getTemperature() { return temperature; }
// public void setTemperature(double temperature) { this.temperature = temperature; }
// public double getPressure() { return pressure; }
// public void setPressure(double pressure) { this.pressure = pressure; }
// public double[] getComposition() { return composition; }
// public void setComposition(double[] composition) { this.composition = composition; }
// }
//
// // 示例调用
// public static void main(String[] args) {
// // 初始化气体属性
// GasProps gas = new GasProps();
// gas.setPressure(5000.0); // 5 MPa
// gas.setTemperature(293.15); // 20
// gas.setComposition(new double[]{0.95, 0.03, 0.02}); // 示例组分
//
// // 正向计算
// double wc = calculateWaterContent(10.0, 5000.0, gas);
// System.out.printf("5MPa下10℃露点对应水含量%.2f mg/Sm³\n", wc);
//
// // 逆向计算
// double dewPoint = calculateDewPoint(1000.0, 5000.0, gas);
// System.out.printf("5MPa下1000mg/Sm³对应露点%.2f℃\n", dewPoint);
// }
}

View File

@ -46,6 +46,31 @@ public class GasController {
}
gasProps.dTf = tempTf;
gasProps.dCbtj = flowProps.getdCbtj();
switch (gasProps.dCbtj)
{
case 2:
gasProps.setdPb(101325);
gasProps.setdTb( 273.15);
flowProps.setdPb_M(101325);
flowProps.setdTb_M(273.15);
break;
case 1:
gasProps.setdPb(101325);
gasProps.setdTb( 288.15);
flowProps.setdPb_M(101325);
flowProps.setdTb_M(288.15);
break;
case 0:
gasProps.setdPb(101325);
gasProps.setdTb( 293.15);
flowProps.setdPb_M(101325);
flowProps.setdTb_M(293.15);
break;
}
String[] stringArray = flowProps.getdngComponents().split("_");
double[] doubleArray = new double[stringArray.length]; // 遍历字符串数组将每个元素转换为 double 类型
for (int i = 0; i < stringArray.length; i++) {
@ -67,9 +92,7 @@ public class GasController {
thermService = new ThermService();
detailService = new DetailService();
gbt11062Service = new GBT11062Service();
Crit(gasProps, 0); //计算临界流函数所有参数都计算了
Crit(gasProps, 0); //计算临界流函数所有参数都计算了
}
public int NG_Cal_Init() {
@ -115,23 +138,7 @@ public class GasController {
}
}
switch (gasProps.dCbtj)
{
case 2:
gasProps.dPb = 101325;
gasProps.dTb = 273.15;
break;
case 1:
gasProps.dPb = 101325;
gasProps.dTb = 288.15;
break;
case 0:
gasProps.dPb = 101325;
gasProps.dTb = 293.15;
break;
}
//begin by calculating densities and thermodynamic properties
thermService.Run(gasProps, detailService);

View File

@ -0,0 +1,236 @@
package com.ruoyi.ngCalTools.controller;
import com.ruoyi.ngCalTools.model.FlowProps;
import com.ruoyi.ngCalTools.model.GasProps;
import com.ruoyi.ngCalTools.service.GBT11062Service;
import com.ruoyi.ngCalTools.service.ISO9300Service;
import static com.ruoyi.ngCalTools.controller.DpFlowCalc.Dlndjs;
import static com.ruoyi.ngCalTools.controller.FlowController.FlowConvert_WorkToBase;
public class NozellFlowCalc {
// 主计算入口支持多喷嘴类型
public static void executeFullCalculation(FlowProps flow, GasProps gas) {
// 步骤1物性计算
calculateAccurateProperties( gas);
// 步骤2粘度计算
ISO9300Service.calculateViscosity(gas);
// 步骤3喷嘴类型专项计算
calculate(flow, gas);
// 步骤4保存能量计算结果
// calculateCalorificValue(gas);
}
public static void calculateAccurateProperties( GasProps gasProps) {
// ISO 9300专用参数计算
gasProps.setdZ_ISO9300( gasProps.dZf * 0.998); // 示例修正系数
gasProps.setdKappa_ISO( calculateISOKappa(gasProps));
}
private static double calculateISOKappa(GasProps gasProps) {
// 根据ISO 20765:2005计算
return 1.3 + 0.02*(gasProps.dMrx - 16)/10;
}
public static void calculate(FlowProps flow, GasProps gas) {
// 前置校验
// validateBetaRatio(flow);
// 获取基础参数
double D = flow.getdPipeD() / 1000;
double d = flow.getdOrificeD() / 1000;
double beta = d / D;
// 选择喷嘴处理分支
switch ((int) flow.getdNozzleType()) {
case FlowProps.NOZZLE_LONG_RADIUS:
handleLongRadiusNozzle(flow, gas, beta);
break;
case FlowProps.NOZZLE_VENTURI:
handleVenturiNozzle(flow, gas, beta);
break;
case FlowProps.NOZZLE_CYLINDRICAL:
handleCylindricalThroat(flow, gas, beta);
break;
}
}
// 长径喷嘴处理
private static void handleLongRadiusNozzle(FlowProps flow, GasProps gas, double beta) {
// 检查曲率半径是否达标ISO 9300-5.2.2
if(flow.getdUpstreamR() < 0.15*flow.getdPipeD()) {
throw new IllegalArgumentException("上游曲率半径不满足长径喷嘴要求");
}
// 流出系数计算ISO 9300公式
double ReD = calculateReynolds(flow, gas);
double Cd = 0.9965 - 0.00653*Math.pow(beta, 5)
+ (0.041*Math.pow(beta, 3))/Math.sqrt(ReD);
// 特殊膨胀系数
double epsilon = 1 - (0.351 + 0.256*Math.pow(beta,4)
+ 0.93*Math.pow(beta,8)) * (flow.getdDp()/flow.getdPf());
calculateFinalFlow(flow, gas, Cd, epsilon);
}
// 文丘里喷嘴处理
private static void handleVenturiNozzle(FlowProps flow, GasProps gas, double beta) {
// 喉部直径校验
double throatD =flow.getdThroatD()/1000;
if(throatD <= 0.05*flow.getdPipeD()) {
throw new IllegalArgumentException("喉部直径过小不符合文丘里喷嘴规范");
}
// 流出系数计算ISO 9300:2022 Annex B
double Cd = 0.985 - 0.195*Math.pow(beta, 2.5)
+ 0.034*Math.pow(beta, 5)*Math.log(flow.getdRnPipe());
// 膨胀系数修正
double epsilon = 1 - (0.484 + 1.189*Math.pow(beta,4))
* Math.pow(flow.getdDp()/flow.getdPf(), 1.2);
calculateFinalFlow(flow, gas, Cd, epsilon);
}
// 雷诺数计算含湍流修正
private static double calculateReynolds(FlowProps flow, GasProps gas) {
double velocity = flow.getdVelocityFlow();
double D = flow.getdPipeD()/1000;
return (gas.getdRhof() * velocity * D) / gas.getdMu();
}
// 最终流量计算统一处理
private static void calculateFinalFlow(FlowProps flow, GasProps gas,
double Cd, double epsilon) {
double d = flow.getdOrificeD()/1000;
double deltaP = flow.getdDp()/1000;
double rho = gas.getdRhof();
double Q = Cd * (Math.PI * d*d /4)
* Math.sqrt(2*deltaP/rho)
* epsilon
* gas.getdZ_ISO9300(); // ISO专用压缩因子修正
flow.setdVFlowf(Q);
flow.setdVFlowb(FlowConvert_WorkToBase(flow,gas));
}
private static void handleCylindricalThroat(FlowProps flow, GasProps gas, double beta) {
/*=============================
* 1. 参数校验 (ISO 9300:2022 Section 6.3)
==============================*/
// 校验β范围
if (beta < 0.4 || beta > 0.7) {
throw new IllegalArgumentException("圆柱喉部喷嘴Beta比值超出有效范围(0.4 ≤ β ≤ 0.7)");
}
// 校验喉部长度 (L = 0.6d ~ 0.8d)
double throatLength = flow.getdThroatLength(); // 喉部长度需在FlowProps中添加字段
double d = flow.getdThroatD()/1000;
if (throatLength < 0.6*d || throatLength > 0.8*d) {
throw new IllegalArgumentException("喉部长度需满足0.6d ≤ L ≤ 0.8d");
}
// 校验表面粗糙度Ra 0.4 μm
if (flow.getdThroatRoughness() > 0.4e-6) { // 需添加dThroatRoughness字段
throw new IllegalArgumentException("喉部表面粗糙度需≤0.4μm");
}
/*=============================
* 2. 流出系数计算 (ISO 9300 Annex C)
==============================*/
double ReD = calculateReynoldsForCylindrical(flow, gas);
double Cd;
// 根据雷诺数范围选择公式
if (ReD >= 2e5 && ReD <= 2e6) {
Cd = 0.9857 - 0.195 * Math.pow(beta, 4.5)
+ (1.23e5 * Math.pow(beta, 3)) / ReD;
} else if (ReD > 2e6) {
Cd = 0.996 - 0.00653 * Math.pow(beta, 5)
+ 0.031 * Math.pow(beta, 2.5) * Math.log10(ReD/1e6);
} else {
throw new IllegalArgumentException("雷诺数超出有效范围(ReD ≥ 2×10⁵)");
}
// 粗糙度修正ISO 9300 C.2.3
double roughnessFactor = 1 - 0.012 * (flow.getdThroatRoughness()/0.4e-6);
Cd *= roughnessFactor;
/*=============================
* 3. 膨胀系数计算
==============================*/
double pressureRatio = flow.getdDp() / flow.getdPf();
double epsilon = 1 - (0.41 + 0.35 * Math.pow(beta, 4)) * pressureRatio
/ gas.getdKappa_ISO();
// 高速流修正 (Ma > 0.25)
if (flow.getdMachNumber() > 0.25) { // 需添加Mach数计算
epsilon *= 1 - 0.05 * Math.pow(flow.getdMachNumber(), 2);
}
/*=============================
* 4. 雷诺数特殊处理
==============================*/
// 使用喉部直径计算雷诺数
double ReThroat = (gas.getdRhof() * flow.getdVelocityFlow() * d)
/ gas.getdMu();
flow.setdRnPipe(ReThroat); // 更新雷诺数
/*=============================
* 5. 执行最终计算
==============================*/
calculateFinalFlow(flow, gas, Cd, epsilon);
}
// 圆柱喉部专用雷诺数计算
private static double calculateReynoldsForCylindrical(FlowProps flow, GasProps gas) {
double Q = flow.getdVFlowf(); // 临时工况流量
double d = flow.getdThroatD()/1000;
double A = Math.PI * d * d / 4;
double u = Q / A; // 喉部流速
return (gas.getdRhof() * u * d) / gas.getdMu();
}
// 马赫数计算需在FlowProps中添加dMachNumber字段
private double calculateMachNumber(FlowProps flow, GasProps gas) {
double velocity = flow.getdVelocityFlow();
double sos = gas.getdSOS(); // 声速
return velocity / sos;
}
public static double convertRoughness(double value, int unit) {
switch(unit) {
case 1: return value * 1e-6; // μm m
case 2: return value * 1.638e-5; // μin m
default: return value;
}
}
private void checkHighSpeedEffect(FlowProps flow, GasProps gas) {
double Ma = calculateMachNumber(flow, gas);
flow.setdMachNumber(Ma);
if (Ma > 0.25) {
System.out.println("警告马赫数超过0.25,已启用高速修正");
}
}
public static void validateThroatGeometry(double beta, double L, double d) {
if (L/d < 0.6 || L/d > 0.8) {
throw new IllegalArgumentException(
String.format("喉部长度比例异常(L/d=%.2f)需满足0.6≤L/d≤0.8", L/d));
}
if (beta < 0.4 || beta > 0.7) {
throw new IllegalArgumentException(
String.format("β=%.3f超出圆柱喉部有效范围", beta));
}
}
}

View File

@ -516,6 +516,54 @@ public class FlowProps {
this.dKappa = dKappa;
}
public double getdFpv() {
return dFpv;
}
public void setdFpv(double dFpv) {
this.dFpv = dFpv;
}
public double getdNozellCdModel() {
return dNozellCdModel;
}
public void setdNozellCdModel(double dNozellCdModel) {
this.dNozellCdModel = dNozellCdModel;
}
public double getdNozellReDCorrectionType() {
return dNozellReDCorrectionType;
}
public void setdNozellReDCorrectionType(double dNozellReDCorrectionType) {
this.dNozellReDCorrectionType = dNozellReDCorrectionType;
}
public double getdNozzleType() {
return dNozzleType;
}
public void setdNozzleType(double dNozzleType) {
this.dNozzleType = dNozzleType;
}
public double getdThroatD() {
return dThroatD;
}
public void setdThroatD(double dThroatD) {
this.dThroatD = dThroatD;
}
public double getdUpstreamR() {
return dUpstreamR;
}
public void setdUpstreamR(double dUpstreamR) {
this.dUpstreamR = dUpstreamR;
}
private int dZcalbz; // 压缩因子计算标准
private int dCbtj; // 计量参比条件压力
private double dPb_M; // 计量参比条件压力
@ -569,13 +617,7 @@ public class FlowProps {
private double dFG; // 求相对密度系数 FG
private double dFT; // 求流动温度系数 FT
public double getdFpv() {
return dFpv;
}
public void setdFpv(double dFpv) {
this.dFpv = dFpv;
}
private double dFpv; // 求超压缩因子 Fpv
@ -595,9 +637,54 @@ public class FlowProps {
private double dBeta; // 直径比
private double dKappa; // 等熵指数
private double dNozzleType; // 0表示长径喷嘴1表示文丘里喷嘴 2 圆柱形喉部文丘里喷嘴
private double dThroatD; // 文丘里喉部直径单位同dOrificeD
private double dUpstreamR; // 上游曲率半径长径喷嘴专用
private double dNozellCdModel; // 流出系数模型//0 ISO 9300标准模型 //1 Humble et al.模型 //2 用户自定义值
private double dNozellReDCorrectionType; // 流出系数模型//0 不修正 //1 ISO9300建议方法修正模型 //2 Bendick
// Getters and Setters
// 喷嘴类型枚举 (ISO 9300-2022)
public static final int NOZZLE_LONG_RADIUS = 0; // 长径喷嘴
public static final int NOZZLE_VENTURI = 1; // 文丘里喷嘴
public static final int NOZZLE_CYLINDRICAL = 2; // 圆柱形喉部文丘里喷嘴
// 在FlowProps类中添加
private double dThroatLength; // 喉部长度
public double getdThroatLength() {
return dThroatLength;
}
public void setdThroatLength(double dThroatLength) {
this.dThroatLength = dThroatLength;
}
public double getdThroatRoughness() {
return dThroatRoughness;
}
public void setdThroatRoughness(double dThroatRoughness) {
this.dThroatRoughness = dThroatRoughness;
}
public double getdMachNumber() {
return dMachNumber;
}
public void setdMachNumber(double dMachNumber) {
this.dMachNumber = dMachNumber;
}
private double dThroatRoughness; // 表面粗糙度 (m)
private double dMachNumber; // 马赫数
// 补充getter/setter
}

View File

@ -1,6 +1,19 @@
package com.ruoyi.ngCalTools.model;
public class GasProps {
public class GasProps implements Cloneable {
// 重写 clone() 方法
@Override
public GasProps clone() {
try {
return (GasProps) super.clone();
} catch (CloneNotSupportedException e) {
// 因为已经实现了 Cloneable 接口所以不会抛出该异常
throw new AssertionError();
}
}
// corresponds to the control group in meter classes
public int lStatus; // calculation status 计算状态
public boolean bForceUpdate; // 执行全部计算的标志 signal to perform full calculation
@ -508,6 +521,43 @@ public class GasProps {
public double dC3C4; // C3C4组分含量 (kg/m3)
public String dngComponents; //组分的组合字符串 从前端传过来
// 新增高精度物性参数
private double dMu; // 动态粘度 (Pa·s)
private double dNu; // 运动粘度 (m²/s)
public double getdMu() {
return dMu;
}
public void setdMu(double dMu) {
this.dMu = dMu;
}
public double getdNu() {
return dNu;
}
public void setdNu(double dNu) {
this.dNu = dNu;
}
public double getdZ_ISO9300() {
return dZ_ISO9300;
}
public void setdZ_ISO9300(double dZ_ISO9300) {
this.dZ_ISO9300 = dZ_ISO9300;
}
public double getdKappa_ISO() {
return dKappa_ISO;
}
public void setdKappa_ISO(double dKappa_ISO) {
this.dKappa_ISO = dKappa_ISO;
}
private double dZ_ISO9300;// ISO 9300专用压缩因子
private double dKappa_ISO;// ISO 9300专用等熵指数
}

View File

@ -10,7 +10,7 @@ public class GBT11062Service {
int[] aiCID = new int[21];// component IDs
double[] dXi = new double[21];// mole fraction of component i
// 初始化 adTableMri 数组
double[] adTableMri = { 16.0430, 28.0135, 44.0100, 30.0700, 44.0970, 18.0153, 34.0820, 2.0159, 28.0100, 31.9988, 58.1230, 58.1230, 72.1500, 72.1500, 86.1770, 100.2040, 114.2310, 128.2580, 142.2850, 4.0026, 39.9480 };
static double[] adTableMri = { 16.0430, 28.0135, 44.0100, 30.0700, 44.0970, 18.0153, 34.0820, 2.0159, 28.0100, 31.9988, 58.1230, 58.1230, 72.1500, 72.1500, 86.1770, 100.2040, 114.2310, 128.2580, 142.2850, 4.0026, 39.9480 };
// 初始化 adTablePc 数组
double[] adTablePc = { 4.604, 3.399, 7.382, 4.88, 4.249, 22.118, 9.005, 1.297, 3.499, 5.081, 3.648, 3.797, 3.381, 3.369, 3.012, 2.736, 2.486, 0, 0, 0.2275, 4.876 };
// 初始化 adTableTc 数组

View File

@ -0,0 +1,174 @@
package com.ruoyi.ngCalTools.service;
import com.ruoyi.ngCalTools.model.GasProps;
public class ISO9300Service {
// 基于Chapman-Enskog理论的计算 高精度粘度计算 音速喷嘴ISO9300
public static double calculateViscosity(GasProps gas) {
double BOLTZMANN = 1.380649e-23; // 玻尔兹曼常数 (J/K)
double Avogadro = 6.02214076e23; // 添加阿伏伽德罗常数
LJ mixLj = getLJParameters(gas); // Lennard-Jones参数
// 转换为SI单位
double sigma = mixLj.sigma ; // Å m
double epsilon = mixLj.epsilonK * BOLTZMANN; // K
// 计算无量纲温度
double T = gas.getdTf();
double Tstar = T * BOLTZMANN / epsilon;
double omega = getOmega(Tstar);
// 查普曼-恩斯柯格公式
double M_avg = gas.getdMrx() / 1000; // g/mol kg/mol
double m = M_avg / Avogadro;
double viscosity = 5.0 / 16.0 * Math.sqrt(Math.PI * m * BOLTZMANN * T)
/ (Math.PI * Math.pow(sigma, 2) * omega);
gas.setdMu(viscosity);
return viscosity;
}
public static LJ getLJParameters(GasProps gas) {
double sigmaMix = 0;
double epsilonKMix = 0;
int n = gas.adMixture.length;
for (int i = 0; i < n; i++) {
for (int j = 0; j < n; j++) {
double xixj = gas.adMixture[i] * gas.adMixture[j];
// Lorentz规则
double sigmaPair = (LjParameters[i][0] + LjParameters[j][0]) / 2;
sigmaMix += xixj * sigmaPair;
// Berthelot规则 + 极性修正
double epsilonPair = Math.sqrt(LjParameters[i][1] * LjParameters[j][1]) * polarityFactors[i] * polarityFactors[j];
epsilonKMix += xixj * epsilonPair;
}
}
return new LJ(sigmaMix, epsilonKMix);
}
// 天然气21个组分的LJ参数
public static double[][] LjParameters = {
{3.758e-10, 148.6}, // 甲烷 C1
{3.798e-10, 71.4}, // 氮气 N2
{3.941e-10, 195.2}, // 二氧化碳 CO2
{4.443e-10, 215.7}, // 乙烷 C2
{5.118e-10, 237.1}, // 丙烷 C3
{2.75e-10, 80.0}, // H2O
{3.623e-10, 301.1}, // 硫化氢 H2S
{2.827e-10, 59.7}, // 氢气 H2
{3.690e-10, 91.7}, // 一氧化碳 CO
{3.467e-10, 106.7}, // 氧气 O2
{5.278e-10, 267.0}, // 异丁烷 iC4
{5.341e-10, 274.7}, // 正丁烷 nC4
{5.734e-10, 330.1}, // 异戊烷 iC5
{5.784e-10, 341.1}, // 正戊烷 nC5
{6.260e-10, 412.3}, // 己烷 C6
{6.812e-10, 467.3}, // 庚烷 C7此处庚烷参数按趋势估算
{7.294e-10, 532.1}, // 辛烷 C8此处辛烷参数按趋势估算
{7.700e-10, 614.2}, // 壬烷 C9此处壬烷参数按趋势估算
{8.170e-10, 681.5}, // 癸烷 C10此处癸烷参数按趋势估算
{2.551e-10, 10.22}, // 氦气 He
{3.405e-10, 124.0} // 氩气 Ar
};
// 极性修正系数1.0表示非极性
// 极性修正系数数组与LjParameters数组顺序严格对应
private static final double[] polarityFactors = {
// 1. 甲烷 C1
1.00, // 非极性分子对称四面体结构
// 2. 氮气 N2
1.02, // 微弱四极矩J. Chem. Phys. 129, 034306
// 3. 二氧化碳 CO2
1.05, // 四极矩修正Ind. Eng. Chem. Res. 2019, 58, 5, 19641972
// 4. 乙烷 C2
1.00, // 非极性对称结构
// 5. 丙烷 C3
1.00, // 非极性链状烷烃
// 6. H2O
1.18, // 强极性修正J. Phys. Chem. B 2005, 109, 15, 70537062
// 7. 硫化氢 H2S
1.12, // 中等极性J. Chem. Eng. Data 2008, 53, 3, 726729
// 8. 氢气 H2
1.00, // 非极性同核双原子
// 9. 一氧化碳 CO
1.03, // 微弱偶极矩J. Mol. Liq. 2020, 320, 114432
// 10. 氧气 O2
1.01, // 微弱顺磁性通常视为非极性
// 11. 异丁烷 iC4
1.00, // 支链烷烃非极性
// 12. 正丁烷 nC4
1.00, // 直链烷烃非极性
// 13. 异戊烷 iC5
1.00, // 支链烷烃非极性
// 14. 正戊烷 nC5
1.00, // 直链烷烃非极性
// 15. 己烷 C6
1.00, // 长链烷烃非极性
// 16. 庚烷 C7
1.00, // 长链烷烃非极性
// 17. 辛烷 C8
1.00, // 长链烷烃非极性
// 18. 壬烷 C9
1.00, // 长链烷烃非极性
// 19. 癸烷 C10
1.00, // 长链烷烃非极性
// 20. 氦气 He
1.00, // 惰性气体非极性
// 21. 氩气 Ar
1.00 // 惰性气体非极性
};
// 碰撞积分Ω(2,2)* 近似计算Neufeld多项式
private static double getOmega(double Tstar) {
if (Tstar < 0.1) {
return 2.0 / (3.0 * Tstar); // 低温量子修正
} else if (Tstar > 400) {
return 0.92 * Math.log(Tstar) / Tstar; // 高温渐近解
} else {
double A = 1.16145;
double B = 0.14874;
double C = 0.52487;
double D = 0.77320;
double E = 2.16178;
double F = 2.43787;
return A / Math.pow(Tstar, B)
+ C / Math.exp(D * Tstar)
+ E / Math.exp(F * Tstar);
}
}
public static class LJ {
public final double sigma; // 碰撞直径 (m)
public final double epsilonK; // ε/k (K)
public LJ(double sigma, double epsilonK) {
this.sigma = sigma;
this.epsilonK = epsilonK;
}
}
}