Modelica/Dymola中加热管道的建模

Question

我目前正在学习化学工程，在我的学士学位论文中，我应该通过热口将两根管道连接在一起，模拟可用于过热器的加热管。尽管我付出了很大的努力来理解如何在Modelica中正确编码，但我的代码仍然无法正常工作，我变得非常绝望。

因此，该模型基本上必须适用于流体水和过热蒸汽，因此在非静止条件下只有单相流。传热应该是对流发生的。此外，我忽略了该模型中由于摩擦引起的压力损失。

以下是我对模型应该如何工作的想法： 我几乎在尝试构建一个像MSL中的模型"Dynamic Pipe"，只是更容易，以便从事相同主题的学生能够快速理解我的代码。所以我把管道分成了若干个节点n，第一个体积是入口状态，所以基本上这个状态并不真正属于管道。之后，平衡方程适用。我不太确定动量方程，因此非常感谢对它们的任何帮助。对流传热由MSL的"对流"模型Thermal.HeatTransfer.Components定义。 当使用流动源、具有固定压力的边界和壁上的固定温度测试模型时，我也收到错误"无法降低 DAE 指数"，我完全不知道这意味着什么。

另外，这是我的代码：

model Pipe_base3
//Import
import Modelica.SIunits.*;
import Modelica.Constants.pi;
replaceable package Medium =
Modelica.Media.Interfaces.PartialTwoPhaseMedium                          annotation (choicesAllMatching = true);
parameter Integer n=2;
parameter Integer np=1;
// Geometry==================================================================//
parameter Diameter d_pipe = 0.05 "Inner diameter of pipe"
annotation (Dialog(tab="Geometry"));
parameter Length L = 1 "Length of unit"
annotation (Dialog(tab="Geometry"));
parameter Area A_hex = pi * d_pipe * L
"Shell surface of pipe for heat exchange"                                                 annotation (Dialog(tab="Geometry"));
parameter Area A_q = (pi/4)*d_pipe^2
annotation (Dialog(tab="Geometry"));
//Initialisation=============================================================//
parameter Medium.Temperature T_start = 403.15 annotation (Dialog(tab="Initialization"));
parameter Medium.SpecificEnthalpy h_start = Medium.specificEnthalpy_pT(p_start, T_start) annotation (Dialog(tab="Initialization"));
parameter AbsolutePressure p_start = Medium.saturationPressure(T_start) annotation (Dialog(tab="Initialization"));
parameter Medium.MassFlowRate m_flow_start = 0.5 annotation (Dialog(tab="Initialization"));
//Temperature, pressure, energy==============================================//
Medium.Temperature T[n+1]( each start=T_start, fixed=false);
Medium.SpecificEnthalpy h[n+1]( each start=h_start, fixed=false);
Medium.AbsolutePressure p[n+1](each start=p_start, fixed=false);
HeatFlowRate Q_flow[n](fixed = false);
Energy U[n](min=0);
Energy KE[n]; //Kinetic Energy
Medium.ThermodynamicState state[n+1];
// Nondimensional Variables + HeatTransfer===================================//
Medium.PrandtlNumber Pr[n](fixed=false);
ReynoldsNumber Re[n](fixed=false);
Real Xi[n];
NusseltNumber Nu[n];
CoefficientOfHeatTransfer alpha[n];
// Thermodynamic properties==================================================//
Medium.SpecificInternalEnergy u[n](fixed=false);
Medium.DynamicViscosity eta[n];
Density rho[n+1];
Medium.SpecificHeatCapacity cp[n];
Medium.ThermalConductivity lambda_fluid[n];
//Segmental properties
Mass ms[n]; //Mass per Segment
MassFlowRate m_flow[n+1]( each start=m_flow_start/np, fixed=false);
Velocity w[n+1](fixed=false);
// Momentum
Force F_p[n];
Momentum I[n];
Force Ib_flow[n];
parameter Boolean init = false;
Modelica.Fluid.Interfaces.FluidPort_a fluidin( redeclare package Medium = Medium, m_flow(start = m_flow_start, min = 0), p(start = p_start))
annotation (Placement(transformation(extent={{-90,-100},{-70,-80}}),
iconTransformation(extent={{-90,-100},{-70,-80}})));
Modelica.Fluid.Interfaces.FluidPort_b fluidout( redeclare package Medium = Medium, m_flow(start = -m_flow_start, max = 0), p(start = p_start), h_outflow(start=h_start))
annotation (Placement(transformation(extent={{70,-100},{90,-80}}),
iconTransformation(extent={{70,-100},{90,-80}})));
Modelica.Thermal.HeatTransfer.Interfaces.HeatPort_a[n] heatport
annotation (Placement(transformation(extent={{-10,60},{10,80}}),
iconTransformation(extent={{-10,60},{10,80}})));
Modelica.Blocks.Interfaces.RealOutput[n] alpha_output annotation (Placement(
transformation(extent={{-100,38},{-140,78}}), iconTransformation(extent={{-100,
38},{-140,78}})));
protected 
parameter Volume vn = (A_q * L) / n; //Volume per segment
parameter Real x[n] = linspace((L/n), L, n);
parameter Length length = L/n;
initial equation 
for i in 1:(n+1) loop
//h[i] = Medium.specificEnthalpy_pTX(p_start, T_start, {1});
p[i] = p_start;
end for;
equation 
//Port equations=============================================================//
fluidout.p = p[n];
//fluidin.p-fluidout.p=p[1]-p[n+1];
fluidout.h_outflow = h[n];
fluidout.m_flow = -m_flow[n+1];
//===========================================================================//
h[1]=inStream(fluidin.h_outflow);
p[1]=fluidin.p;
state[1]=Medium.setState_ph(p[1],h[1]);
T[1]=Medium.temperature(state[1]);
rho[1]=Medium.density(state[1]);
m_flow[1]=fluidin.m_flow/np;
m_flow[1]=A_q*rho[1]*w[1];
for i in 1:(n) loop
// Heatport equations======================================================//
T[i] = heatport[i].T;
Q_flow[i] = heatport[i].Q_flow;
// Momentum Balance =======================================================//
der(I[i]) = Ib_flow[i] - F_p[i];
I[i]=m_flow[i]*length;
Ib_flow[i] = (p[i+1]*w[i+1]*w[i+1] - p[i]*w[i]*w[i])*A_q*np;
F_p[i] = (A_q*p[i+1]-A_q*p[i]);
// Energy Balance=========================================================//
U[i] = ms[i] * u[i];
KE[i] = 0.5*ms[i]*w[i+1]*w[i+1];
der(U[i]+KE[i])=m_flow[i]*(h[i]+0.5*w[i]) - m_flow[i+1]*(h[i+1]+0.5*w[i+1]) + Q_flow[i];
der(rho[i+1])= -((rho[i+1]-rho[i])*w[i+1] + (w[i+1]-w[i])*rho[i+1]); //Konti

