二维、结构网格单相流SIMPLE 算法的C语言实现

来源：互联网发布：淘宝消保在哪里交编辑：程序博客网时间：2024/05/03 05:38

下载：SIMPLE_C源代码 exam1 exam6 exam10 exam6.mpg(46K)

SIMPLE_C说明文件(PDF 1.7M)

与我联系：qzhen@flotu.org

第一节前言

SIMPLE算法的全称是“Semi-Implicit Method for Pressure-Linked Equations”，意思是求解压力耦合的质量/动量/能量传递方程的半隐式方法。自从S.V.Patankar(帕坦卡)和D.B.Spalding(斯伯丁)于1972年提出SIMPLE算法后，它实际上已经成为许多工程流动、传热以及反应体系的数值模拟的最重要的方法。许多商业CFD软件，如cfx与fluent,其核心也都基于SIMPLE算法。

了解SIMPLE算法是学习差分数值计算以及CFD（计算流体力学）一个非常好的起点。而且，SIMPLE算法对于很多工程问题也足够用了。介绍SIMPLE算法最好的书自然是帕坦卡自己所著的《传热与流体流动的数值计算》[1]。不过，这本优秀、简洁的著作没有附带源代码。

我在学习SIMPLE的过程中，了解到一份可靠、简单、说明文档完整、可方便地编译执行的源代码对于初学者的重要性。于是我用C语言改写了SIMPLE算法(称之为SIMPLE_C)，并给出了详细的文档，希望能够缩短算法与代码之间的距离。算法源自水木清华BBS上orange01与moran两位同学所给予的FORTRAN SIMPLE代码。同时清华大学化工系反应工程教研组刘飞同学的笔记也给了我很大的帮助，在此一并谢过。

目前，SIMPLE_C只能解决二维直角坐标、结构网格的问题。取决于表征问题的微分方程，SIMPLE_C在这个几何框架下能够模拟瞬态或者时均的流动、传热以及反应耦合的单相流体系。应该说明，下面的一些问题可能或者根本不在SIMPLE_C的能力范围之内：高精度的流场预测；多相流；三维流动；复杂流场或者几何结构；边界层行为预测；马赫数接近或超过1的高速流动。

关于orange01与moran两位同学的SIMPLE FORTRAN代码的版权问题，我在源代码上看见这么一段，在此原封不动地摘录如下：

Cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
Cccc This computer program was copied from the graduate student course program
Cccc of the University of Minnesota. Part of it was re-formulated to meet the
Cccc personal computer environment. Some inappropriate expressions were also
Cccc corrected. The program is used only for teaching purpose. No part of it
Cccc may be published. You may use it as a frame to re-develop your own code
Cccc for research purpose. XJTU Instructor, 1995.11

Cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc

　

本文的工作是非常粗浅的。本质上它是帕坦卡在70年代就已解决的问题的一个简单的复述。不过，希望它能起到抛砖引玉的作用――期望更多的研究者能够将手头基于SIMPLE的先进算法的源代码整理并共享出来。下面一些特性是特别欢迎的：

l 三维网格；

l 非结构/自适应网格；

l 高精度差分格式；

l 欧拉+欧拉模式的多相流；

l 欧拉+拉格朗日模式的多相流；

l 基于PC CLUSTER(个人电脑集群，MIMD体系)的并行算法。

我个人在计算方法的兴趣主要集中欧拉+拉格朗日模式的多相流模拟。要解决的主要物理问题是连续相与离散相的耦合作用。

第四节计算实例以及源程序编译指南

4.1 例子的编译

http://www.flotu.org/~qzhen/cfd/simplec/index.htm 给出了三个计算例子的ZIP文件。它们都是VC6或者VC7的项目文件，对应后面的例子4.2~4.4。相信使用GCC的同学能够轻松地写出make文件。编译过程中需要注意的matlab相关部分的编译。下面我以Windows平台介绍一下。关于使用matlab的engine接口的编程的详细介绍可以浏览：

[5] http://www.mathworks.com/access/helpdesk/help/techdoc/matlab_external/ch06eng2.shtml#using_the_matlab_engine。

SIMPLE_C的matlab engine的C语言接口是matlab 5.3中所带的。编译过程的几个简要的技巧叙述如下：

1) 程序中只有exam*.c文件中的user_plot()用到了matlab engine的C语言接口。如果你不想使用它或者编译遇到麻烦，直接将user_plot()以及其它matlab相关语句注释调掉即可。这些语句只在exam*.c出现，包括：

exam*.c中的#include "engine.h"与#include "matrix.h";

exam*.c中的give_matlab_a_double_matrix()与give_matlab_a_interger()的定义;

exam*.c中的user_plot()中的语句;

global_var.h中的give_matlab_a_double_matrix()与give_matlab_a_interger()的声明。

下面的(2)~(4)是需要加入matlab支持的步骤。

2) 把matlab的头文件路径”%matlab%/extern/include”加入到你的编译环境中。%matlab%是你的matlab安装主目录。

3) 产生.lib文件。”%matlab%/extern/include”有产生.lib文件的.def文件。使用”%vc%/bin/lib.exe /def[:deffilename] /out[:libfilename]”将.def编译为.lib文件。%vc%是你的Visual C++的安装主目录。需要产生的.lib文件包括libeng.lib, libmat.lib与libmx.lib。它们所使用的.def文件分别是libeng.def, libmat.def与libmx.def。不同版本的Visual C++使用的.lib文件好像并不兼容。

4) 将libeng.lib, libmat.lib与libmx.lib加入到你的项目中。VC7与VC6的设置见下面的图：

图23. 在VC7中添加matlab engine库文件

图24. 在VC6中添加matlab engine库文件

另外，关于制作matlab 5.3下动画制作的问题

可以用下面的方法

mpgwrite :matlab user contributed function
功能，把一系列matlab的图像帧转化为mpg文件
平台matlab5 :win32,solaris,sgi,hpux,linux
具体可以参考

http://fluid.stanford.edu/~fringer/movies/gif_howto/gif_howto.html
http://www.rpi.edu/~holmes/NumDiffEqs/Demos/matlab.movie.info.pdf

需要从www.mathworks.com下载一个文件,我已经上载到 http://www.flotu.org/~qzhen/cfd/simplec/mpgwrite.tar 2M

4.1 例子exam1：二维稳态传热

边壁上的温度。

扩散系数。扩散系数在user_gamsor( )中设定。

求解域x方向宽度＝1.0 。求解域y方向宽度＝2.0 。

计算结果：

图25. exam1的计算结果－二维温度场的分布

4.2 例子exam10：湍流中被动标量扩散（基于时均的 k –ε方程）

exam10中没有考虑温度对于密度以及粘度的影响，所以可以把温度当成只是一个被动标量。

入口左半边的速度为5m/s，温度为100摄氏度；右半边的速度为10m/s，温度为400摄氏度。

计算结果：

图26. exam10的计算结果－温度与速度的分布

4.3 例子exam6：常粘度系数湍流的动态模拟

一股热流体（速度100m/s，温度500摄氏度）从右下角注入系统（壁面温度为300摄氏度），出口的流速分布设为均匀。

微分方程包括时均的连续性方程、动量、能量的输运方程。所用的粘度系数是1。热传导系数的量级与粘度系数相仿。

图27. exam6的计算结果－温度与速度的分布

动画文件可以从http://www.flotu.org/~qzhen/cfd/simplec/exam6.mpg 下载。

本文原地址：http://www.flotu.org/~qzhen/cfd/simplec/