-1

我正在使用与代数方程混合的抛物线 PDE,并且所有这些方程都是耦合的。我在 OM(1.9.1)("residualFcn[some number]") 中使用了欧拉方法(Dassl 太慢)和大容差(用于快速模拟)和接收错误(两种类型)。问题是求解器不能t 求解非线性系统(数学上,系统是正确的)。第一个问题是什么类型的方法在 OM 中使用欧拉积分法(显式或隐式或 Crank-Nicholson ......或......)?所以我试图用数字解决它(显式欧拉方法(在“new [N]”下面的代码中)(也许问题可能是CFL条件。),但我有问题(特定采样时间的样本重建)。所以,第二个问题是指重现特定采样时间的值?!在下面的代码中有数组“a [3]”。想法是针对每个“

还有一件事,如果 delta(t)/(delta(x))^2 >= 0.5 (delta(t) 定义用户,并指代方程部分,则 delta(x) 如下面的代码所示,并使用空间离散化在方程部分(经典前馈方法)),是否满足数值稳定性?同样的例子,但算法部分?问候

这是代码:

model Euler1D
import Modelica.Utilities.*;
parameter Integer N=10; //50
parameter Real Lp=1e-6;
parameter Real deltax=1/(N-1)*Lp;
Real a[3];
Real old[N];
Real new[N];
Real b;
equation 
a[1]=if 
   (time>5) then 0 else time+5;
a[2]=time;
a[3]=2;
when 
(sample(0,1)) then
d=b;
end when;
algorithm 
// IN t=ts;
 when (sample(0,1)) then
 for i in 1:2 loop
  b:=a[i];
 Streams.print(String(time)+"   "+String(a[i])+ "     "+String(b), "C:/Some_Path/text.txt");
end for;
end when;

// Another problem
old[1]:=10;
old [N]:=0;
new[1]:=10;
new [N]:=0;
// Boundary
for i in 2:N-1 loop
old [i]:=10;
new[i]:=10;
end for;

for dx in deltax:deltax:Lp-deltax loop  // spatial discretization

for i in 2:N-1 loop
(new[i]):=(old[i]+0.5*(old[i + 1] +old[i-1]- 2*old[i]));
//def:=def+abs(new[i]-old[i]);
end for;
 for i in 2:N-1 loop
 old[i]:=new[i]; // switch the values
 end for;

 for i in 1:N loop
 Streams.print(String(time)+"   "+ String(new[2]), "C:/Some_Path/Anel_Nodes.txt");
 end for;

 annotation (uses(Modelica(version="3.2")));
 end Euler1D;
4

1 回答 1

3

鉴于您的评论,这就是我构建模型的方式:

model Euler1D 
  parameter Integer N=10 "Spatial discretization";
  parameter Modelica.SIunits.Length L=1;
protected 
  parameter Modelica.SIunits.Length dx=L/(N-1) "Segment size";
  Real c[N] "Solution variable c";
  Real J[N] "Solution variable J";
  Real dc_dx "Spatial derivative at x==L";
  Real d2c_dx2[N] "Second spatial derivative of c";
initial equation 
  der(c[2:N-1]) = zeros(N-2);
equation 
  // Equations for spatial derivatives
  d2c_dx2[2:N-1] = { (c[i+1]-2*c[i]+c[i-1])/(dx^2) for i in 2:N-1};
  dc_dx = (c[N]-c[N-1])/dx;

  // Boundary conditions
  c[1] = 0 "Dirichlet B.C.";
  dc_dx = 1 "Neumann B.C.";

  // PDE
  J = c .* c "J = c^2";
  der(c) = d2c_dx2+J "PDE";
end Euler1D;

初始方程避免了模拟开始时的大瞬态。如此有效,我提出的模型将为您提供稳态解决方案。但是您可以使事情随时间变化(例如边界条件)。

这仍然可能存在一些错误,但我很难说,因为我对这个系统没有太多的直觉,而且我也不知道我使用的值是否在物理上有意义。

但希望它概述了如何构建这样一个模型。该模型在 Dymola 中运行。我没有在 OpenModelica 中尝试过。

话虽如此,你总是用 PDE 来解决这个问题,即存在时间积分方案无法明确看到的物理产生的隐式约束(例如 CFL 条件)。一般来说,可变时间步长算法具有误差估计,它们可以间接看到违反这些隐式约束所引入的不稳定性。但他们必须“摸索”通过它们。我不确定上述模型在多大程度上会出现问题。

如果你想给你的解增加一点刺激,你可以改变方程

  dc_dx = 1 "Neumann B.C.";

...至...

  dc_dx = 0.7+0.2*sin(time) "Neumann B.C.";

...你会得到一个很好的时变解决方案。我注意到,在玩这个的时候,如果你把末端的梯度设置得太陡,模型就会变得不稳定。由于我上面提到的问题,这是对基础数学方程的真正解决方案还是数值伪影,我不知道。因此,如果您想要一个稳定的解决方案,您可以在此模型中添加的内容是有限的。

于 2013-02-27T15:16:05.860 回答