Gekko结合CoolProp求解微分方程遇错及双变量插值问题求助

1. 微分方程定义错误 #

原因 #

Gekko无法直接对符号变量（如T作为m.Var）调用CoolProp函数。错误源于尝试在ODE方程中直接使用CoolProp，导致生成无效的符号表达式。与Scipy（在每个时间步数值计算）不同，Gekko需要符号关系或插值函数才能与求解器兼容。

解决方案 #

预计算温度范围内所需的热物性参数（密度d、比热c），然后使用Gekko的插值函数为ODE创建符号关系。

示例代码：

import numpy as np
import CoolProp.CoolProp as CP
from gekko import GEKKO
import matplotlib.pyplot as plt

# 初始化Gekko
m = GEKKO(remote=False)

# 时间数组
m.time = np.linspace(0, 100, 101)

# 温度变量
T = m.Var(value=78+273)

# 预计算温度范围内的物性参数
T_range = np.linspace(350, 375, 26)  # 覆盖预期温度范围
d_vals = [CP.PropsSI('D', 'T', t, 'Q', 1, "Water") for t in T_range]
c_vals = [CP.PropsSI('C', 'T', t, 'Q', 1, "Water") for t in T_range]  # 修正原代码中的'O'为比热参数'C'

# 创建密度和比热的样条插值
d = m.Var()
c = m.Var()
m.spline(T, d, T_range, d_vals)
m.spline(T, c, T_range, c_vals)

# 定义ODE方程
m.Equation(T.dt() == (50 + 0.1/d*T)/(1000*c))

# 设置求解模式并求解
m.options.IMODE = 4  # 动态仿真模式
m.solve(disp=True)

# 绘制结果
plt.plot(m.time, T.value)
plt.xlabel('时间 (s)')
plt.ylabel('温度 (K)')
plt.title('Gekko温度仿真')
plt.grid(True)
plt.show()

2. B样条配置错误（水压力插值） #

原因 #

Gekko的二维B样条需要指定每个维度的样条阶数（kx、ky）和平滑因子（sf）。当未提供这些参数时，Gekko会尝试从未配置的信息文件中读取，导致错误。

解决方案 #

在m.bspline调用中添加kx、ky和sf参数。通常使用kx=3、ky=3（三次样条）和sf=0（无平滑，完全拟合数据）。

修正后的代码片段：

# 数据生成代码保持不变
m.bspline(T1, D1, p1, T_data, D_data, p_data, kx=3, ky=3, sf=0)

3. 氮介质B样条插值失败 #

原因 #

同样是缺少样条配置参数（kx、ky、sf）。此外需确保输入数据维度正确（h_data为匹配D_data×u_data的二维数组）。

解决方案 #

在m.bspline调用中添加样条参数，修正后的完整代码：

import numpy as np
import CoolProp.CoolProp as CP
from gekko import GEKKO

m = GEKKO(remote=False)
u1 = m.Param(value=-81e3)
D1 = m.Param(value=808)
h1 = m.Var()

D_data = np.array([807.5, 808.0, 808.5, 809.0, 809.5])
u_data = np.array([-121000, -111000, -101000, -91000, -81000])
h_data = [] 
for n in D_data:
    row = []     
    for i in u_data:         
        row.append(CP.PropsSI('H','U',i,'D',n,"Nitrogen"))
    h_data.append(row) 

# 定义带显式参数的B样条
m.bspline(D1, u1, h1, D_data, u_data, h_data, kx=3, ky=3, sf=0)

m.options.SOLVER = 3
m.options.DIAGLEVEL = 10
m.solve(disp=True)

print(f"插值焓值: {h1.value[0]}")

Gekko与Scipy的核心差异要点 #

• 符号 vs 数值: Gekko使用符号方程进行优化/仿真，因此直接对Gekko变量调用数值函数（如CoolProp）无效，需通过插值创建符号关系。
• IMODE设置: 务必正确设置IMODE（如4用于动态仿真，3用于稳态求解）。
• 插值参数: 二维B样条需显式指定kx、ky和sf，避免配置错误。

内容的提问来源于stack exchange，提问作者Visitor