汇川PLC ,流量计算库，流量计算机

NGFlowMeter/CppProject/NGTools/Detail.h

@@ -0,0 +1,120 @@
+/*************************************************************************
+* <EFBFBD>ļ<EFBFBD>:    detail.h  
+**************************************************************************/
+
+#ifndef _DETAIL_H
+#define _DETAIL_H
+
+#include "NGCal.h"
+
+
+typedef struct Detail {
+
+    int iNCC;
+
+    int aiCID[21];
+
+    double dOldMixID;
+    double dOldPb;
+    double dOldTb;
+    double dOldPf;
+    double dOldTf;
+
+
+    double adAn[58];
+    double adUn[58];
+
+
+    double dMri[21];
+    double dEi[21];
+    double dKi[21];
+    double dGi[21];
+    double dQi[21];
+    double dFi[21];
+    double dSi[21];
+    double dWi[21];
+
+    double dEij[21][21];
+    double dUij[21][21];
+    double dKij[21][21];
+    double dGij[21][21];
+
+    double adTable6Eij[21][21];
+    double adTable6Uij[21][21];
+    double adTable6Kij[21][21];
+    double adTable6Gij[21][21];
+
+    double adTable5Qi[21];
+    double adTable5Fi[21];
+    double adTable5Si[21];
+    double adTable5Wi[21];
+
+    double dXi[21];
+    double dPCalc;
+    double dT;
+    double dP;
+    double dRhoTP;
+    double dB;
+    double adBcoef[18];
+    double adFn[58];
+    double fx[58];
+    double dU;
+    double dKp3;
+    double dW;
+    double dQp2;
+    double dF;
+    double dRho;
+    double dRhoL;
+    double dRhoH;
+    double dPRhoL;
+    double dPRhoH;
+
+
+    double dZ;
+    double ddZdT;
+    double dd2ZdT2;
+    double ddZdD;
+    double ddBdT;
+    double dd2BdT2;
+} Detail;
+
+
+Detail *Detail_Construct(void);
+
+void Detail_Destroy(Detail *pDetail);
+
+
+int Detail_compositionchange(Detail *pDetail, const NGParSTRUCT *pAGA10);
+
+int Detail_table(Detail *pDetail);
+
+void Detail_paramdl(Detail *pDetail);
+
+void Detail_chardl(Detail *pDetail, NGParSTRUCT *pAGA10);
+
+void Detail_dhvMol(Detail *pDetail, NGParSTRUCT *pAGA10);
+
+void Detail_bvir(Detail *pDetail);
+
+void Detail_temp(Detail *pDetail);
+
+void Detail_braket(Detail *pDetail, NGParSTRUCT *pAGA10);
+
+void Detail_pdetail(Detail *pDetail, double dRho);
+
+void Detail_ddetail(Detail *pDetail, NGParSTRUCT *pAGA10);
+
+void Detail_relativedensity(const Detail *pDetail, NGParSTRUCT *pAGA10);
+
+
+double Detail_zdetail(Detail *pDetail, double dRho);
+
+double Detail_dZdT(Detail *pDetail, double dRho);
+
+double Detail_d2ZdT2(Detail *pDetail, double dRho);
+
+double Detail_dZdD(Detail *pDetail, double dRho);
+
+void Detail_Run(Detail *pDetail, NGParSTRUCT *ptNGPar);
+
+#endif

NGFlowMeter/CppProject/NGTools/FlowCal.c

