commit eb2717f0994b25d00a490ad294eb7258961a13d1
Author: ljz <425868052@qq.com>
Date: Fri May 8 15:15:12 2026 +0800
first commit
diff --git a/README.md b/README.md
new file mode 100644
index 0000000..21a6d42
--- /dev/null
+++ b/README.md
@@ -0,0 +1,38 @@
+### 3 分钟了解如何进入开发
+
+欢迎使用云效代码管理 Codeup，通过阅读以下内容，你可以快速熟悉 Codeup ，并立即开始今天的工作。
+
+### 提交**文件**
+
+Codeup 支持两种方式进行代码提交：网页端提交，以及本地 Git 客户端提交。
+
+* 如需体验本地命令行操作，请先安装 Git 工具，安装方法参见[安装Git](https://help.aliyun.com/document_detail/153800.html)。
+
+* 如需体验 SSH 方式克隆和提交代码，请先在平台账号内配置 SSH 公钥，配置方法参见[配置 SSH 密钥](https://help.aliyun.com/document_detail/153709.html)。
+
+* 如需体验 HTTP 方式克隆和提交代码，请先在平台账号内配置克隆账密，配置方法参见[配置 HTTPS 克隆账号密码](https://help.aliyun.com/document_detail/153710.html)。
+
+现在，你可以在 Codeup 中提交代码文件了，跟着文档「[__提交第一行代码__](https://help.aliyun.com/document_detail/153707.html?spm=a2c4g.153710.0.0.3c213774PFSMIV#6a5dbb1063ai5)」一起操作试试看吧。
+
+
+
+
+### 进行代码检测
+
+开发过程中，为了更好的维护你的代码质量，你可以开启 Codeup 内置开箱即用的「[代码检测服务](https://help.aliyun.com/document_detail/434321.html)」，开启后提交或合并请求的变更将自动触发检测，识别代码编写规范和安全漏洞问题，并及时提供结果报表和修复建议。
+
+
+
+### 开展代码评审
+
+功能开发完毕后，通常你需要发起「[代码评审并执行合并](https://help.aliyun.com/document_detail/153872.html)」，Codeup 支持多人协作的代码评审服务，你可以通过「[保护分支设置合并规则](https://help.aliyun.com/document_detail/153873.html?spm=a2c4g.203108.0.0.430765d1l9tTRR#p-4on-aep-l5q)」策略及「[__合并请求设置__](https://help.aliyun.com/document_detail/153874.html?spm=a2c4g.153871.0.0.3d38686cJpcdJI)」对合并过程进行流程化管控，同时提供在线代码评审及冲突解决能力，让评审过程更加流畅。
+
+
+
+### 成员协作
+
+是时候邀请成员一起编写卓越的代码工程了，请点击左下角「成员」邀请你的小伙伴开始协作吧！
+
+### 更多
+
+Git 使用教学、高级功能指引等更多说明，参见[Codeup帮助文档](https://help.aliyun.com/document_detail/153402.html)。
\ No newline at end of file
diff --git a/brayton_cycle_test.py b/brayton_cycle_test.py
new file mode 100644
index 0000000..9de7a8c
--- /dev/null
+++ b/brayton_cycle_test.py
@@ -0,0 +1,980 @@
+# -*- coding: utf-8 -*-
+"""
+Created on Mon Dec 22 15:35:07 2025
+
+"""
+import ctREFPROP.ctREFPROP as ct
+from scipy.optimize import minimize_scalar
+import numpy as np
+import matplotlib.pyplot as plt
+
+class CO2PropertyCalculator():
+ """二氧化碳物性计算"""
+ def __init__(self, refprop_path = None):
+ """初始化库"""
+ self.rp = ct.REFPROPFunctionLibrary(refprop_path)
+ # 设置流体文件
+ self.rp.SETUPdll(1, 'CO2.FLD', 'HMX.BNC', 'DEF')
+ # 设置单位
+ self.rp.SETUPdll(2, 'SI', 'SI', 'DEF')
+ self.z = [1.0]
+ self.mw = self.rp.WMOLdll(self.z)
+ def calculate_properties(self, T=None, P=None, h=None, s=None):
+ """计算二氧化碳物性"""
+ if T is not None and P is not None:
+ # 已知Tp
+ result = self.rp.TPFLSHdll(T, P, self.z)
+ properties = {
+ 'T': T,
+ 'P': P,
+ 'h': result.h,
+ 's': result.s,
+ }
+ elif P is not None and h is not None:
+ # 已知Ph
+ result = self.rp.PHFLSHdll(P, h, self.z)
+ properties = {
+ 'T': result.T,
+ 'P': P,
+ 'h': h,
+ 's': result.s,
+ }
+ elif T is not None and h is not None:
+ # 已知Th
+ result = self.rp.THFLSHdll(T, h, self.z)
+ properties = {
+ 'T': T,
+ 'P': result.P,
+ 'h': h,
+ 's': result.s,
+ }
+ elif P is not None and s is not None:
+ # 已知Ps
+ result = self.rp.PSFLSHdll(P, s, self.z)
+ properties = {
+ 'T': result.T,
+ 'P': P,
+ 'h': result.h,
+ 's': s,
+ }
+ else:
+ raise ValueError("提供的参数不足")
+ if result.ierr > 0:
+ raise ValueError(f"REFPROP计算错误:{result.ierr}")
+
+ # 补充提取的物性，这里由于后续还要使用，不进行参数变换
+ properties['D'] = result.D
+ properties['cp'] = result.Cp
+ properties['cv'] = result.Cv,
+ return properties
+
+class Compressor():
+ """压缩机类"""
+ def __init__(self, name, eff):
+ """
+ 初始化参数
+
+ name: 名称
+ eff: 等熵效率
+ Wc: 压缩功
+ inlet_state: 入口参数
+ outlet_state: 出口参数
+ outlet_state_is: 等熵状态下出口参数
+ """
+
+ self.name = name
+ self.eff = eff
+ self.variables = None
+
+ def calculator(self, p_in, T_in, p_out, property_calculator):
+ # 先计算熵值
+ inlet_state = property_calculator.calculate_properties(T=T_in, P=p_in)
+ mw = property_calculator.mw
+ s = inlet_state['s']
+ h_in = inlet_state['h']
+ outlet_state_is = property_calculator.calculate_properties(P=p_out, s=s)
+ h_out_is = outlet_state_is['h']
+ h_out = h_in + (h_out_is - h_in) / self.eff
+ outlet_state = property_calculator.calculate_properties(P=p_out, h=h_out)
+ Wc = h_out - h_in
+ T_out = outlet_state['T']
+ # 计算结果
+ self.variables = {
+ 'name': self.name,
+ 'inlet_state':{
+ 'P': p_in, # 压强(kPa)
+ 'T': T_in,
+ 'h': h_in/mw, # 比焓(J/mol)->(J/kg)
+ 's': s/mw, # 比熵(J/mol.K)->(J/kg.K)
+ },
+ 'outlet_state':{
+ 'P': p_out, # 压强(kPa)
+ 'T': T_out,
+ 'h': h_out/mw, # 比焓(kJ/mol)->(kJ/kg)
+ 's': s/mw, # 比熵(kJ/mol.K)->(kJ/kg.K)
+ },
+ 'eff': self.eff,
+ 'Wc': Wc/mw, # 压缩功(kJ/mol)->(kJ/kg)
+ 'pi': p_out / p_in
+ }
+
+class Turbine():
+ """透平类"""
+ def __init__(self, name, eff):
+ """
+ 初始化参数
+ name: 名称
+ eff: 透平效率
+ """
+ self.name = name
+ self.eff = eff
+ self.variables = None
+
+ def calculator(self, p_in, T_in, p_out, property_calculator):
+ """涡轮参数计算"""
+ inlet_state = property_calculator.calculate_properties(P=p_in, T=T_in)
+ mw = property_calculator.mw
+ s = inlet_state['s']
+ h_in = inlet_state['h']
+ outlet_state_is = property_calculator.calculate_properties(P=p_out, s=s)
+ h_out_is = outlet_state_is['h']
+ h_out = h_in - (h_in - h_out_is) * self.eff
+ outlet_state = property_calculator.calculate_properties(P=p_out, h=h_out)
+ T_out = outlet_state['T']
+ Wt = h_in - h_out
+ self.variables = {
+ 'name': self.name,
+ 'inlet_state':{
+ 'P': p_in, # 压强(kPa)
+ 'T': T_in,
+ 'h': h_in/mw, # 比焓(J/mol)->(J/kg)
+ 's': s/mw, # 比熵(J/mol.K)->(J/kg.K)
+ },
+ 'outlet_state':{
+ 'P': p_out, # 压强(kPa)
+ 'T': T_out,
+ 'h': h_out/mw, # 比焓(J/mol)->(J/kg)
+ 's': s/mw, # 比熵(J/mol.K)->(J/kg.K)
+ },
+ 'eff': self.eff,
+ 'Wt': Wt/mw, # 透平做功(J/mol)->(J/kg)
+ 'pi': p_in / p_out
+ }
+
+class Recuperator():
+ """换热器类"""
+ def __init__(self, name, eff, x=0):
+ """
+ 初始化参数
+ name: 名称
+ eff: 换热效率(基于焓的计算方法)
+ """
+ self.name = name
+ self.eff = eff
+ self.Q_ex = None
+ self.variables = None
+ self.x = x
+
+ def calculator(self, cold_inlet_state, hot_inlet_state, bypass_info, ploss, property_calculator):
+ """
+ 计算换热器两侧参数，默认逆流
+ bypass_info: 是否存在分流，0为不存在，1为存在
+
+ 下标含义
+ ----------
+ 1: 换热器冷端入口
+ 2: 换热器冷端出口
+ 3: 换热器热端入口
+ 4: 换热器热端出口
+ ass: 迭代中间变量，假设值
+ """
+
+ # 计算两入口参数, 这里单位是kg
+ p1 = cold_inlet_state['P']
+ T1 = cold_inlet_state['T']
+ h1 = cold_inlet_state['h']
+ p2 = p1 * (1 - ploss)
+ p3 = hot_inlet_state['P']
+ T3 = hot_inlet_state['T']
+ h3 = hot_inlet_state['h']
+ p4 = p3 * (1 - ploss)
+ mw = property_calculator.mw
+
+ # 设定质量流量
+ m_cold = (1-self.x) if bypass_info == 1 else 1.0
+ m_hot = 1.0
+
+# =============================================================================
+# # 采用焓差效能的方式来计算换热器进出口参数
+# # 假设最大温差发生在冷端, 计算冷端出口温度和比焓
+# cold_outlet_state = property_calculator.calculate_properties(T=T3, P=p2)
+# h2 = cold_outlet_state['h'] / mw
+# Q_ass_cold = m_cold * abs(h2 - h1)
+#
+# # 假设最大温差发生在热端, 计算热端出口温度和比焓
+# hot_outlet_state = property_calculator.calculate_properties(T=T1, P=p4)
+# h4 = hot_outlet_state['h'] / mw
+# Q_ass_hot = m_hot * abs(h3 - h4)
+#
+# # 比较两个可能的Q，取最小值与焓差效能的乘积作为实际换热量
+# self.Q_ex = min(Q_ass_hot, Q_ass_cold) * self.eff
+# # 由实际换热量计算出口焓和出口状态
+# h2 = h1 + self.Q_ex / m_cold
+# h4 = h3 - self.Q_ex / m_hot
+# cold_outlet_state = property_calculator.calculate_properties(P=p2, h=h2*mw)
+# hot_outlet_state = property_calculator.calculate_properties(P=p4, h=h4*mw)
+# T2 = cold_outlet_state['T']
+# T4 = hot_outlet_state['T']
+# =============================================================================
+
+ # 采用温差效能的方式来计算换热器进出口参数
+ T_ass_max = abs(T1 - T3)
+ # 假设最大温差发生在冷端, 计算冷端出口温度和比焓
+ T2 = T1 + T_ass_max * self.eff
+ cold_outlet_state = property_calculator.calculate_properties(T=T2, P=p2)
+ h2 = cold_outlet_state['h'] / mw
+ Q_ass_cold = m_cold * abs(h2 - h1)
+
+ # 假设最大温差发生在热端, 计算热端出口温度和比焓
+ T4 = T3 - T_ass_max * self.eff
+ hot_outlet_state = property_calculator.calculate_properties(T=T4, P=p4)
+ h4 = hot_outlet_state['h'] / mw
+ Q_ass_hot = m_hot * abs(h3 - h4)
+
+ # 比较两个可能的Q，取最小值与焓差效能的乘积作为实际换热量
+ self.Q_ex = min(Q_ass_hot, Q_ass_cold)
+ # 由实际换热量计算出口焓和出口状态
+ h2 = h1 + self.Q_ex / m_cold
+ h4 = h3 - self.Q_ex / m_hot
+ cold_outlet_state = property_calculator.calculate_properties(P=p2, h=h2*mw)
+ hot_outlet_state = property_calculator.calculate_properties(P=p4, h=h4*mw)
+ T2 = cold_outlet_state['T']
+ T4 = hot_outlet_state['T']
+
+ # 拼装变量
+ self.variables = {
+ 'name': self.name,
+ 'cold_inlet_state':{
+ 'P': p1,
+ 'T': T1,
+ 'h': h1,
+ 's': cold_inlet_state['s']
+ },
+ 'cold_outlet_state':{
+ 'P': p2,
+ 'T': cold_outlet_state['T'],
+ 'h': h2,
+ 's': cold_outlet_state['s']/mw
+ },
+ 'hot_inlet_state':{
+ 'P': p3,
+ 'T': T3,
+ 'h': h3,
+ 's': hot_inlet_state['s']
+ },
+ 'hot_outlet_state':{
+ 'P': p4,
+ 'T': T4,
+ 'h': h4,
+ 's': hot_outlet_state['s']/mw
+ },
+ 'eff': self.eff,
+ 'Q_exchange': self.Q_ex
+ }
+ def check_pinch_point(self, property_calculator, num_segments=20):
+ """
+ 换热器内部夹点校验
+ 将换热量均分为 num_segments 段，检查内部每个微元的冷热流体温度
+ """
+ h_cold_in = self.variables['cold_inlet_state']['h']
+ h_hot_in = self.variables['hot_inlet_state']['h']
+ p_cold = self.variables['cold_inlet_state']['P']
+ p_hot = self.variables['hot_inlet_state']['P']
+
+ m_cold = (1 - self.x) if self.name == "Low Temperature recuprerator" else 1.0
+ m_hot = 1.0
+
+ dQ = self.Q_ex / num_segments
+
+ # 沿冷流体流动方向步进检查
+ for i in range(num_segments + 1):
+ q_current = i * dQ
+
+ # 当前微元截面的焓值
+ h_cold_local = h_cold_in + q_current / m_cold
+ h_hot_local = (h_hot_in - self.Q_ex / m_hot) + q_current / m_hot
+
+ # 查温度
+ T_cold_local = property_calculator.calculate_properties(P=p_cold, h=h_cold_local * property_calculator.mw)['T']
+ T_hot_local = property_calculator.calculate_properties(P=p_hot, h=h_hot_local * property_calculator.mw)['T']
+
+ # 如果热流体温度低于等于冷流体温度 (设定一个 0.1K 的最小逼近温差容差)
+ if T_hot_local - T_cold_local < 0.1:
+ return False # 发生温度交叉，物理不可行！
+
+ return True # 校验通过
+
+class Heater():
+ """加热器类"""
+ def __init__(self, name):
+ self.name = name
+ self.variables = None
+ def calculator(self, inlet_state, outlet_state):
+ h_in = inlet_state['h']
+ h_out = outlet_state['h']
+ Q_input = h_out - h_in
+ self.variables = {
+ 'name': self.name,
+ 'inlet_state':{
+ 'P': inlet_state['P'], # 压强(kPa)
+ 'T': inlet_state['T'],
+ 'h': inlet_state['h'],
+ 's': inlet_state['s'],
+ },
+ 'outlet_state':{
+ 'P': outlet_state['P'], # 压强(kPa)
+ 'T': outlet_state['T'],
+ 'h': outlet_state['h'],
+ 's': outlet_state['s'],
+ },
+ 'Q_in': Q_input,
+ }
+
+class Condenser():
+ """冷凝器类"""
+ def __init__(self, name):
+ self.name = name
+ self.variables = None
+ def calculator(self, inlet_state, outlet_state):
+ h_in = inlet_state['h']
+ h_out = outlet_state['h']
+ Q_output = h_in - h_out
+ self.variables = {
+ 'name': self.name,
+ 'inlet_state':{
+ 'P': inlet_state['P'], # 压强(kPa)
+ 'T': inlet_state['T'],
+ 'h': inlet_state['h'],
+ 's': inlet_state['s'],
+ },
+ 'outlet_state':{
+ 'P': outlet_state['P'], # 压强(kPa)
+ 'T': outlet_state['T'],
+ 'h': outlet_state['h'],
+ 's': outlet_state['s'],
+ },
+ 'Q_out': Q_output,
+ }
+class Concentrator():
+ """汇流组件"""
+ def __init__(self, name):
+ self.name = name
+ self.variables = None
+ def calculator(self, inlet_state_bypass, inlet_state_mroad, x, property_calculator):
+ h_in_bypass = inlet_state_bypass['h']
+ h_in_mroad = inlet_state_mroad['h']
+ p_in = inlet_state_bypass['P']
+ mw = property_calculator.mw
+ h_out = x * h_in_bypass + (1-x) * h_in_mroad
+ outlet_state = property_calculator.calculate_properties(P=p_in, h=h_out*mw)
+ self.variables = {
+ 'name': self.name,
+ 'inlet_state_bypass': inlet_state_bypass,
+ 'inlet_state_mroad': inlet_state_mroad,
+ 'outlet_state':{
+ 'P': p_in,
+ 'T': outlet_state['T'],
+ 'h': h_out,
+ 's': outlet_state['s']/mw
+ }
+ }
+class BraytonCycle():
+ """循环计算"""
+ def __init__(self, name, refprop_path = "C:/Program Files (x86)/REFPROP 10.0+/REFPROP"):
+ self.name = name
+ self.property_calculator = None
+ self.compressor = None
+ self.turbine = None
+ self.recuperator = None
+ self.heater = None
+ self.condenser = None
+ self.concentrator = None
+ self.refprop_path = refprop_path
+ self.property_calculator = CO2PropertyCalculator(self.refprop_path)
+
+ def cycle_eff_calculator(self):
+ if type(self.compressor) == list:
+ Wc = (self.compressor[0].variables['Wc'] +
+ self.compressor[1].variables['Wc'])
+ else:
+ Wc = self.compressor.variables['Wc']
+ Wt = self.turbine.variables['Wt']
+ Q_input = self.heater.variables['Q_in']
+ cycle_eff = (Wt - Wc) / Q_input
+ return cycle_eff
+
+ def SC(self, T_low, T_high, p_low, p_high, param = None):
+ if param == None:
+ param = {
+ 'compressor_eff': 0.98,
+ 'turbine_eff': 0.95
+ }
+ # 定义循环组件
+ self.heater = Heater(name = "Main heater")
+ self.condenser = Condenser(name = "Main condenser")
+ self.compressor = Compressor(name = "Main compressor", eff = param['compressor_eff'])
+ self.turbine = Turbine(name = "Turbine", eff = param['turbine_eff'])
+
+ # 计算循环参数
+ self.compressor.calculator(p_low, T_low, p_high, self.property_calculator)
+ self.turbine.calculator(p_high, T_high, p_low, self.property_calculator)
+ self.heater.calculator(self.compressor.variables['outlet_state'], self.turbine.variables['inlet_state'])
+ self.condenser.calculator(self.turbine.variables['outlet_state'], self.compressor.variables['inlet_state'])
+
+ # 计算循环效率
+ cycle_eff = self.cycle_eff_calculator()
+ return cycle_eff
+
+ def SRC(self, T_low, T_high, p_low, p_high, param = None):
+ if ('recuperator_eff' not in param) or param == None:
+ param = {
+ 'compressor_eff': 0.98,
+ 'turbine_eff': 0.95,
+ 'recuperator_eff': 0.85
+ }
+
+ # 定义循环组件
+ self.heater = Heater(name = "Main heater")
+ self.condenser = Condenser(name = "Main condenser")
+ self.compressor = Compressor(name = "Main compressor", eff = param['compressor_eff'])
+ self.turbine = Turbine(name = "Turbine", eff = param['turbine_eff'])
+ self.recuperator = Recuperator(name = "Recuperator", eff = param['recuperator_eff'])
+
+ # 计算循环参数
+ self.compressor.calculator(p_low, T_low, p_high, self.property_calculator)
+ self.turbine.calculator(p_high, T_high, p_low, self.property_calculator)
+ self.recuperator.calculator(self.compressor.variables['outlet_state'],
+ self.turbine.variables['outlet_state'], 0,
+ self.property_calculator)
+ self.condenser.calculator(self.recuperator.variables['hot_outlet_state'],
+ self.compressor.variables['inlet_state'])
+ self.heater.calculator(self.recuperator.variables['cold_outlet_state'],
+ self.turbine.variables['inlet_state'])
+
+ # 计算循环效率
+ cycle_eff = self.cycle_eff_calculator()
+ return cycle_eff
+
+ def RC(self, T_low, T_high, p_low, p_high, ploss, param = None):
+ if ('recompressor_eff' not in param or
+ 'highT_recuperator_eff' not in param or
+ 'x' not in param or
+ param == None):
+ param = {
+ 'compressor_eff': 0.98,
+ 'turbine_eff': 0.95,
+ 'recuperator_eff': 0.85,
+ 'recompressor_eff': 0.98,
+ 'highT_recuperator_eff': 0.85,
+ 'x': 0.9
+ }
+ self.x = param['x']
+ # 定义循环组件
+ self.heater = Heater(name = "Main heater")
+ self.condenser = Condenser(name = "Main_condenser")
+ self.turbine = Turbine(name = "Main Turbine",
+ eff = param['turbine_eff'])
+ main_compressor = Compressor(name = "Main compressor",
+ eff = param['compressor_eff'])
+ recompressor = Compressor(name = "Recompressor",
+ eff = param['recompressor_eff'])
+ lrecuperator = Recuperator(name = "Low Temperature recuprerator",
+ eff = param['recuperator_eff'],
+ x = self.x)
+ hrecuperator = Recuperator(name = "High Temperature recuperator",
+ eff = param['highT_recuperator_eff'])
+ self.concentrator = Concentrator(name = "Concentrator")
+ self.recuperator = []
+ self.compressor = []
+
+
+ # 计算组件进出口参数
+ # 先算压缩机和涡轮
+ main_compressor.calculator(p_low, T_low, p_high, self.property_calculator)
+ self.turbine.calculator(p_high*(1-ploss)**3,
+ T_high, p_low*(1-ploss)**(-3),
+ self.property_calculator)
+
+ # 假定低温回热器出口参数并进行迭代
+ mw = self.property_calculator.mw
+ T_hr_inlet = ((self.turbine.variables['outlet_state']['T']
+ + main_compressor.variables['outlet_state']['T'])/2)
+ # T_lr_inlet = self.turbine.variables['outlet_state']['T'] * 1.01
+ max_iter = 300
+ relax_fac = 0.4
+ for i in range(max_iter):
+ hr_inlet_state_hot_mol = self.property_calculator.calculate_properties(P=p_high*(1-ploss), T=T_hr_inlet)
+ # 单位转换成mol
+ hr_inlet_state_hot = {
+ 'P': hr_inlet_state_hot_mol['P'],
+ 'T': hr_inlet_state_hot_mol['T'],
+ 'h': hr_inlet_state_hot_mol['h']/mw,
+ 's': hr_inlet_state_hot_mol['s']/mw
+ }
+
+ hrecuperator.calculator(hr_inlet_state_hot,
+ self.turbine.variables['outlet_state'], 0,
+ ploss, self.property_calculator)
+ lrecuperator.calculator(main_compressor.variables['outlet_state'],
+ hrecuperator.variables['hot_outlet_state'], 1,
+ ploss, self.property_calculator)
+ recompressor.calculator(p_in = p_low/(1-ploss), p_out = p_high*(1-ploss),
+ T_in = lrecuperator.variables['hot_outlet_state']['T'],
+ property_calculator = self.property_calculator)
+ self.concentrator.calculator(recompressor.variables['outlet_state'],
+ lrecuperator.variables['cold_outlet_state'],
+ self.x, self.property_calculator)
+
+ Tc_outlet = self.concentrator.variables['outlet_state']['T']
+ err = abs(Tc_outlet - T_hr_inlet)
+ if err <= 1e-5:
+ break
+ if i == max_iter - 1:
+ raise ValueError(f"迭代次数超过范围，当前误差{err:.4f}")
+
+ T_hr_inlet = relax_fac * Tc_outlet + (1 - relax_fac) * T_hr_inlet
+
+ hrecuperator.calculator(self.concentrator.variables['outlet_state'],
+ self.turbine.variables['outlet_state'], 0,
+ ploss, self.property_calculator)
+ lrecuperator.calculator(main_compressor.variables['outlet_state'],
+ hrecuperator.variables['hot_outlet_state'], 1,
+ ploss, self.property_calculator)
+
+ main_compressor.variables['Wc'] *= (1-self.x)
+ recompressor.variables['Wc'] *= self.x
+ self.recuperator.append(lrecuperator)
+ self.recuperator.append(hrecuperator)
+ self.compressor.append(main_compressor)
+ self.compressor.append(recompressor)
+ self.condenser.calculator(lrecuperator.variables['hot_outlet_state'],
+ main_compressor.variables['inlet_state'])
+ self.condenser.variables['Q_out'] *= (1-self.x)
+ self.heater.calculator(hrecuperator.variables['cold_outlet_state'],
+ self.turbine.variables['inlet_state'])
+ # 计算循环效率
+ cycle_eff = self.cycle_eff_calculator()
+ return cycle_eff
+
+ def base_params_single_optimize(self, fixed_var, base_params, target_var_name, bounds):
+ """
+ 组件性能单变量优化器
+ : 固定边界条件
+ : 默认参数
+ : 优化变量名称
+ : 变量范围
+ """
+
+ if target_var_name not in base_params:
+ raise ValueError(f"参数{target_var_name}不在参数字典中")
+
+ def opt_fun(opt_var):
+ # 复制变量字典
+ opt_param = base_params.copy()
+ # 修改要优化的变量为参数
+ opt_param[target_var_name] = opt_var
+ # 带入循环参数计算
+ try:
+ eff = self.RC(
+ T_low = fixed_var['T_low'],
+ T_high= fixed_var['T_high'],
+ p_low = fixed_var['p_low'],
+ p_high = fixed_var['p_high'],
+ ploss = fixed_var['ploss'],
+ param = opt_param)
+ return -eff
+ except Exception as e:
+ return 0.0 # 遇到物性计算崩溃时返回极差值
+ res = minimize_scalar(opt_fun, bounds=bounds, method='bounded')
+ if res.success:
+ print(f"✅ 优化完成！")
+ print(f"👉 最佳 {target_var_name} = {res.x:.4f}")
+ print(f"👉 此时系统最高效率 = {-res.fun:.2%}\n")
+ else:
+ print("❌ 优化失败。")
+
+ return res
+
+ def fixed_params_single_optimize(self, fixed_params, target_var_name, params, bounds):
+ """
+ 边界条件单变量优化器
+ : 固定边界条件
+ : 默认参数
+ : 优化变量名称
+ : 变量范围
+ """
+ if target_var_name not in fixed_params:
+ raise ValueError(f"参数{target_var_name}不在参数字典中")
+
+ def opt_fun(opt_var):
+ # 复制变量
+ opt_params = fixed_params.copy()
+ # 变量替换
+ opt_params[target_var_name] = opt_var
+ # 带入循环
+ try:
+ eff = self.RC(T_low = opt_params['T_low'],
+ T_high = opt_params['T_high'],
+ p_low = opt_params['p_low'],
+ p_high = opt_params['p_high'],
+ ploss = opt_params['ploss'],
+ param = params)
+ return -eff
+ except Exception as e:
+ return 0.0
+ res = minimize_scalar(opt_fun, bounds=bounds, method='bounded')
+ if res.success:
+ print(f"✅ 优化完成！")
+ print(f"👉 最佳 {target_var_name} = {res.x:.4f}")
+ print(f"👉 此时系统最高效率 = {-res.fun:.2%}\n")
+ else:
+ print("❌ 优化失败。")
+ return res
+
+ def plot_optimization_landscape(self, fixed_params, params, target_var_name, bounds, res, num_points=50):
+ """
+ 绘制单变量优化地形图
+ :param target_var_name: 要扫描和优化的变量名（如 'x'）
+ :param bounds: 扫描和优化的范围 (min, max)
+ :param num_points: 扫描的采样点数量，越大曲线越平滑，但计算越慢
+ """
+ print(f"开始对【{target_var_name}】进行区间扫描，共计算 {num_points} 个点...")
+ plt.rcParams['font.sans-serif'] = ['SimHei'] # Windows 用黑体
+ plt.rcParams['axes.unicode_minus'] = False # 正常显示负号
+ # 1. 生成扫描数组
+ x_vals = np.linspace(bounds[0], bounds[1], num_points)
+ eff_vals = []
+ valid_x = [] # 记录那些没有报错的 x
+ if target_var_name not in fixed_params:
+ # 2. 遍历计算曲线上的点
+ for val in x_vals:
+ current_param = params.copy()
+ current_param[target_var_name] = val
+ try:
+ # 调用你的黑盒物理模型（注意：这里取正效率用于画图）
+ eff = self.RC(
+ T_low=fixed_params['T_low'],
+ T_high=fixed_params['T_high'],
+ p_low=fixed_params['p_low'],
+ p_high=fixed_params['p_high'],
+ ploss=fixed_params['ploss'],
+ param=current_param
+ )
+ # 如果系统加了夹点校验且没通过，可能会返回 None 或者抛异常
+ # 这里确保只有成功的点才画上去
+ eff_vals.append(eff * 100) # 乘以 100 转换为百分比
+ valid_x.append(val)
+ except Exception as e:
+ # 如果某个 x 导致计算崩溃，我们跳过这个点，不画它
+ pass
+ else:
+ for val in x_vals:
+ current_param = fixed_params.copy()
+ current_param[target_var_name] = val
+ try:
+ # 调用你的黑盒物理模型（注意：这里取正效率用于画图）
+ eff = self.RC(
+ T_low=current_param['T_low'],
+ T_high=current_param['T_high'],
+ p_low=current_param['p_low'],
+ p_high=current_param['p_high'],
+ ploss=current_param['ploss'],
+ param=params
+ )
+ # 如果系统加了夹点校验且没通过，可能会返回 None 或者抛异常
+ # 这里确保只有成功的点才画上去
+ eff_vals.append(eff * 100) # 乘以 100 转换为百分比
+ valid_x.append(val)
+ except Exception as e:
+ # 如果某个 x 导致计算崩溃，我们跳过这个点，不画它
+ pass
+ print("扫描完成！正在使用优化器寻找精确最高点...")
+
+ # 4. 开始绘图
+ plt.figure(figsize=(8, 6), dpi=120) # 设置画布大小和清晰度
+
+ # 画出目标函数曲线
+ plt.plot(valid_x, eff_vals, linestyle='-', color='#1f77b4', linewidth=2, label='系统热效率曲线')
+
+ # 如果优化成功，用醒目的红星标出最优点
+ if res.success:
+ best_x = res.x
+ best_eff = -res.fun * 100
+ plt.scatter(best_x, best_eff, color='red', marker='*', s=200, zorder=5, label=f'最优点 ({best_x:.4f}, {best_eff:.2f}%)')
+
+ # 画辅助虚线对齐坐标轴
+ plt.axvline(x=best_x, color='gray', linestyle='--', alpha=0.6)
+ plt.axhline(y=best_eff, color='gray', linestyle='--', alpha=0.6)
+ else:
+ print("优化结果未输入!")
+ # 设置图表装饰
+ plt.title(f'系统热效率随 {target_var_name} 的变化趋势', fontsize=14)
+ plt.xlabel(f'优化变量: {target_var_name}', fontsize=12)
+ plt.ylabel('循环热效率 η (%)', fontsize=12)
+ plt.grid(True, linestyle=':', alpha=0.7)
+ plt.legend(fontsize=11)
+
+ # 显示图像
+ plt.tight_layout()
+ plt.show()
+
+ # 调用测试
+ # plot_optimization_landscape('x', bounds=(0.6, 0.95), num_points=40)
+
+
+ def sweep_and_optimize(self, fixed_params, params, sweep_var, sweep_bounds, opt_var, opt_bounds, num_points=50):
+ """
+ 带内部动态优化的单变量扫描器 (极度通用版)
+
+ :param fixed_params: 固定的边界条件字典
+ :param params: 组件性能参数字典
+ :param sweep_var: 你要扫描/遍历的变量名 (例如 'T_low')
+ :param sweep_bounds: 扫描变量的范围 (min, max)
+ :param opt_var: 在每个扫描点下，你需要动态寻找最优值的变量名 (例如 'x')
+ :param opt_bounds: 优化变量的搜索范围 (min, max)
+ :param num_points: 扫描点数
+ """
+ print(f"\n🚀 开始执行嵌套扫描：")
+ print(f" - 扫描变量 (X轴): 【{sweep_var}】 范围 {sweep_bounds}")
+ print(f" - 内部动态优化变量: 【{opt_var}】 范围 {opt_bounds}")
+
+ # 1. 生成扫描节点
+ sweep_vals = np.linspace(sweep_bounds[0], sweep_bounds[1], num_points)
+
+ # 记录数据的列表
+ valid_sweep_vals = []
+ best_effs = []
+ best_opt_vals = []
+
+ # 2. 开始逐点扫描
+ for s_val in sweep_vals:
+
+ # 【核心1：每次必须使用干净的字典副本】
+ current_fixed = fixed_params.copy()
+ current_param = params.copy()
+
+ # 判断扫描变量是属于 fixed_params 还是 params，并赋值
+ if sweep_var in current_fixed:
+ current_fixed[sweep_var] = s_val
+ elif sweep_var in current_param:
+ current_param[sweep_var] = s_val
+ else:
+ raise ValueError(f"找不到扫描变量: {sweep_var}")
+
+ # 3. 定义内部优化目标函数 (闭包)
+ def inner_objective(guess_val):
+ # 将优化器猜的值赋给 opt_var
+ current_param[opt_var] = guess_val
+
+ try:
+ # 调用黑盒计算
+ eff = self.RC(
+ T_low=current_fixed['T_low'],
+ T_high=current_fixed['T_high'],
+ p_low=current_fixed['p_low'],
+ p_high=current_fixed['p_high'],
+ ploss=current_fixed['ploss'],
+ param=current_param
+ )
+ # 【预留口：此处可加入换热器内部夹点校验】
+ # 取出两个换热器进行夹点校验
+ for rec in self.recuperator:
+ if not rec.check_pinch_point(self.property_calculator):
+ return 0.0 # 核心！如果交叉了，直接返回 0 效率，强迫优化器换参数
+ return -eff
+ except Exception:
+ return 0.0
+
+ # 4. 调用一维优化器
+ res = minimize_scalar(inner_objective, bounds=opt_bounds, method='bounded')
+
+ # 5. 结果校验与存储
+ if res.success and -res.fun > 0:
+ best_eff = -res.fun * 100
+ best_opt_val = res.x
+
+ valid_sweep_vals.append(s_val)
+ best_effs.append(best_eff)
+ best_opt_vals.append(best_opt_val)
+
+ print(f"✔️ {sweep_var} = {s_val:.2f} | 寻得最优 {opt_var} = {best_opt_val:.4f} | 最高效率 = {best_eff:.2f}%")
+ else:
+ print(f"❌ {sweep_var} = {s_val:.2f} | 优化失败或物理无解，已跳过")
+
+ # 6. 调用画图方法 (将画图剥离，保持代码干净)
+ self._plot_results(sweep_var, valid_sweep_vals, best_effs, opt_var, best_opt_vals)
+
+ return valid_sweep_vals, best_effs, best_opt_vals
+
+ def _plot_results(self, sweep_var, x_data, y_eff_data, opt_var, y_opt_data):
+ """专门用来画图的内部方法，支持双Y轴"""
+ if not x_data:
+ print("没有有效数据可供绘制！")
+ return
+
+ plt.rcParams['font.sans-serif'] = ['SimHei']
+ plt.rcParams['axes.unicode_minus'] = False
+
+ fig, ax1 = plt.subplots(figsize=(9, 6), dpi=120)
+
+ # 画左Y轴：最高效率曲线
+ color1 = '#1f77b4'
+ ax1.set_xlabel(f'扫描变量: {sweep_var}', fontsize=12)
+ ax1.set_ylabel('最优循环热效率 η (%)', color=color1, fontsize=12)
+ ax1.plot(x_data, y_eff_data, color=color1, linewidth=2.5, label='系统热效率')
+ ax1.tick_params(axis='y', labelcolor=color1)
+ ax1.grid(True, linestyle=':', alpha=0.6)
+
+ # 画右Y轴：对应的最优分流量走势
+ ax2 = ax1.twinx()
+ color2 = '#d62728'
+ ax2.set_ylabel(f'匹配的最优动态变量: {opt_var}', color=color2, fontsize=12)
+ ax2.plot(x_data, y_opt_data, color=color2, linestyle='--', linewidth=2, label=f'最优 {opt_var} 值')
+ ax2.tick_params(axis='y', labelcolor=color2)
+
+ plt.title(f'系统最高效率及对应的最优 {opt_var} 随 {sweep_var} 的变化', fontsize=14)
+ fig.tight_layout()
+ plt.show()
+
+if __name__ == "__main__":
+ brayton1 = BraytonCycle(name = "simple brayton cycle test")
+ T_high = 650+273.15
+ T_low = 42+273.15
+ p_high = 20.0e6
+ p_low = 9.09e6
+ fixed_var = {
+ 'T_high': 650+273.15,
+ 'T_low': 42+273.15,
+ 'p_high': 20.0e3,
+ 'p_low': 8.16e3,
+ 'ploss': 0.01
+ }
+ # param = {
+ # 'compressor_eff': 0.9,
+ # 'turbine_eff': 0.93,
+ # 'recuprator_eff': 0.96
+ # }
+
+ # =========================================================================
+ # 焓差效能参数定义
+ # param = {
+ # 'x': 0.279,
+ # 'compressor_eff': 0.9,
+ # 'recompressor_eff': 0.855,
+ # 'turbine_eff': 0.93,
+ # 'recuperator_eff': 0.9354,
+ # 'highT_recuperator_eff': 0.8857
+ # }
+ # =========================================================================
+ # eff = brayton1.simple_brayton_cycle(T_low, T_high, p_low/1e3, p_high/1e3, param)
+ #
+ param = {
+ 'x': 0.279,
+ 'compressor_eff': 0.9,
+ 'recompressor_eff': 0.9,
+ 'turbine_eff': 0.93,
+ 'recuperator_eff': 0.94,
+ 'highT_recuperator_eff': 0.96
+ }
+ ploss = 0.01
+ eff = brayton1.RC(T_low, T_high, p_low/1e3, p_high/1e3, ploss, param)
+ c = brayton1.compressor
+ h = brayton1.heater
+ cond = brayton1.condenser
+ conc = brayton1.concentrator.variables
+ t = brayton1.turbine
+ r = brayton1.recuperator
+ vl = brayton1.recuperator[0].variables
+ vh = brayton1.recuperator[1].variables
+ Wc = brayton1.compressor[0].variables['Wc'] + brayton1.compressor[1].variables['Wc']
+ # Wc = brayton1.compressor.variables['Wc']
+ Wt = brayton1.turbine.variables['Wt']
+ Q_input = brayton1.heater.variables['Q_in']
+ Q_output = brayton1.condenser.variables['Q_out']
+ print(Q_input+Wc-(Q_output+Wt))
+ print((Wt - Wc)/Q_input)
+
+ # res_x = brayton1.base_params_single_optimize(fixed_var = fixed_var,
+ # base_params = param,
+ # target_var_name = 'x',
+ # bounds = (0.1,0.5))
+ # res_min_T = brayton1.fixed_params_single_optimize(fixed_params = fixed_var,
+ # target_var_name = 'T_low',
+ # params = param,
+ # bounds = (30+273.15, 60+273.15))
+
+ x_vals, eff_vals, opt_vals = brayton1.sweep_and_optimize(
+ fixed_params=fixed_var,
+ params=param,
+ sweep_var='T_low', # 扫描变量
+ sweep_bounds=(30+273.15, 60+273.15),
+ opt_var='x', # 动态优化变量
+ opt_bounds=(0.1, 0.5),
+ num_points=30)
+
+ # refprop_path = "C:/Program Files (x86)/REFPROP 10.0+/REFPROP"
+ # recuprator = Recuperator(name='Main Recuprator', eff=param['recuprator_eff'])
+ # calculator = CO2PropertyCalculator(refprop_path)
+ # cold_inlet_state = brayton1.compressor.variables['outlet_state']
+ # hot_inlet_state = brayton1.turbine.variables['outlet_state']
+ # recuprator.calculator(hot_inlet_state, cold_inlet_state, calculator)
+ # print(recuprator.variables)
+
+
+
+
+ # self.cycle_eff = None
+ # self.compressors = []
+ # self.turbines = []
+ # self.recuperators = []
+ # self.heaters = []
+ # self.condensers = []
+
+ # def add_compressor(self, name, eff):
+ # """添加压缩机"""
+ # compressor = Compressor(name, eff)
+ # self.compressors.append(compressor)
+ # return compressor
+
+ # def add_turbine(self, name, eff):
+ # """添加透平"""
+ # turbine = Turbine(name, eff)
+ # self.turbines.append(turbine)
+ # return turbine
+
+ # def add_heater(self, name):
+ # heater = Heater(name)
+ # self.heaters.append(heater)
+ # return heater
+
+ # def add_condenser(self,name):
+ # condenser = Condenser(name)
+ # self.condensers.append(condenser)
+ # return condenser
+
+ # def efficiency_calculation(self):
+ # cycle_Wc = 0.0
+ # for compressor in self.compressors:
+ # cycle_Wc += compressor.Wc
+ # cycle_Wt = 0.0
+ # for turbine in self.turbines:
+ # cycle_Wt += turbine.Wt
+ # cycle_Q_in = 0.0
+ # for heater in self.heaters:
+ # cycle_Q_in += self.heater.Q_input
+ # cycle_Q_out = 0.0
+ # for condenser in self.condensers:
+ # cycle_Q_out += self.condenser.Q_output
+ # cycle_eff = (cycle_Wt - cycle_Wc - cycle_Q_out) / cycle_Q_in
+ # return cycle_eff
+
\ No newline at end of file
diff --git a/ctREFTest.py b/ctREFTest.py
new file mode 100644
index 0000000..958986e
--- /dev/null
+++ b/ctREFTest.py
@@ -0,0 +1,33 @@
+# -*- coding: utf-8 -*-
+"""
+Created on Tue Dec 16 15:46:51 2025
+
+@author: LJZ
+"""
+
+import ctREFPROP.ctREFPROP as ct
+
+path = "C:/Program Files (x86)/REFPROP 10.0+/REFPROP"
+R = ct.REFPROPFunctionLibrary(path)
+R.SETPATHdll(path)
+fluid = "CO2.FLD"
+hfm = "HMX.BNC"
+hrf = "DEF"
+z = [1.0]
+ierr, herr = R.SETUPdll(1, fluid, hfm, hrf)
+R.SETUPdll(2, 'SI', 'SI', 'DEF')
+R.SETREFdll("DEF",1,[1.0],0,0,0,0)
+mw = R.WMOLdll(z)
+
+p = 20000 # kpa
+T = 391.05100419252705 # K
+
+
+result = R.TPFLSHdll(T,p,z)
+s_SI = result.s / mw
+s = result.s
+result2 = R.PSFLSHdll(p, s, z).T
+rho = result.D * mw
+h = result.h / mw
+
+x = min(1, 50)
\ No newline at end of file
diff --git a/refdata.xlsx b/refdata.xlsx
new file mode 100644
index 0000000..afe0d28
Binary files /dev/null and b/refdata.xlsx differ