ms[i]=vn*rho[i+1];
T[i+1]=Medium.temperature(state[i+1]);
state[i+1] = Medium.setState_ph(p[i+1], h[i+1], 1); //Sets thermodynamic state from which other properties can be determined
u[i] = Medium.specificInternalEnergy(state[i+1]);
cp[i] = Medium.specificHeatCapacityCp(state[i+1]);
rho[i+1] = Medium.density(state[i+1]);
eta[i] = Medium.dynamicViscosity(state[i+1]);
lambda_fluid[i] = Medium.thermalConductivity(state[i+1]);

Re[i] * eta[i] = (rho[i+1] * abs(w[i+1]) * d_pipe);
Pr[i] *lambda_fluid[i] = (eta[i] * cp[i]);
Xi[i] = (1.8 * log10(abs(Re[i])+1) - 1.5)^(-2);
Nu[i] = ((Xi[i]/8)*Re[i]*Pr[i])/(1+12.7*sqrt(Xi[i]/8)*((Pr[i])^(2/3)-1))*(1+(1/3)*(d_pipe/x[i])^(2/3));
Nu[i] = Modelica.Fluid.Pipes.BaseClasses.CharacteristicNumbers.NusseltNumber(alpha[i], d_pipe, lambda_fluid[i]);
alpha_output[i] = alpha[i] * (A_hex/n);
m_flow[i+1] = A_q * w[i+1] * rho[i+1];
// der(p[i]) = - w[i]*der(w[i]) * rho[i];
// 0 = m_flow[i-1] - m_flow[i];
// der(rho[i]) = -((rho[i]-rho[i-1])*w[i] + (v[i]-v[i-1])*rho[i]);
//m_flow[i] = A_q * w[i] * rho[i]; //Calculation of flow velocity
//ms[i] = vn * rho[i]; //Mass per segment
//Calculation of thermodynamic properties for each segment=================//

//Heat Transfer============================================================//

end for;
fluidin.h_outflow = h[1]; //
annotation (Icon(coordinateSystem(preserveAspectRatio=false, extent={{-100,-100},
{100,100}}), graphics={Line(
points={{-80,-80},{-80,94},{-80,100},{0,20},{80,100},{80,-80}},
color={0,0,255},
smooth=Smooth.None), Line(
points={{-60,-60},{-60,-48},{-60,0},{60,0},{60,-60},{48,-40},{72,-40},
{60,-60}},
color={0,0,255},
smooth=Smooth.None)}), __Dymola_selections);
end Pipe_base3;

提前非常感谢！

Answer 1