@@ -0,0 +1,501 @@
+//
+// Created by ldeyu on 2025/7/7.
+//
+#include "NGCal.h"
+#include "FlowCal.h"
+#include "math.h"
+
+
+double format_double(double value, int digits) {
+    // 处理默认位数（4位）
+    if (digits == 0) {
+        digits = 4;
+    }
+
+    // 验证位数有效性
+    if (digits < 1 || digits > 5) {
+        fprintf(stderr, "Error: Invalid digit value (must be 1-5 or 0 for default)\n");
+        return NAN; // 返回 NaN 表示错误
+    }
+
+    char format_str[10];
+    char buffer[50];
+
+    // 生成动态格式字符串（如 "%.4f"）
+    snprintf(format_str, sizeof(format_str), "%%.%df", digits);
+
+    // 格式化数值到字符串
+    snprintf(buffer, sizeof(buffer), format_str, value);
+
+    // 转换回 double
+    double result = strtod(buffer, NULL);
+    return result;
+}
+
+void OFlowCal(FlowParSTRUCT *ptFlowPar, NGParSTRUCT *ptNGPar) {
+    double tempPatm = ptFlowPar->dPatm * 1000000;
+
+    double tempPf = ptFlowPar->dPf * 1000000;
+
+
+    double tempDP = ptFlowPar->dDp * 1000;
+
+    double tempTf = ptFlowPar->dTf + 273.15;
+
+    if (ptFlowPar->dPfType == 0) {
+        ptFlowPar->dPf = tempPatm + tempPf;
+        ptNGPar->dPf = tempPatm + tempPf;
+    } else {
+        ptFlowPar->dPf = tempPf;
+        ptNGPar->dPf = tempPf;
+    }
+
+    ptFlowPar->dDp = tempDP;
+    ptFlowPar->dTf = tempTf;
+    ptNGPar->dTf = tempTf;
+    ptNGPar->dCbtj = ptFlowPar->dCbtj;
+
+    switch (ptNGPar->dCbtj) {
+        case 2:
+            ptNGPar->dPb = 101325;
+            ptNGPar->dTb = 273.15;
+            ptFlowPar->dPb_M = (101325);
+            ptFlowPar->dTb_M = (273.15);
+            break;
+        case 1:
+
+            ptNGPar->dPb = (101325);
+            ptNGPar->dTb = (288.15);
+            ptFlowPar->dPb_M = (101325);
+            ptFlowPar->dTb_M = (288.15);
+            break;
+        case 0:
+
+            ptNGPar->dPb = (101325);
+            ptNGPar->dTb = (293.15);
+            ptFlowPar->dPb_M = (101325);
+            ptFlowPar->dTb_M = (293.15);
+            break;
+        default: ;
+    }
+
+    double ngArray[NUMBEROFCOMPONENTS];
+int i=0;
+    for ( i = 0; i < NUMBEROFCOMPONENTS; i++) {
+        ngArray[i] = ptFlowPar->dNG_Compents[i] / 100;
+        ptNGPar->adMixture[i] = ngArray[i];
+    }
+    Crit(ptNGPar, 0);
+
+
+    ptFlowPar->dFpv = format_double(ptNGPar->dFpv, 4);
+
+
+    ptFlowPar->dOrificeD = format_double(ptFlowPar->dOrificeD * (
+                                             1 + 0.000001 * CaiLiaoPzxs(ptFlowPar->dOrificeMaterial) * (
+                                                 ptFlowPar->dTf - 293.15)), 2);
+    ptFlowPar->dPipeD = format_double(ptFlowPar->dPipeD * (1 + 0.000001 * CaiLiaoPzxs(ptFlowPar->dPipeMaterial) * (
+                                                               ptFlowPar->dTf - 293.15)), 2);
+    ptFlowPar->dBeta = format_double(ptFlowPar->dOrificeD / ptFlowPar->dPipeD, 4);
+
+    ptFlowPar->dE = format_double(calculateE(ptFlowPar->dBeta), 4);
+    ptFlowPar->dFG = format_double(calculateFG(ptNGPar->dRD_Real), 4);
+    ptFlowPar->dFT = format_double(calculateFT(ptFlowPar->dTb_M, ptFlowPar->dTf), 4);
+
+    ptFlowPar->dKappa = format_double(calculateKappa(ptNGPar->dZf), 4);
+    ptFlowPar->dDViscosity = format_double(Dlndjs(ptFlowPar->dPf / 1e6, ptFlowPar->dTf), 5);
+    ptFlowPar->dDExpCoefficient = format_double(calculateEpsilon(ptFlowPar->dPf, ptFlowPar->dDp,
+                                                                 ptFlowPar->dBeta, ptFlowPar->dKappa), 4);
+
+    double D = ptFlowPar->dPipeD / 1000.0;
+    double d = ptFlowPar->dOrificeD / 1000.0;
+    double beta = ptFlowPar->dBeta;
+    double P1 = ptFlowPar->dPf;
+    double deltaP = ptFlowPar->dDp;
+    double Tf = ptFlowPar->dTf;
+
+
+    double C_initial = 0.6;
+    double Qf_initial = (C_initial * ptFlowPar->dE * ptFlowPar->dDExpCoefficient * M_PI * pow(d, 2) / 4)
+                        * sqrt(2 * deltaP / (ptNGPar->dRhof * (1 - pow(beta, 4))));
+    ptFlowPar->dVFlowf = Qf_initial;
+    double tolerance = 1e-6;
+    double currentC = C_initial;
+    double currentReD = calculateReD(Qf_initial, D, ptNGPar->dRhof, ptFlowPar->dDViscosity);
+    int iter = 0;
+    double prevC = 0;
+    do {
+        int maxIter = 100;
+        prevC = currentC;
+        currentC = calculateCd(beta, currentReD, ptFlowPar->dPipeD, ptFlowPar->dPtmode);
+        double Qf = (currentC * ptFlowPar->dDExpCoefficient * M_PI * pow(d, 2) / 4)
+                    * sqrt(2 * deltaP / (ptNGPar->dRhof * (1 - pow(beta, 4))));
+        ptFlowPar->dVFlowf = Qf;
+        currentReD = calculateReD(Qf, D, ptNGPar->dRhof, ptFlowPar->dDViscosity);
+        iter++;
+        if (iter > maxIter) {
+            fprintf(stderr, "\n");
+        }
+    } while (fabs(currentC - prevC) / currentC > tolerance);
+
+    double K = calculateK(ptFlowPar->dPipeType);
+    double G_me = calculateRoughnessFactor(ptFlowPar->dPipeD, K, currentC);
+    double C_corrected = currentC * G_me;
+
+    ptFlowPar->dCd = C_corrected;
+    ptFlowPar->dRoughNessPipe = G_me;
+    ptFlowPar->dRnPipe = currentReD;
+    double Qn = ptFlowPar->dVFlowf * (ptFlowPar->dFpv * ptFlowPar->dFpv * P1 / ptFlowPar->dPb_M)
+                * (ptFlowPar->dTb_M) / Tf;
+    ptFlowPar->dVFlowb = Qn;
+    ptFlowPar->dMFlowb = ptFlowPar->dVFlowb * ptNGPar->dRhob;
+    ptFlowPar->dEFlowb = ptFlowPar->dVFlowb * ptNGPar->dHhvMol * ptFlowPar->dPb_M * 1e-6 / RGASKJ / ptFlowPar->dTb_M;
+    ptFlowPar->dVelocityFlow = ptFlowPar->dVFlowf / (M_PI * pow((ptFlowPar->dPipeD / 2000), 2));
+}
+
+
+double CaiLiaoPzxs(const int tempCaiLiao) {
+    double CaiLiaoPzxs = 0;
+    switch (tempCaiLiao) {
+        case 0:
+            CaiLiaoPzxs = 11.75;
+            break;
+
+        case 1:
+            CaiLiaoPzxs = 11.6;
+            break;
+
+        case 2:
+            CaiLiaoPzxs = 11.16;
+            break;
+
+        case 3:
+            CaiLiaoPzxs = 11.59;
+            break;
+
+        case 4:
+            CaiLiaoPzxs = 10.5;
+            break;
+
+        case 5:
+            CaiLiaoPzxs = 10.0;
+            break;
+
+        case 6:
+            CaiLiaoPzxs = 10.2;
+            break;
+
+        case 7:
+            CaiLiaoPzxs = 15.5;
+            break;
+
+        case 8:
+            CaiLiaoPzxs = 11.5;
+            break;
+
+        case 9:
+            CaiLiaoPzxs = 10.8;
+            break;
+
+        case 10:
+            CaiLiaoPzxs = 16.6;
+            break;
+
+        case 11:
+            CaiLiaoPzxs = 11.4;
+            break;
+
+        case 12:
+            CaiLiaoPzxs = 16.55;
+            break;
+
+        case 13:
+            CaiLiaoPzxs = 17.8;
+            break;
+
+        case 14:
+            CaiLiaoPzxs = 17.2;
+            break;
+        default: ;
+    }
+    return CaiLiaoPzxs;
+}
+
+
+double calculateK(int dPipeType) {
+    double Jdccd;
+    switch (dPipeType) {
+        case 0:
+            Jdccd = 0.029;
+            break;
+        case 1:
+        case 2:
+        case 3:
+            Jdccd = 0.075;
+            break;
+        case 4:
+            Jdccd = 0.1;
+            break;
+        case 5:
+            Jdccd = 0.15;
+            break;
+        case 6:
+            Jdccd = 1;
+            break;
+        case 7:
+            Jdccd = 2.1;
+            break;
+        case 8:
+            Jdccd = 0.04;
+            break;
+        case 9:
+            Jdccd = 0.15;
+            break;
+        case 10:
+            Jdccd = 0.13;
+            break;
+        case 11:
+            Jdccd = 0.25;
+            break;
+        default:
+
+            fprintf(stderr, "δ֪<EFBFBD>Ĺܵ<EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>: %d\n", dPipeType);
+            return FLOW_CALC_ERROR;
+    }
+    return Jdccd;
+}
+
+
+double calculateRoughnessFactor(double D_pipe, double K, double C) {
+    double K_over_D = K / D_pipe;
+
+
+    if (K_over_D <= 0.0004) {
+        return 1.0000;
+    }
+
+
+    double term = (K_over_D * 1e6) - 400;
+    if (term < 0) {
+        fprintf(
+            stderr,
+            "K/D <20><><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>ʽ<EFBFBD><CABD><EFBFBD>÷<EFBFBD>Χ\n");
+        return FLOW_CALC_ERROR;
+    }
+
+    double G_me = 1 + (0.011 / C) * sqrt(term);
+    return G_me;
+}
+
+
+double calculateE(double beta) {
+    return 1 / sqrt(1 - pow(beta, 4));
+}
+
+
+double calculateFG(double dRD_Real) {
+    return 1 / sqrt(dRD_Real);
+}
+
+
+double calculateFT(double dTb_M, double dTf) {
+    return sqrt(dTb_M / dTf);
+}
+
+
+double calculateEpsilon(double dPf, double dDp, double beta, double dKappa) {
+    double tau = (dPf - dDp) / dPf;
+    double epsilon = 1 - (0.351 + 0.256 * pow(beta, 4) + 0.93 * pow(beta, 8)) * (1 - pow(tau, 1 / dKappa));
+    return epsilon;
+}
+
+
+double calculateKappa(double dZf) {
+    double gamma = 1.3;
+    double Z = dZf;
+
+
+    double kappa = gamma / (1 - (gamma - 1) * (1 / Z - 1));
+    return kappa;
+}
+
+
+double calculateReD(double Qf, double D, double rho, double mu) {
+    return (4 * Qf * rho * 1000) / (M_PI * D * mu);
+}
+
+
+double calculateCd(double beta, double ReD, double D_mm, int ptMode) {
+    double L1, L2;
+
+    switch (ptMode) {
+        case 1:
+            L1 = L2 = 0;
+            break;
+        case 0:
+            L1 = L2 = 25.4 / D_mm;
+            break;
+        case 2:
+            L1 = 1.0;
+            L2 = 0.47;
+            break;
+        default:
+            fprintf(stderr, "<EFBFBD><EFBFBD>֧<EFBFBD>ֵ<EFBFBD>ȡѹ<EFBFBD><EFBFBD>ʽ: %d\n", ptMode);
+            return FLOW_CALC_ERROR;
+    }
+
+    double term1 = 0.5961 + 0.0261 * pow(beta, 2) - 0.216 * pow(beta, 8);
+    double term2 = 0.000521 * pow(1e6 * beta / ReD, 0.7);
+    double A = pow(19000 * beta / ReD, 0.8);
+    double term3 = (0.0188 + 0.0063 * A) * pow(beta, 3.5) * pow(1e6 / ReD, 0.3);
+    double term4 = (0.043 + 0.08 * exp(-10 * L1) - 0.123 * exp(-7 * L1))
+                   * (1 - 0.11 * A) * pow(beta, 4) / (1 - pow(beta, 4));
+    double term5 = -0.031 * (2 * L2 / (1 - beta) - 0.8 * pow(2 * L2 / (1 - beta), 1.1))
+                   * pow(beta, 1.3);
+
+    double Cd = term1 + term2 + term3 + term4 + term5;
+
+
+    if (D_mm < 71.12) {
+        Cd += 0.011 * (0.75 - beta) * (2.8 - D_mm / 25.4);
+    }
+    return Cd;
+}
+
+
+double Dlndjs(double tempP_jy, double tempT) {
+    double Dlndjs_Dlnd_Data[8][11] = {
+        {976, 991, 1014, 1044, 1073, 1114, 1156, 1207, 1261, 1331, 1405},
+        {1027, 1040, 1063, 1091, 1118, 1151, 1185, 1230, 1276, 1331, 1389},
+        {1071, 1082, 1106, 1127, 1149, 1180, 1211, 1250, 1289, 1335, 1383},
+        {1123, 1135, 1153, 1174, 1195, 1224, 1253, 1289, 1324, 1366, 1409},
+        {1167, 1178, 1196, 1216, 1236, 1261, 1287, 1318, 1350, 1385, 1421},
+        {1213, 1224, 1239, 1257, 1275, 1297, 1320, 1346, 1373, 1403, 1435},
+        {1260, 1270, 1281, 1297, 1313, 1333, 1352, 1374, 1396, 1424, 1451},
+        {1303, 1312, 1323, 1338, 1352, 1372, 1391, 1412, 1432, 1456, 1482}
+    };
+
+    double Dlndjs_Dlnd_T[8] = {
+        -15 + 273.15, 0 + 273.15, 15 + 273.15, 30 + 273.15,
+        45 + 273.15, 60 + 273.15, 75 + 273.15, 90 + 273.15
+    };
+
+    double Dlndjs_Dlnd_P[11] = {0.1, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10};
+
+    double ky, kx;
+    int i, m = 0, n = 0;
+
+    if (tempT < Dlndjs_Dlnd_T[0]) {
+        tempT = Dlndjs_Dlnd_T[0];
+    }
+    if (tempT > Dlndjs_Dlnd_T[7]) {
+        tempT = Dlndjs_Dlnd_T[7];
+    }
+    if (tempP_jy < Dlndjs_Dlnd_P[0]) {
+        tempP_jy = Dlndjs_Dlnd_P[0];
+    }
+    if (tempP_jy > Dlndjs_Dlnd_P[10]) {
+        tempP_jy = Dlndjs_Dlnd_P[10];
+    }
+
+    for (i = 0; i <= 6; i++) {
+        if (tempT >= Dlndjs_Dlnd_T[i] && tempT <= Dlndjs_Dlnd_T[i + 1]) {
+            m = i;
+            break;
+        }
+    }
+
+    for (i = 0; i <= 9; i++) {
+        if (tempP_jy >= Dlndjs_Dlnd_P[i] && tempP_jy <= Dlndjs_Dlnd_P[i + 1]) {
+            n = i;
+            break;
+        }
+    }
+
+    if (Dlndjs_Dlnd_P[n + 1] - Dlndjs_Dlnd_P[n] != 0) {
+        ky = (tempP_jy - Dlndjs_Dlnd_P[n]) / (Dlndjs_Dlnd_P[n + 1] - Dlndjs_Dlnd_P[n]);
+    } else {
+        ky = 0;
+    }
+
+    if (Dlndjs_Dlnd_T[m + 1] - Dlndjs_Dlnd_T[m] != 0) {
+        kx = (tempT - Dlndjs_Dlnd_T[m]) / (Dlndjs_Dlnd_T[m + 1] - Dlndjs_Dlnd_T[m]);
+    } else {
+        kx = 0;
+    }
+
+    double s1 = Dlndjs_Dlnd_Data[m][n] + (Dlndjs_Dlnd_Data[m][n + 1] - Dlndjs_Dlnd_Data[m][n]) * ky;
+    double s2 = Dlndjs_Dlnd_Data[m + 1][n] + (Dlndjs_Dlnd_Data[m + 1][n + 1] - Dlndjs_Dlnd_Data[m + 1][n]) * ky;
+    return (s1 + (s2 - s1) * kx) / 100000.0;
+}
+
+// 压力损失计算
+double YaLiSunShi(double tempLiuChuXiShu, double tempZjb, double tempDp, int JieLiuZhuangZhi) {
+    double ylss = 0;
+    switch (JieLiuZhuangZhi) {
+        case 0:
+        case 1:
+        case 2:
+            ylss = tempDp * (sqrt(1 - tempZjb) - tempLiuChuXiShu * pow(tempZjb, 2))
+                   / (sqrt(1 - tempZjb) + tempLiuChuXiShu * pow(tempZjb, 2));
+            break;
+        default: ;
+    }
+    return ylss;
+}
+
+
+// 标况转工况流量转换
+double FlowConvert_BaseToWork(double Pf, double Tf, double Zb, double Zf, double FlowBase, int Cbtj) {
+    double tempPn = 0;
+    double tempTn = 0;
+
+
+    switch (Cbtj) {
+        case 2:
+            tempPn = 101325;
+            tempTn = 273.15;
+            break;
+        case 1:
+            tempPn = 101325;
+            tempTn = 288.15;
+            break;
+        case 0:
+            tempPn = 101325;
+            tempTn = 293.15;
+            break;
+        case 3:
+            tempPn = 10155981;
+            tempTn = 288.7055555;
+            break;
+        default: ;
+    }
+    return FlowBase * tempPn * Tf * Zf / Pf / tempTn / Zb;
+}
+
+// 工况转标况流量转换
+double FlowConvert_WorkToBase(double Pf, double Tf, double Zb, double Zf, double FlowWork, int Cbtj) {
+    double tempPn = 0;
+    double tempTn = 0;
+
+    switch (Cbtj) {
+        case 2:
+            tempPn = 101325;
+            tempTn = 273.15;
+            break;
+        case 1:
+            tempPn = 101325;
+            tempTn = 288.15;
+            break;
+        case 0:
+            tempPn = 101325;
+            tempTn = 293.15;
+            break;
+        case 3:
+            tempPn = 10155981;
+            tempTn = 288.7055555;
+            break;
+        default: ;
+    }
+    return FlowWork * Pf * tempTn * Zb / tempPn / Tf / Zf;
+}

NGFlowMeter/CppProject/NGTools/FlowCal.h

@@ -0,0 +1,87 @@
+#ifndef FLOWCAL_H
+#define FLOWCAL_H
+#include "NGCal.h"
+typedef struct FlowParSTRUCT {
+
+
+	int dCbtj;
+	double dPb_M;
+	double dTb_M;
+	double dPb_E;
+	double dTb_E;
+	double dPatm;
+	double dNG_Compents[21];
+
+	int dMeterType;
+	int dCoreType;
+	int dPtmode;
+	int dPipeType;
+	double dPipeD;
+	int dPipeMaterial;
+
+	double dOrificeD;
+	int dOrificeMaterial;
+
+	double dPf;
+	int dPfType;
+	double dTf;
+	double dDp;
+
+	double dMeterFactor;
+	double dPulseNum;
+
+
+
+	double dE;
+	double dFG;
+	double dFT;
+	double dDViscosity;
+	double dDExpCoefficient;
+	double dRnPipe;
+	double dBk;
+	double dRoughNessPipe;
+	double dCd;
+	double dCdCorrect;
+	double dCdNozell;
+	double dVFlowb;
+	double dVFlowf;
+	double dMFlowb;
+	double dEFlowb;
+	double dVelocityFlow;
+	double dPressLost;
+	double dBeta;
+	double dKappa;
+	double dFpv;
+} FlowParSTRUCT;
+
+double CaiLiaoPzxs(int tempCaiLiao);
+
+double calculateK(int dPipeType);
+
+double calculateRoughnessFactor(double D_pipe, double K, double C);
+
+void thermalExpansionCorrection(double dOrificeMaterial, double dOrificeD,
+								double dPipeMaterial, double dPipeD,
+								double dTf, double correctedValues[3]);
+
+double calculateE(double beta);
+
+double calculateFG(double dRD_Real);
+
+double calculateFT(double dTb_M, double dTf);
+
+double calculateEpsilon(double dPf, double dDp, double beta, double dKappa);
+
+double calculateKappa(double dZf);
+
+double calculateReD(double Qf, double D, double rho, double mu);
+
+double calculateCd(double beta, double ReD, double D_mm, int ptMode);
+
+
+double Dlndjs(double tempP_jy, double tempT);
+
+void OFlowCal(FlowParSTRUCT *ptFlowPar, NGParSTRUCT *ptNGPar);
+
+#endif
+