当我开始使用Modelica时，我遇到了同样的情况：我想要Modelica.Fluid.Pipes.DynamicPipe的功能，但复杂性较低(我希望代码更具可读性，层次结构更少)。所以，和你一样，我开始从头开始构建自己的管道模型。但是，由于我希望能够替换压降和传热相关性并具有很大的灵活性，因此我最终得到了一个与Modelica.Fluid.Pipes.DynamicPipe几乎相同复杂性的模型。

我给你的建议是

1. 构建您自己的简单动态管道模型，无需任何复杂的 特征。这只能用于教育目的(例如 让其他学生了解您的编码原理)
2. 了解如何使用Modelica.Fluid.Pipes.DynamicPipe来解决需要不同模型复杂性(段数、可更换压降和传热方法等)的问题。Modelica.Fluid.Examples.HeatExchanger是一个示例，说明如何使用Modelica.Fluid.Pipes.DynamicPipe对像您请求的热交换器一样进行建模。

在这里，我分享了一个非常简单的动态管道的例子，可以用作热交换器。管道由n管段组成，并利用了您可以实例化组件数组并在for循环中连接元素的事实。

至于动量平衡，正确/完整的方法是通过对作用在每个控制体积上的所有力求和来解释动量的变化(牛顿第二定律)。然而，在大多数集总模型中，稳态动量平衡就足够了，这将方程简化为质量流量和压降之间的线性或二次关系。Modelica.Fluid.Pipes.DynamicPipe有许多不同的压力/流量相关性可供选择。

此致敬意

雷内·贾斯特尼尔森

Answer 2

我已经构建了一个使用您的模型的小示例/测试。它应该是模型的一个非常简单的应用程序。不幸的是，我收到相同的错误消息：

找不到微分功能： Modelica.Media.Water.IF97_Utilities.waterBaseProp_ph(boundary1.p， pipe_base3_1.h[2]， 0， 1) 关于时间

指数减少基本上意味着模型包含没有未知数的方程。这是通过这些方程相对于时间的微分来解决的(可以发生多次)。有关更多信息，您可以查看
https://www.inf.ethz.ch/personal/cellier/Lect/NSDS/Ppt/nsds_ppt_engl.html 尤其是第16讲，可能还有之前的那些:)

因此，Modelica工具必须知道如何进行这种区分。对于方程式，这通常是自动完成的，但对于函数，它必须由开发人员指定。似乎没有为

Modelica.Media.Water.IF97_Utilities.waterBaseProp_ph()

执行此操作，这就是您收到错误消息的原因。 解决这个问题基本上有两种可能：

1. 您更改模型以摆脱或修改约束方程(没有未知数的方程)。它应该是错误消息中显示的那个：

   der(pipe_base3_1.rho[2]) = ...

2. 您将微分功能添加到介质中(我对流体/介质不太感兴趣，所以我不知道这有多复杂，所以我会尝试先使用 1.如何做到这一点 https://modelica.org/documents/ModelicaSpec33Revision1.pdf 在第 12.7 节中显示

下面是示例的代码：

model PipeTest
Pipe_base3 pipe_base3_1(redeclare package Medium = Modelica.Media.Water.WaterIF97_R1pT)
annotation (Placement(transformation(extent={{-10,-10},{10,10}})));
Modelica.Fluid.Sources.FixedBoundary boundary(
nPorts=1,
p=100000,
redeclare package Medium = Modelica.Media.Water.WaterIF97_R1pT)
annotation (Placement(transformation(extent={{-60,-40},{-40,-20}})));
Modelica.Fluid.Sources.FixedBoundary boundary1(
nPorts=1,
p=100000,
redeclare package Medium = Modelica.Media.Water.WaterIF97_R1pT)
annotation (Placement(transformation(extent={{60,-40},{40,-20}})));
Modelica.Thermal.HeatTransfer.Sources.FixedHeatFlow fixedHeatFlow[2](Q_flow={0,0})
annotation (Placement(transformation(extent={{-40,20},{-20,40}})));
equation 
connect(boundary.ports[1], pipe_base3_1.fluidin) annotation (Line(points={{-40,-30},{-8,-30},{-8,-9}}, color={0,127,255}));
connect(boundary1.ports[1], pipe_base3_1.fluidout) annotation (Line(points={{40,-30},{8,-30},{8,-9}}, color={0,127,255}));
connect(fixedHeatFlow.port, pipe_base3_1.heatport) annotation (Line(points={{-20,30},{0,30},{0,7}}, color={191,0,0}));
annotation (
Icon(coordinateSystem(preserveAspectRatio=false)),
Diagram(coordinateSystem(preserveAspectRatio=false)),
uses(Modelica(version="3.2.2")));
end PipeTest;

希望这有帮助...