NGFlowMeter/CppProject/NGTools/NGCal.c

@@ -0,0 +1,87 @@
+#include "NGCal.h"
+#include "Therm.h"
+#include "Detail.h"
+#include "math.h"
+static Therm *ptTherm;
+static Detail *ptDetail;
+
+int NGCal_Init(NGParSTRUCT *ptNGPar) {
+
+	ptDetail = Detail_Construct();
+	if (NULL == ptDetail) {
+		return MEMORY_ALLOCATION_ERROR;
+	}
+
+
+	ptTherm = (Therm *) malloc(sizeof(Therm));
+	Therm_Init(ptTherm);
+	if (NULL == ptTherm) {
+		return MEMORY_ALLOCATION_ERROR;
+	}
+
+	return NGCal_NGCal;
+}
+
+int NGCal_UnInit(void) {
+	if (ptDetail) free(ptDetail);
+	if (ptTherm) free(ptTherm);
+	return 0;
+}
+
+double SOS(NGParSTRUCT *ptNGPar) {
+	if (NULL == ptDetail || NULL == ptTherm) {
+		NGCal_UnInit();
+		NGCal_Init(ptNGPar);
+	}
+	Therm_Run(ptTherm, ptNGPar, ptDetail);
+	ptNGPar->dCstar = 0.0;
+	return ptNGPar->dSOS;
+}
+
+double Crit(NGParSTRUCT *ptNGPar, double dPlenumVelocity) {
+	if (NULL == ptDetail || NULL == ptTherm) {
+		NGCal_UnInit();
+		if (NGCal_NGCal != NGCal_Init(ptNGPar)) {
+			ptNGPar->lStatus = MEMORY_ALLOCATION_ERROR;
+			return 0.0;
+		}
+	}
+
+
+	Therm_Run(ptTherm, ptNGPar, ptDetail);
+
+	double DH = (ptNGPar->dSOS * ptNGPar->dSOS - dPlenumVelocity * dPlenumVelocity) / 2.0;
+	double S = ptNGPar->dS;
+	double H = ptNGPar->dH;
+	double R = ptNGPar->dRhof;
+	double P = ptNGPar->dPf;
+	double Z = ptNGPar->dZf;
+	double T = ptNGPar->dTf;
+	// DDH = 10.0;
+	int i=0;
+	for ( i = 1; i < MAX_NUM_OF_ITERATIONS; i++) {
+		double tolerance = 1.0;
+		Therm_HS_Mode(ptTherm, ptNGPar, ptDetail, H - DH, S, 1);
+		Therm_Run(ptTherm, ptNGPar, ptDetail);
+		double DDH = DH;
+		DH = (ptNGPar->dSOS * ptNGPar->dSOS - dPlenumVelocity * dPlenumVelocity) / 2.0;
+		if (fabs(DDH - DH) < tolerance) break;
+	}
+	ptNGPar->dCstar = (ptNGPar->dRhof * ptNGPar->dSOS) / sqrt(R * P * Z);
+	ptNGPar->dPf = P;
+	ptNGPar->dTf = T;
+	Therm_Run(ptTherm, ptNGPar, ptDetail);
+	Detail_dhvMol(ptDetail,ptNGPar);
+	return ptNGPar->dCstar;
+}
+
+double Cperf(const NGParSTRUCT *ptNGPar) {
+	double k = ptNGPar->dKappa;
+	double root = 2.0 / (k + 1.0);
+	double exponent = (k + 1.0) / (k - 1.0);
+	return (sqrt(k * pow(root, exponent)));
+}
+
+double CRi(const NGParSTRUCT *ptNGPar) {
+	return (Cperf(ptNGPar) / sqrt(ptNGPar->dZf));
+}

NGFlowMeter/CppProject/NGTools/NGCal.h

@@ -0,0 +1,80 @@
+/*************************************************************************
+* <EFBFBD>ļ<EFBFBD>:	NGCal.h
+**************************************************************************/
+#ifndef _NGCal_H
+#define _NGCal_H
+#include <stdio.h>
+#include <stdlib.h>
+#define NORMAL              9000
+#define NGCal_NGCal         9001
+#define MEMORY_ALLOCATION_ERROR      9002
+#define GENERAL_CALCULATION_FAILURE  9003
+#define MAX_NUM_OF_ITERATIONS_EXCEEDED 9004
+#define NEGATIVE_DENSITY_DERIVATIVE 9005
+#define MAX_DENSITY_IN_BRAKET_EXCEEDED 9006
+#define FLOW_CALC_ERROR 9007
+#define FLOW_CALC_DIEDAI_ERROR 9008
+#define NUMBEROFCOMPONENTS 21
+#define M_PI 3.1415926535897932
+#define MAX_NUM_OF_ITERATIONS 100
+#define P_CHG_TOL 0.001
+#define T_CHG_TOL 0.001
+#define P_MAX    1.379e8
+#define P_MIN    0.0
+#define T_MAX    473.15
+#define T_MIN    143.0
+#define RGASKJ  8.314510e-3
+#define RGAS    8.314510
+typedef struct tagNGParSTRUCT
+{
+	long lStatus;
+	int bForceUpdate;
+	double adMixture[21];
+	int dCbtj;
+	double dPb;
+	double dTb;
+	double dPf;
+	double dTf;
+	double dMrx;
+	double dZb;
+	double dZf;
+	double dFpv;
+	double dDb;
+	double dDf;
+	double dRhob;
+	double dRhof;
+	double dRD_Ideal;
+	double dRD_Real;
+	double dHo;
+	double dH;
+	double dS;
+	double dCpi;
+	double dCp;
+	double dCv;
+	double dk;
+	double dKappa;
+	double dSOS;
+	double dCstar;
+	double dHhvMol;
+	double dLhvMol;
+} NGParSTRUCT;
+
+
+enum gascomp {
+    XiC1=0, XiN2, XiCO2, XiC2, XiC3,
+    XiH2O, XiH2S, XiH2, XiCO, XiO2,
+    XiIC4, XiNC4, XiIC5, XiNC5, XiNC6,
+    XiNC7, XiNC8, XiNC9, XiNC10, XiHe, XiAr
+};
+
+
+ int NGCal_Init(NGParSTRUCT * ptNGPar);
+ int NGCal_UnInit(void);
+
+
+ double SOS(NGParSTRUCT *);
+
+
+ double Crit(NGParSTRUCT *, double);
+
+#endif

NGFlowMeter/CppProject/NGTools/Therm.c

@@ -0,0 +1,413 @@
+#include "therm.h"
+#include <math.h>
+#include "Detail.h" 
+
+
+void Therm_Init(Therm *therm) {
+    therm->CAL_TH = 4.1840;
+    therm->coefA = 0;
+    therm->coefB = 1;
+    therm->coefC = 2;
+    therm->coefD = 3;
+    therm->coefE = 4;
+    therm->coefF = 5;
+    therm->coefG = 6;
+    therm->coefH = 7;
+    therm->coefI = 8;
+    therm->coefJ = 9;
+    therm->coefK = 10;
+
+    therm->dPdD = 0.0;
+    therm->dPdT = 0.0;
+    therm->dSi = 0.0;
+    therm->dTold = 0.0;
+    therm->dMrxold = 0.0;
+
+    therm->GK_points = 5;
+
+
+    therm->GK_root[0] = 0.14887433898163121088;
+    therm->GK_root[1] = 0.43339539412924719080;
+    therm->GK_root[2] = 0.67940956829902440263;
+    therm->GK_root[3] = 0.86506336668898451073;
+    therm->GK_root[4] = 0.97390652851717172008;
+
+
+    therm->GK_weight[0] = 0.29552422471475286217;
+    therm->GK_weight[1] = 0.26926671930999634918;
+    therm->GK_weight[2] = 0.21908636251598204295;
+    therm->GK_weight[3] = 0.14945134915058059038;
+    therm->GK_weight[4] = 0.066671344308688137179;
+
+
+    double thermConstants[21][11] = {
+        {-29776.4, 7.95454, 43.9417, 1037.09, 1.56373, 813.205, -24.9027, 1019.98, -10.1601, 1070.14, -20.0615},
+        {-3495.34, 6.95587, 0.272892, 662.738, -0.291318, -680.562, 1.78980, 1740.06, 0.0, 100.0, 4.49823},
+        {20.7307, 6.96237, 2.68645, 500.371, -2.56429, -530.443, 3.91921, 500.198, 2.13290, 2197.22, 5.81381},
+        {-37524.4, 7.98139, 24.3668, 752.320, 3.53990, 272.846, 8.44724, 1020.13, -13.2732, 869.510, -22.4010},
+        {-56072.1, 8.14319, 37.0629, 735.402, 9.38159, 247.190, 13.4556, 1454.78, -11.7342, 984.518, -24.0426},
+        {-13773.1, 7.97183, 6.27078, 2572.63, 2.05010, 1156.72, 0.0, 100.0, 0.0, 100.0, -3.24989},
+        {-10085.4, 7.94680, -0.08380, 433.801, 2.85539, 843.792, 6.31595, 1481.43, -2.88457, 1102.23, -0.51551},
+        {-5565.60, 6.66789, 2.33458, 2584.98, 0.749019, 559.656, 0.0, 100.0, 0.0, 100.0, -7.94821},
+        {-2753.49, 6.95854, 2.02441, 1541.22, 0.096774, 3674.81, 0.0, 100.0, 0.0, 100.0, 6.23387},
+        {-3497.45, 6.96302, 2.40013, 2522.05, 2.21752, 1154.15, 0.0, 100.0, 0.0, 100.0, 9.19749},
+        {-72387.0, 17.8143, 58.2062, 1787.39, 40.7621, 808.645, 0.0, 100.0, 0.0, 100.0, -44.1341},
+        {-72674.8, 18.6383, 57.4178, 1792.73, 38.6599, 814.151, 0.0, 100.0, 0.0, 100.0, -46.1938},
+        {-91505.5, 21.3861, 74.3410, 1701.58, 47.0587, 775.899, 0.0, 100.0, 0.0, 100.0, -60.2474},
+        {-83845.2, 22.5012, 69.5789, 1719.58, 46.2164, 802.174, 0.0, 100.0, 0.0, 100.0, -62.2197},
+        {-94982.5, 26.6225, 80.3819, 1718.49, 55.6598, 802.069, 0.0, 100.0, 0.0, 100.0, -77.5366},
+        {-103353.0, 30.4029, 90.6941, 1669.32, 63.2028, 786.001, 0.0, 100.0, 0.0, 100.0, -92.0164},
+        {-109674.0, 34.0847, 100.253, 1611.55, 69.7675, 768.847, 0.0, 100.0, 0.0, 100.0, -106.149},
+        {-122599.0, 38.5014, 111.446, 1646.48, 80.5015, 781.588, 0.0, 100.0, 0.0, 100.0, -122.444},
+        {-133564.0, 42.7143, 122.173, 1654.85, 90.2255, 785.564, 0.0, 100.0, 0.0, 100.0, -138.006},
+        {0.0, 4.9680, 0.0, 100.0, 0.0, 100.0, 0.0, 100.0, 0.0, 100.0, 0.0},
+        {0.0, 4.9680, 0.0, 100.0, 0.0, 100.0, 0.0, 100.0, 0.0, 100.0, 0.0}
+    };
+
+int i = 0;
+int j=0;
+    for (i = 0; i < 21; i++) {
+        for ( j = 0; j < 11; j++) {
+            therm->ThermConstants[i][j] = thermConstants[i][j];
+        }
+    }
+}
+
+
+void Therm_Run(Therm *therm, NGParSTRUCT *ptNGPar, Detail *detail) {
+
+    double c, x, y, z;
+    Detail_Run(detail, ptNGPar);
+    Detail_dZdD(detail, ptNGPar->dDf);
+    Therm_CprCvrHS(therm, ptNGPar, detail);
+    ptNGPar->dk = ptNGPar->dCp / ptNGPar->dCv;
+    x = ptNGPar->dk * RGAS * 1000.0 * ptNGPar->dTf;
+    y = ptNGPar->dMrx;
+    z = ptNGPar->dZf + ptNGPar->dDf * detail->ddZdD;
+    c = (x / y) * z;
+    ptNGPar->dSOS = sqrt(c);
+    ptNGPar->dKappa = (c * ptNGPar->dRhof) / ptNGPar->dPf;
+}
+
+
+double Therm_CpiMolar(Therm *therm, NGParSTRUCT *ptNGPar) {
+    double cp = 0.0;
+    double Cpx;
+    double DT, FT, HT, JT;
+    double Dx, Fx, Hx, Jx;
+    double T;
+
+    T = ptNGPar->dTf;
+int i = 0;
+    for ( i = 0; i < NUMBEROFCOMPONENTS; i++) {
+        if (ptNGPar->adMixture[i] > 0) {
+            Cpx = 0.0;
+            DT = therm->ThermConstants[i][therm->coefD] / T;
+            FT = therm->ThermConstants[i][therm->coefF] / T;
+            HT = therm->ThermConstants[i][therm->coefH] / T;
+            JT = therm->ThermConstants[i][therm->coefJ] / T;
+            Dx = DT / sinh(DT);
+            Fx = FT / cosh(FT);
+            Hx = HT / sinh(HT);
+            Jx = JT / cosh(JT);
+            Cpx += therm->ThermConstants[i][therm->coefB];
+            Cpx += therm->ThermConstants[i][therm->coefC] * Dx * Dx;
+            Cpx += therm->ThermConstants[i][therm->coefE] * Fx * Fx;
+            Cpx += therm->ThermConstants[i][therm->coefG] * Hx * Hx;
+            Cpx += therm->ThermConstants[i][therm->coefI] * Jx * Jx;
+            Cpx *= ptNGPar->adMixture[i];
+            cp += Cpx;
+        }
+    }
+    cp *= therm->CAL_TH;
+    return cp;
+}
+
+
+double Therm_coth(double x) {
+    return 1.0 / tanh(x);
+}
+
+
+double Therm_Ho(Therm *therm, NGParSTRUCT *ptNGPar) {
+    double H = 0.0;
+    double Hx;
+    double DT, FT, HT, JT;
+    double cothDT, tanhFT, cothHT, tanhJT;
+    double T = ptNGPar->dTf;
+int i = 0;
+    for (i = 0; i < NUMBEROFCOMPONENTS; i++) {
+        if (ptNGPar->adMixture[i] <= 0.0) continue;
+        Hx = 0.0;
+        DT = therm->ThermConstants[i][therm->coefD] / T;
+        FT = therm->ThermConstants[i][therm->coefF] / T;
+        HT = therm->ThermConstants[i][therm->coefH] / T;
+        JT = therm->ThermConstants[i][therm->coefJ] / T;
+        cothDT = Therm_coth(DT);
+        tanhFT = tanh(FT);
+        cothHT = Therm_coth(HT);
+        tanhJT = tanh(JT);
+        Hx += therm->ThermConstants[i][therm->coefA];
+        Hx += therm->ThermConstants[i][therm->coefB] * T;
+        Hx += therm->ThermConstants[i][therm->coefC] * therm->ThermConstants[i][therm->coefD] * cothDT;
+        Hx -= therm->ThermConstants[i][therm->coefE] * therm->ThermConstants[i][therm->coefF] * tanhFT;
+        Hx += therm->ThermConstants[i][therm->coefG] * therm->ThermConstants[i][therm->coefH] * cothHT;
+        Hx -= therm->ThermConstants[i][therm->coefI] * therm->ThermConstants[i][therm->coefJ] * tanhJT;
+        Hx *= ptNGPar->adMixture[i];
+        H += Hx;
+    }
+
+    H *= therm->CAL_TH;
+    H /= ptNGPar->dMrx;
+    return H * 1000.0;
+}
+
+
+double Therm_So(Therm *therm, NGParSTRUCT *ptNGPar) {
+    double S = 0.0;
+    double Sx;
+    double DT, FT, HT, JT;
+    double cothDT, tanhFT, cothHT, tanhJT;
+    double sinhDT, coshFT, sinhHT, coshJT;
+    double T = ptNGPar->dTf;
+    int i = 0;
+    for (i = 0; i < NUMBEROFCOMPONENTS; i++) {
+        if (ptNGPar->adMixture[i] <= 0.0) continue;
+        Sx = 0.0;
+        DT = therm->ThermConstants[i][therm->coefD] / T;
+        FT = therm->ThermConstants[i][therm->coefF] / T;
+        HT = therm->ThermConstants[i][therm->coefH] / T;
+        JT = therm->ThermConstants[i][therm->coefJ] / T;
+        cothDT = Therm_coth(DT);
+        tanhFT = tanh(FT);
+        cothHT = Therm_coth(HT);
+        tanhJT = tanh(JT);
+        sinhDT = sinh(DT);
+        coshFT = cosh(FT);
+        sinhHT = sinh(HT);
+        coshJT = cosh(JT);
+        Sx += therm->ThermConstants[i][therm->coefK];
+        Sx += therm->ThermConstants[i][therm->coefB] * log(T);
+        Sx += therm->ThermConstants[i][therm->coefC] * (DT * cothDT - log(sinhDT));
+        Sx -= therm->ThermConstants[i][therm->coefE] * (FT * tanhFT - log(coshFT));
+        Sx += therm->ThermConstants[i][therm->coefG] * (HT * cothHT - log(sinhHT));
+        Sx -= therm->ThermConstants[i][therm->coefI] * (JT * tanhJT - log(coshJT));
+
+        Sx *= ptNGPar->adMixture[i];
+        S += Sx;
+    }
+
+    S *= therm->CAL_TH;
+    S /= ptNGPar->dMrx;
+    return S * 1000.0;
+}
+
+
+void Therm_CprCvrHS(Therm *therm, NGParSTRUCT *ptNGPar, Detail *detail) {
+    double Cvinc = 0.0;
+    double Cvr, Cpr;
+    double Hinc = 0.0;
+    double Sinc = 0.0;
+    double Smixing = 0.0;
+   double Si;
+
+    double Cp = Therm_CpiMolar(therm, ptNGPar);
+    ptNGPar->dHo = Therm_Ho(therm, ptNGPar);
+    Si = Therm_So(therm, ptNGPar);
+    ptNGPar->dCpi = (Cp * 1000.0) / ptNGPar->dMrx;
+     int i = 0;
+    for ( i = 0; i < therm->GK_points; i++) {
+        double x = ptNGPar->dDf * (1.0 + therm->GK_root[i]) / 2.0;
+        Detail_zdetail(detail, x);
+        Detail_dZdT(detail, x);
+        Detail_d2ZdT2(detail, x);
+        Hinc += therm->GK_weight[i] * detail->ddZdT / x;
+        Cvinc += therm->GK_weight[i] * (2.0 * detail->ddZdT + ptNGPar->dTf * detail->dd2ZdT2) / x;
+        Sinc += therm->GK_weight[i] * (detail->dZ + ptNGPar->dTf * detail->ddZdT - 1.0) / x;
+        x = ptNGPar->dDf * (1.0 - therm->GK_root[i]) / 2.0;
+        Detail_zdetail(detail, x);
+        Detail_dZdT(detail, x);
+        Detail_d2ZdT2(detail, x);
+        Hinc += therm->GK_weight[i] * detail->ddZdT / x;
+        Cvinc += therm->GK_weight[i] * (2.0 * detail->ddZdT + ptNGPar->dTf * detail->dd2ZdT2) / x;
+        Sinc += therm->GK_weight[i] * (detail->dZ + ptNGPar->dTf * detail->ddZdT - 1.0) / x;
+    }
+
+    Detail_zdetail(detail, ptNGPar->dDf);
+    Detail_dZdT(detail, ptNGPar->dDf);
+    Detail_d2ZdT2(detail, ptNGPar->dDf);
+    Cvr = Cp - RGAS * (1.0 + ptNGPar->dTf * Cvinc * 0.5 * ptNGPar->dDf);
+    double a = (ptNGPar->dZf + ptNGPar->dTf * detail->ddZdT);
+    double b = (ptNGPar->dZf + ptNGPar->dDf * detail->ddZdD);
+    double dPdT = RGAS * ptNGPar->dDf * a;
+    double dPdD = RGAS * ptNGPar->dTf * b;
+    Cpr = Cvr + RGAS * ((a * a) / b);
+
+    Cpr /= ptNGPar->dMrx;
+    Cvr /= ptNGPar->dMrx;
+
+    ptNGPar->dCv = Cvr * 1000.0;
+    ptNGPar->dCp = Cpr * 1000.0;
+
+    ptNGPar->dH = ptNGPar->dHo + 1000.0 * RGAS * ptNGPar->dTf * (
+                      ptNGPar->dZf - 1.0 - ptNGPar->dTf * Hinc * 0.5 * ptNGPar->dDf) / ptNGPar->dMrx;
+    for ( i = 0; i < NUMBEROFCOMPONENTS; i++) {
+        if (ptNGPar->adMixture[i] != 0) Smixing += ptNGPar->adMixture[i] * log(ptNGPar->adMixture[i]);
+    }
+    Smixing *= RGAS;
+
+    ptNGPar->dS = Si - Smixing - 1000.0 * RGAS * (
+                      log(ptNGPar->dPf / 101325.0) - log(ptNGPar->dZf) + Sinc * 0.5 * ptNGPar->dDf) / ptNGPar->dMrx;
+}
+
+
+void Therm_HS_Mode(Therm *therm, NGParSTRUCT *ptNGPar, Detail *detail, double H, double S, int bGuess) {
+    double s0 = S;
+    double h0 = H;
+    double t1, t2, tmin, tmax;
+    double p1, p2, px, pmin, pmax;
+    double delta1, delta2;
+    double tolerance = 0.001;
+
+    if (bGuess) {
+        t1 = ptNGPar->dTf;
+        px = ptNGPar->dPf;
+        pmax = px * 2.0;
+        pmin = px * 0.1;
+        tmax = t1 * 1.5;
+        tmin = t1 * 0.67;
+    } else {
+        t1 = 273.15;
+        px = 1013250.0;
+        pmax = P_MAX;
+        pmin = 10000.0;
+        tmax = T_MAX;
+        tmin = T_MIN;
+    }
+
+    t2 = t1 + 10.0;
+
+    Detail_Run(detail, ptNGPar);
+    double h1 = Therm_H(therm, ptNGPar, detail) - h0;
+int i = 0;
+
+    for ( i = 0; i < MAX_NUM_OF_ITERATIONS; i++) {
+        ptNGPar->dTf = t2;
+        p1 = px;
+        p2 = px * 0.1;
+        ptNGPar->dPf = p1;
+        Detail_Run(detail, ptNGPar);
+        double s1 = Therm_S(therm, ptNGPar, detail) - s0;
+int j=0;
+        for ( j = 0; j < MAX_NUM_OF_ITERATIONS; j++) {
+            ptNGPar->dPf = p2;
+            Detail_Run(detail, ptNGPar);
+            double s2 = Therm_S(therm, ptNGPar, detail) - s0;
+            delta2 = fabs(s1 - s2) / s0;
+            if (delta2 < tolerance) break;
+            double p0 = p2;
+            p2 = (p1 * s2 - p2 * s1) / (s2 - s1);
+            if (p2 <= pmin) p2 = pmin;
+            if (p2 >= pmax) p2 = pmax;
+            p1 = p0;
+            s1 = s2;
+        }
+        if (ptNGPar->lStatus == MAX_NUM_OF_ITERATIONS_EXCEEDED) break;
+        double h2 = Therm_H(therm, ptNGPar, detail) - h0;
+        delta1 = fabs(h1 - h2) / h0;
+        if (delta1 < tolerance && i > 0) break;
+        double t0 = t2;
+        t2 = (t1 * h2 - t2 * h1) / (h2 - h1);
+        if (t2 >= tmax) t2 = tmax;
+        if (t2 <= tmin) {
+            t2 = t0 + 10.0;
+            ptNGPar->dTf = t2;
+            Detail_Run(detail, ptNGPar);
+            h2 = Therm_H(therm, ptNGPar, detail) - h0;
+        }
+
+        t1 = t0;
+        h1 = h2;
+    }
+
+    if (ptNGPar->lStatus == MAX_NUM_OF_ITERATIONS_EXCEEDED) {
+        ptNGPar->lStatus = MAX_NUM_OF_ITERATIONS_EXCEEDED;
+    }
+}
+
+
+double Therm_H(Therm *therm, NGParSTRUCT *ptNGPar, Detail *detail) {
+    double Hinc = 0.0;
+    ptNGPar->dHo = Therm_Ho(therm, ptNGPar);
+int i = 0;
+    for ( i = 0; i < therm->GK_points; i++) {
+        double x = ptNGPar->dDf * (1.0 + therm->GK_root[i]) / 2.0;
+        Detail_zdetail(detail, x);
+        Detail_dZdT(detail, x);
+        Detail_d2ZdT2(detail, x);
+        Hinc += therm->GK_weight[i] * detail->ddZdT / x;
+        if (i == 10) break;
+        x = ptNGPar->dDf * (1.0 - therm->GK_root[i]) / 2.0;
+        Detail_zdetail(detail, x);
+        Detail_dZdT(detail, x);
+        Detail_d2ZdT2(detail, x);
+
+        Hinc += therm->GK_weight[i] * detail->ddZdT / x;
+    }
+
+    Detail_zdetail(detail, ptNGPar->dDf);
+    Detail_dZdT(detail, ptNGPar->dDf);
+    Detail_d2ZdT2(detail, ptNGPar->dDf);
+
+    ptNGPar->dH = ptNGPar->dHo + 1000.0 * RGAS * ptNGPar->dTf * (
+                      ptNGPar->dZf - 1.0 - ptNGPar->dTf * Hinc * 0.5 * ptNGPar->dDf) / ptNGPar->dMrx;
+    return ptNGPar->dH;
+}
+
+
+double Therm_S(Therm *therm, NGParSTRUCT *ptNGPar, Detail *detail) {
+    double Sinc;
+    double Smixing;
+    double x;
+    int i;
+    Sinc = 0.0;
+    Smixing = 0.0;
+    for (i = 0; i < therm->GK_points; i++) {
+
+        x = ptNGPar->dDf * (1.0 + therm->GK_root[i]) / 2.0;
+        Detail_zdetail(detail, x);
+        Detail_dZdT(detail, x);
+        Detail_d2ZdT2(detail, x);
+        Sinc += therm->GK_weight[i] * (detail->dZ + ptNGPar->dTf * detail->ddZdT - 1.0) / x;
+        if (i == 10) break;
+
+        x = ptNGPar->dDf * (1.0 - therm->GK_root[i]) / 2.0;
+        Detail_zdetail(detail, x);
+        Detail_dZdT(detail, x);
+        Detail_d2ZdT2(detail, x);
+        Sinc += therm->GK_weight[i] * (detail->dZ + ptNGPar->dTf * detail->ddZdT - 1.0) / x;
+    }
+
+
+    Detail_zdetail(detail, ptNGPar->dDf);
+    Detail_dZdT(detail, ptNGPar->dDf);
+    Detail_d2ZdT2(detail, ptNGPar->dDf);
+
+
+    if (ptNGPar->dTf != therm->dTold || ptNGPar->dMrx != therm->dMrxold) {
+        therm->dSi = Therm_So(therm, ptNGPar);
+        therm->dTold = ptNGPar->dTf;
+        therm->dMrxold = ptNGPar->dMrx;
+    }
+
+
+    for (i = 0; i < NUMBEROFCOMPONENTS; i++) {
+        if (ptNGPar->adMixture[i] != 0) Smixing += ptNGPar->adMixture[i] * log(ptNGPar->adMixture[i]);
+    }
+    Smixing *= RGAS;
+
+
+    ptNGPar->dS = therm->dSi - Smixing - 1000.0 * RGAS * (
+                      log(ptNGPar->dPf / 101325.0) - log(ptNGPar->dZf) + Sinc * 0.5 * ptNGPar->dDf) / ptNGPar->dMrx;
+    return (ptNGPar->dS);
+}

NGFlowMeter/CppProject/NGTools/Therm.h

@@ -0,0 +1,60 @@
+/*************************************************************************
+* <EFBFBD>ļ<EFBFBD>: therm.h
+* <EFBFBD><EFBFBD><EFBFBD><EFBFBD>: Therm<EFBFBD><EFBFBD><EFBFBD>ͷ<EFBFBD>ļ<EFBFBD>
+**************************************************************************/
+
+#ifndef _THERM_H
+#define _THERM_H
+
+#include "NGCal.h"
+#include "Detail.h"
+
+typedef struct Therm {
+    double CAL_TH;
+    int coefA;
+    int coefB;
+    int coefC;
+    int coefD;
+    int coefE;
+    int coefF;
+    int coefG;
+    int coefH;
+    int coefI;
+    int coefJ;
+    int coefK;
+
+    double dPdD;
+    double dPdT;
+    double dSi;
+    double dTold;
+    double dMrxold;
+
+    int GK_points;
+    double GK_root[5];
+    double GK_weight[5];
+
+    double ThermConstants[21][11];
+} Therm;
+
+
+void Therm_Init(Therm *therm);
+
+void Therm_Run(Therm *therm, NGParSTRUCT *ptNGPar, Detail *detail);
+
+double Therm_CpiMolar(Therm *therm, NGParSTRUCT *ptNGPar);
+
+double Therm_coth(double x);
+
+double Therm_Ho(Therm *therm, NGParSTRUCT *ptNGPar);
+
+double Therm_So(Therm *therm, NGParSTRUCT *ptNGPar);
+
+void Therm_CprCvrHS(Therm *therm, NGParSTRUCT *ptNGPar, Detail *detail);
+
+void Therm_HS_Mode(Therm *therm, NGParSTRUCT *ptNGPar, Detail *detail, double H, double S, int bGuess);
+
+double Therm_H(Therm *therm, NGParSTRUCT *ptNGPar, Detail *detail);
+
+double Therm_S(Therm *therm, NGParSTRUCT *ptNGPar, Detail *detail);
+
+#endif

NGFlowMeter/CppProject/NGTools/dll.h

@@ -0,0 +1,17 @@
+#ifndef _DLL_H_
+#define _DLL_H_
+
+#if BUILDING_DLL
+#define DLLIMPORT __declspec(dllexport)
+#else
+#define DLLIMPORT __declspec(dllimport)
+#endif
+
+DLLIMPORT int __stdcall NGToolsCal(void);
+
+
+DLLIMPORT int __stdcall Test(void);
+
+
+
+#endif

NGFlowMeter/CppProject/NGTools/dllmain.c

@@ -0,0 +1,185 @@
+/* Replace "dll.h" with the name of your header */
+#include "dll.h"
+#include <windows.h>
+#include <stdio.h>
+#include "NGCal.h"
+#include "FlowCal.h"
+DLLIMPORT int  __stdcall Test(void)
+{
+	return 33;
+}
+DLLIMPORT int __stdcall  NGToolsCal(void)
+{
+	FlowParSTRUCT flowParams = {0};
+    NGParSTRUCT ngParams = {0};
+
+    // 设置基本参数
+    flowParams.dPatm = 0.0981; // 标准大气压(bar)
+    flowParams.dPf = 1.48; // 压力(MPa)
+    flowParams.dPfType = 0; // 0=表压，1=绝压
+    flowParams.dDp = 12.50; // 差压(kPa)
+    flowParams.dTf = 15; // 温度(°C)
+    flowParams.dCbtj = 0; // 参比条件类型(0=标准状态)
+
+    // 设置管道参数
+    flowParams.dPipeD = 259.38; // 管道内径(mm)
+    flowParams.dOrificeD = 150.25; // 孔板孔径(mm)
+    flowParams.dPipeType = 0; // 管道类型
+    flowParams.dPtmode = 0; // 取压方式(0=法兰取压，1=角接取压)
+
+    // 设置材料参数
+    flowParams.dPipeMaterial = 2; // 20号钢
+    flowParams.dOrificeMaterial = 9; // 镍铬合金
+
+    // 设置天然气组分(示例: 95%甲烷，5%其他)
+    // 初始化天然气组分数组(GB/T 21446-2008 典型示例组成)
+    int i=0;
+    for ( i = 0; i < NUMBEROFCOMPONENTS; i++) {
+        flowParams.dNG_Compents[i] = 0.0; // 先全部初始化为0
+    }
+
+	flowParams.dNG_Compents[0] = 88.36; // 甲烷(CH4)
+    flowParams.dNG_Compents[1] = 0.68; // 氮气(N2)
+    flowParams.dNG_Compents[2] = 1.57; // 二氧化碳(CO2)
+    flowParams.dNG_Compents[3] =6.25; // 乙烷(C2H6)
+    flowParams.dNG_Compents[4] = 2.4; // 丙烷(C3H8)
+    flowParams.dNG_Compents[5] = 0.00; // 水(H2O)
+    flowParams.dNG_Compents[6] = 0.00; // 硫化氢(H2S)
+    flowParams.dNG_Compents[7] = 0.04; // 氢气(H2)
+    flowParams.dNG_Compents[8] = 0.00; // 一氧化碳(CO)
+    flowParams.dNG_Compents[9] = 0.00; // 氧气(O2)
+    flowParams.dNG_Compents[10] = 0.15; // 异丁烷(i-C4H10)
+    flowParams.dNG_Compents[11] = 0.35; // 正丁烷(n-C4H10)
+    flowParams.dNG_Compents[12] = 0.05; // 异戊烷(i-C5H12)
+    flowParams.dNG_Compents[13] = 0.1; // 正戊烷(n-C5H12)
+    flowParams.dNG_Compents[14] = 0.01; // 己烷(C6H14)
+    flowParams.dNG_Compents[15] = 0.0; // 庚烷(C7H16)
+    flowParams.dNG_Compents[16] = 0.0; // 辛烷(C8H18)
+    flowParams.dNG_Compents[17] = 0.0; // 壬烷(C9H20)
+    flowParams.dNG_Compents[18] = 0.0; // 癸烷(C10H22)
+    flowParams.dNG_Compents[19] = 0.04; // 氦气(He)
+    flowParams.dNG_Compents[20] = 0.0; // 其他组分
+    
+       // 调用流量计算函数
+    OFlowCal(&flowParams, &ngParams);
+    // 打印计算结果
+    printf("工况条件信息:\n");
+    printf("标准参比条件: %d\n", flowParams.dCbtj);
+    printf("计量参比压力: %.2f\n", flowParams.dPb_M);
+    printf("计量参比温度: %.2f\n", flowParams.dTb_M);
+    printf("能量参比压力: %.2f\n", flowParams.dPb_E);
+    printf("能量参比温度: %.2f\n", flowParams.dTb_E);
+    printf("大气压力: %.2f Pa\n", flowParams.dPatm);
+    printf("天然气组分:\n");
+    
+    for ( i = 0; i < 21; i++) {
+        printf("  组分 %d: %.6f\n", i, flowParams.dNG_Compents[i]);
+    }
+
+    printf("\n仪表参数:\n");
+    printf("仪表类型: %d\n", flowParams.dMeterType);
+    printf("核心类型: %d\n", flowParams.dCoreType);
+    printf("取压方式: %d\n", flowParams.dPtmode);
+    printf("管道类型: %d\n", flowParams.dPipeType);
+    printf("管道内径: %.2f mm\n", flowParams.dPipeD);
+    printf("管道材质: %d\n", flowParams.dPipeMaterial);
+    printf("孔板直径: %.2f mm\n", flowParams.dOrificeD);
+    printf("孔板材质: %d\n", flowParams.dOrificeMaterial);
+
+    printf("\n测量值:\n");
+    printf("压力: %.2f Pa\n", flowParams.dPf);
+    printf("压力类型: %d\n", flowParams.dPfType);
+    printf("温度: %.2f K\n", flowParams.dTf);
+    printf("差压: %.2f Pa\n", flowParams.dDp);
+    printf("仪表系数: %.6f\n", flowParams.dMeterFactor);
+    printf("脉冲数: %.2f\n", flowParams.dPulseNum);
+
+    printf("\n计算结果:\n");
+    printf("膨胀系数: %.6f\n", flowParams.dE);
+    printf("相对密度系数: %.6f\n", flowParams.dFG);
+    printf("超压缩系数: %.6f\n", flowParams.dFT);
+    printf("动力粘度: %.6f\n", flowParams.dDViscosity);
+    printf("热膨胀系数: %.6f\n", flowParams.dDExpCoefficient);
+    printf("管道雷诺数: %.2f\n", flowParams.dRnPipe);
+    printf("孔板弯曲系数: %.6f\n", flowParams.dBk);
+    printf("管道粗糙度: %.6f\n", flowParams.dRoughNessPipe);
+    printf("流出系数: %.6f\n", flowParams.dCd);
+    printf("流出系数修正: %.6f\n", flowParams.dCdCorrect);
+    printf("喷嘴流出系数: %.6f\n", flowParams.dCdNozell);
+    printf("标况体积流量: %.6f Nm3/s\n", flowParams.dVFlowb);
+    printf("工况体积流量: %.6f m3/s\n", flowParams.dVFlowf);
+    printf("质量流量: %.6f t/s\n", flowParams.dMFlowb);
+    printf("能量流量: %.6f MJ/s\n", flowParams.dEFlowb);
+    printf("流速: %.6f m/s\n", flowParams.dVelocityFlow);
+    printf("压力损失: %.6f\n", flowParams.dPressLost);
+    printf("直径比: %.6f\n", flowParams.dBeta);
+    printf("等熵指数: %.6f\n", flowParams.dKappa);
+    printf("压缩因子: %.6f\n", flowParams.dFpv);
+
+    printf("状态: %ld\n", ngParams.lStatus);
+    printf("强制更新标志: %d\n", ngParams.bForceUpdate);
+    printf("混合比:\n");
+    for ( i = 0; i < 21; i++) {
+        printf("  组分 %d: %.6f\n", i, ngParams.adMixture[i]);
+    }
+    printf("参比条件: %d\n", ngParams.dCbtj);
+    printf("标准压力: %.2f Pa\n", ngParams.dPb);
+    printf("标准温度: %.2f K\n", ngParams.dTb);
+    printf("工作压力: %.2f Pa\n", ngParams.dPf);
+    printf("工作温度: %.2f K\n", ngParams.dTf);
+
+    printf("\nAGA 8 详细计算结果:\n");
+    printf("平均分子量: %.6f\n", ngParams.dMrx);
+    printf("标准条件下压缩因子: %.6f\n", ngParams.dZb);
+    printf("工作条件下压缩因子: %.6f\n", ngParams.dZf);
+    printf("超压缩因子: %.6f\n", ngParams.dFpv);
+    printf("标准条件下摩尔密度: %.6f moles/dm3\n", ngParams.dDb);
+    printf("工作条件下摩尔密度: %.6f moles/dm3\n", ngParams.dDf);
+    printf("标准条件下密度: %.6f kg/m3\n", ngParams.dRhob);
+    printf("工作条件下密度: %.6f kg/m3\n", ngParams.dRhof);
+    printf("理想相对密度: %.6f\n", ngParams.dRD_Ideal);
+    printf("实际相对密度: %.6f\n", ngParams.dRD_Real);
+
+    printf("\n热力学性质:\n");
+    printf("理想焓: %.6f\n", ngParams.dHo);
+    printf("实际焓: %.6f J/kg\n", ngParams.dH);
+    printf("实际熵: %.6f J/kg-mol.K\n", ngParams.dS);
+    printf("理想定压比热: %.6f J/kg-mol.K\n", ngParams.dCpi);
+    printf("实际定压比热: %.6f J/kg-mol.K\n", ngParams.dCp);
+    printf("实际定容比热: %.6f J/kg-mol.K\n", ngParams.dCv);
+    printf("比热比: %.6f\n", ngParams.dk);
+    printf("等熵指数: %.6f\n", ngParams.dKappa);
+    printf("声速: %.6f m/s\n", ngParams.dSOS);
+    printf("临界流函数: %.6f\n", ngParams.dCstar);
+
+    printf("\n单位摩尔高热值: %.6f\n", ngParams.dHhvMol);
+    printf("单位摩尔低热值: %.6f\n", ngParams.dLhvMol);
+    return 1;
+	
+}
+
+BOOL WINAPI DllMain(HINSTANCE hinstDLL,DWORD fdwReason,LPVOID lpvReserved)
+{
+	switch(fdwReason)
+	{
+		case DLL_PROCESS_ATTACH:
+		{
+			break;
+		}
+		case DLL_PROCESS_DETACH:
+		{
+			break;
+		}
+		case DLL_THREAD_ATTACH:
+		{
+			break;
+		}
+		case DLL_THREAD_DETACH:
+		{
+			break;
+		}
+	}
+	
+	/* Return TRUE on success, FALSE on failure */
+	return TRUE;
+}

NGToolsLibrary/NGToolsLibrary.library.~u

@@ -0,0 +1,4 @@
+Administrator
+XTZJ-20241119ZU
+20200
+638977308197029107