SUNTANS模型学习(3)——学习cylinder算例

(38) 2024-05-06 08:01:01

学习cavity算例

  • 简介
  • 网格配置
  • 参数配置
    • Input file for SUNTANS部分
    • Grid Files部分
    • Output Data Files和Input Data Files部分
    • Use output to calculate Strouhal number部分
    • 其余参数
  • 模型运行
  • 模拟结果

简介

Cylinder Flow(圆柱绕流)是流体数值模型的经典算例。对于平面对称布置的圆柱绕流工况,在理想状态下,圆柱两侧绕流流线呈对称状。但随着雷诺数的增大,圆柱上下游的流线会逐渐失去对称性,而在圆柱后侧的流场中产生卡门涡街。卡门涡街是流体力学中重要的现象,对于卡门涡街的模拟也常被用于模型检验。
在suntans模型代码包中,内置有一圆柱绕流算例cylinder。模拟的配置文件位于 /examples/cylinder 文件夹中。在本算例中,模型预设了两套网格配置(三角形网格、三角形和矩形混合的网格),以及垂向平均流动与三维流动的参数配置。模拟工况汇总如下:

工况 网格配置 二维浅水模拟or三维流动模拟
1 三角形网格 二维浅水模拟
2 三角形网格 三维流动模拟
3 三角形、矩形混合网格 二维浅水模拟
4 三角形、矩形混合网格 三维流动模拟

通过模拟本算例,主要学习参数文件suntans.dat中各主要参数的设定及其物理含义。

网格配置

两套网格配置如下图所示(上图为三角形网格,下图为三角形、矩形混合网格):
SUNTANS模型学习(3)——学习cylinder算例 (https://mushiming.com/)  第1张
其中的红点为各个网格的Voronoi点。
如果模拟三维流动,则需要设定垂向网格配置。在本算例中,我们设定垂向网格总数为Nkmax=20;垂向网格上疏下密(垂向网格由suntans.dat文件中的rstretch参数控制,这个下一节中细讲)。

参数配置

对于本算例,模型参数由suntans.dat文件给出;文件中的参数共包含以下5个部分:

  1. Input file for SUNTANS:主要是网格与物理参数
  2. Grid Files:网格文件名称
  3. Output Data Files:输出文件的名称
  4. Input Data Files:输入文件名称
  5. Use output to calculate Strouhal number in mfiles/strouhal.m:用于计算Strouhal数的输出结果

Input file for SUNTANS部分

这是这个参数文件的主体部分,这里我们以三角形网格的三维流动模拟为例,看看各个控制参数。以下就是各参数,输入方式为【参数名】+【空格】+【参数值】;#号后侧的为注释;即

variablename variablevalue # Description of variable

以下是suntans.dat中第一部分的内容,在此附上中文说明:

Nkmax			20	# 最大的垂向网格数字,以此为例,垂向网格最多的地方有20层网格。
stairstep		0	# 该值可取1或0;1表示每层垂向网格间距在水平方向上不变,不过这也有可能导致某处最低网格底面低于实际深度;
                    #0表示在底面附近网格使用partial-stepping,使单元的底面与深度的实际值相一致。
rstretch		-1.05	# 垂向网格拉伸系数,取值范围(-1.1, 1)或(1, 1.1);当 r<0,底部网格得到了加密,反之,则加密顶上的网格
CorrectVoronoi		-1	# 关于是否启用三角网格修正,包含两种方法,即该值有两个取值:1或-1。
                        # 当取1时,对Voronoi点距离过近的网格进行校正;当取-1时,将钝角三角形的Voronoi点移动到它的形心。
VoronoiRatio		0	# Voronoi点的校正因子。当不启用CorrectVoronoi时,该值取0;
                        # 当CorrectVoronoi=1,该值取[0,1];当CorrectVoronoi=-1,该值取角度[0,90];(详见手册)
vertgridcorrect 	0 	#  对最小水深D_min<Δz的情况进行校正;当该值为1,遇此情况时有Δz=D_min,从而使Nkmax = D_max /D_min
IntDepth 		0	# 如果取1,Voronoi点的水深由depth文件给出(空间插值);否则,水深在initialization.c中设定。
dzsmall			0.1	# 校正前网格最小间距比【这个参数不再使用,不必设置】
scaledepth 		0 	# 对水深的缩放(一般用于测试);当取值为1,则该功能启动
scaledepthfactor 	0 	# 水深缩放比例因子,在scaledepth=1的时候生效
thetaramptime	        1  	# ramp time(详见手册),给初始化时一个缓慢的物理量变化过程,来抑制瞬态的正压效应;在本例子中,该值取1s。
beta			0	# 盐度扩散系数(本例中取0)
theta			0.55	# 控制时间数值积分的显式/隐式/半隐式的系数(对于水动力控制方程)
# θ = 0: 显式(不稳定);θ = 0.5: Crank-Nicolson ;θ = 1.0: 隐式;θ = 0.55: 常用设置。
thetaS			0.55	# 控制时间数值积分格式的系数(对于对流-扩散方程)
thetaB			0.55	# 【不再使用的参数】
kappa_s                 0       # 盐度的垂直扩散系数
kappa_sH                0       # 盐度的水平扩散系数(详见手册)
gamma 			1 	# 温度的扩散系数
kappa_T                 0       # 温度的垂直扩散系数
kappa_TH                0       #  温度的水平扩散系数(详见手册)
nu 			5e-4   	# 水的垂直方向粘度 (m^2 s^-1) (对于层流)
nu_H 			5e-4 	# 水的水平方向粘度 (m^2 s^-1) (对于层流)
tau_T			0 	# 风应力
z0T	                0.0 	# 水面边界的粗糙度 (单位:m),用于计算CdT
z0B	                0.0	# 底面边界的粗糙度 (单位:m),用于计算CdB
CdT	                0.0 	# 如果z0T=0,则顶部的阻力系数由这个常数给出;否则,这个参数不生效。如果用free-slip条件,则取-1。
CdB	                -1	# 如果z0B=0,则顶部的阻力系数由这个常数给出;否则,这个参数不生效。如果用free-slip条件,则取-1。
CdW			0.0	# 侧向边界阻力系数;只有当壁面没有延伸到自由表面时,才计算侧壁上的阻力(例如漫滩情况下,不计算侧向阻力)
turbmodel		0	# 紊流模型 (0 - 不使用, 1 - MY25)
AB          3 # Adam-Bashforth法的阶数(2或3,对应2阶或3阶A-B方法)
dt 			2.5e-3	# 时间步长
grav			9.81e10 # 重力加速度(默认为9.81,此算例中采用了一个较大值以控制雷诺数)
Cmax 			1 	# Courant数的最大允许值
nsteps 			1600    # 时间步总数
ntout   		80  	# 输出文件的步频,此处表示每80个dt输出一次
ntprog   		5 	# 报告运算进度的频率 (%)
ntconserve 		1	# 输出 conserved data 的步频(例如质量、体积、势能等),详见phys.c中函数OutputData
nonhydrostatic		1	# 0 - 应用静水压假定, 1 = 非静压模拟
cgsolver		1	# 线性方程组的求解方法 0 = GS, 1 = CG 【这个参数已不再使用】
maxiters		1000	# CG 法(共轭梯度法)迭代的最大次数
qmaxiters		2000	# CG 法迭代的最大次数 (泊松方程)
qprecond		2	# 是否对泊松方程进行预处理(0 - 不处理;1 - Diagonal预处理;2 - Block-diagonal预处理)
epsilon			1e-10 	# CG法中收敛性判别所允许的最小值
qepsilon		1e-5	# CG法中收敛性判别所允许的最小值(压力泊松方程)
resnorm			0	# 在计算CG求解器的公差准则时是否对残差进行归一化(0和1对应了不同的误差形式)
relax			1	# CG求解器的松弛系数
amp 			0.8	# 用于指定boundaries.c中的参数prop->amp (振幅)
omega 			0	# 用于指定boundaries.c中的参数prop->omega (圆频率)
flux 			0	# 用于指定boundaries.c中的参数prop->flux (通量或流量)
timescale		0	#  用于指定boundaries.c中的参数prop->timescale(应用于开边界的时间尺度)
volcheck		0	# 是否检查体积守恒(0-不检查;1-检查)
masscheck		0	# 是否检查质量守恒(0-不检查;1-检查)
interp			1	# 插值方法(0-经典Perot方法;1-二次插值)
nonlinear		2	# 动量方程对流项格式(0-无对流;1-迎风格式;2-中心差分格式)
laxWendroff		0       # Lax Wendroff扩散系数(0-不应用Lax Wendroff法;1-当nonlinear=2时,涡粘系数通过Lax Wendroff法给出)
newcells		0	# 1 if adjust momentum in surface cells as the volume changes, 0 otherwise
wetdry			0       # 是否采用干-湿网格设置 (1-使用, 0-不使用)
Coriolis_f              0	# 科氏力参数 f=2*Omega*sin(phi)
sponge_distance	        0 	# 海绵边界层衰减距离(详见phys.c)
sponge_decay	        0	# 海绵边界层衰减时间系数 
readSalinity		0	# 0-垂直盐度分布在initilization.c中设定;1-从InitSalinityFile中读取垂直盐度分布,文件必须为二进制格式,且必须按照从上到下的顺序包含N_kmax个双精度值。
readTemperature		0	# 与readSalinity类似,当取1时,从InitTemperatureFile中读取垂直温度分布

Grid Files部分

网格文件均是ASCII的二进制文件。suntans.dat文件中需要设定它们的文件名。

pslg   full_circle.dat 	# 平面直线图形文件;详见手册第3章。 (input) 【本例中不采用】
points points.dat	# 网格节点文件 (input)
edges edges.dat		# 网格边配置文件 (input)
cells cells.dat		# 网格连接关系文件 (input)
depth depth.dat		# 水深文件 (若INTERPDEPTH=1,则应用该文件) (input)
nodes nodes.dat   # 节点的拓扑关系 (output) 
celldata celldata.dat	# 每个线程的网格中心数据 (output)
edgedata edgedata.dat	# 每个线程的网格边界数据t (output)
vertspace vertspace.dat	# 垂向网格间距 (output)
topology topology.dat	# 每个线程的网格拓扑关系

各个文件的格式详见手册第7章。

Output Data Files和Input Data Files部分

########################################################################
#  Output Data Files
########################################################################
FreeSurfaceFile   	fs.dat   #自由水面高程
HorizontalVelocityFile 	u.dat   #水平流速
VerticalVelocityFile 	w.dat   #垂向流速
SalinityFile 		s.dat  #盐度场数据
BGSalinityFile 		s0.dat  #
TemperatureFile		T.dat  #温度
PressureFile		q.dat  #动水压力
VerticalGridFile 	g.dat  #垂向网格
ConserveFile		e.dat  #守恒性检查
ProgressFile	        step.dat  #模拟进度
StoreFile		store.dat  #保存数据,用于热启动的文件
StartFile		start.dat  #用于热启动的文件
EddyViscosityFile	nut.dat  #涡粘系数
ScalarDiffusivityFile	kappat.dat  #标量扩散系数
########################################################################
# Input Data Files 【本例中不采用】
########################################################################
InitSalinityFile	sinit.dat  #初始盐度场
InitTemperatureFile	Tinit.dat  #初始温度场

Use output to calculate Strouhal number部分

ProfileVariables	u	         # 输出的变量值
DataLocations		dataxy.dat       # 变量所在位置的x,y坐标
ProfileDataFile		profdata.dat     # 垂向剖面信息
ntoutProfs		10               # 每次输出配置文件数据所间隔的时间步;此处含义为每10个时间步输出配置文件数据
NkmaxProfs		1                # 输出文件的最顶层
numInterpPoints		1                # 是否”仅输出最邻近的点“

注:这个部分在手册中也并未详细描述;在此先略过。

其余参数

在目录/examples/cylinder下输入以下指令以运行三角形网格浅水流动:

make test

指令框中将会显示一些参数未被设定,但模型运行时自动采用了其默认值。详细信息如下:

Warning in MPI_GetValue...retrieved default value of maxFaces=3.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of ntoutStore=0.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of TVDsalt=4.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of TVDtemp=4.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of TVDturb=0.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of TVDmomentum=3.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of conserveMomentum=1.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of thetaM=-1.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of mergeArrays=1.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of computeSediments=0.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of calcage=0.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of agemethod=1.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of calcaverage=0.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of latitude=29.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of gmtoffset=0.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of metmodel=0.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of varmodel=1.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of nugget=0.10 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of sill=0.90 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of range=100000.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of outputNetcdf=0.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of netcdfBdy=0.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of readinitialnc=0.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of Lsw=2.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of Cda=1100.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of Ce=1400.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of Ch=1400.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of hprecond=1.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of prettyplot=0.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of linearFS=0.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of fixdzz=1.00 which was not specified in data/suntans.dat.
Warning in MPI_GetValue...retrieved default value of smoothbot=0.00 which was not specified in data/suntans.dat.

这里面的参数包括:

maxFaces=3.00 # 网格的最大边界数,默认应用三角形网格,故取3
ntoutStore=0.00 # 
TVDsalt=4.00 # 求解盐度控制方程所选用的TVD格式
TVDtemp=4.00 # 求解温度控制方程所选用的TVD格式
TVDturb=0.00 # 求解紊流模型控制方程所选用的TVD格式
TVDmomentum=3.00 # 求解动量方程所选用的TVD格式
conserveMomentum=1.00 # 是否进行动量的守恒性检查
thetaM=-1.00 # 求解动量方程所选用的TVD格式
mergeArrays=1.00 # 求解垂向对流项的格式;当大于0.5表示隐式,当等于-1表示在UPredictor()中使用了原来的保守方案
computeSediments=0.00 # 是否求解泥沙模块
calcage=0.00 # 是否模拟水龄
agemethod=1.00 # 水龄模拟的方式 1-使用河流边界,2-使用内部源项
calcaverage=0.00 # 是否计算平均水龄

##以下几个参数用于热量和太阳辐射模块##
latitude=29.00 # 纬度,在计算太阳辐射时候使用
gmtoffset=0.00 # 用于修正太阳辐射的参数
metmodel=0.00 # 关于太阳辐射热量模型 0-无气象数据输入,1-COARE3.0,计算辐射短波和长波
varmodel=1.00 # 太阳辐射模型中所用插值方法 0-反距离比例插值,2-克里金法
nugget=0.10 # 克里金插值中所需给定的参数
sill=0.90 # 克里金插值中所需给定的参数
range=100000.00 # 变差函数参数范围
outputNetcdf=0.00 # 输出文件格式 0-二进制格式,1-netcdf格式
netcdfBdy=0.00 # netcdf文件的输出步频
readinitialnc=0.00 # 是否采用netcdf格式的初始条件
Lsw=2.00 # Light extinction depth (m)
Cda=1100.00 # 阻力系数
Ce=1400.00 # 传热系数
Ch=1400.00 # 传热系数
##########################

hprecond=1.00 # 对于自由水面方程求解食肉采用预处理求解器 0-不使用,1-Jacobi预处理求解器
prettyplot=0.00 # 是否使用二次插值得出的输出值,若为1则使用;否则,请指定
linearFS=0.00 # 对自由水面求解是否采用线性方法
fixdzz=1.00 # 对底网格高度的处理 0-不调整,1-假设底部单元高度必须大于dzsmall*dz[Nkmax-1],2-假设底部单元高度必须大于dzsmall
smoothbot=0.00 # 在SetFluxHeight和ComputeVelocityVector中是否进行平滑处理

模型运行

考虑到有四个模拟工况,模型代码中预载了下面四个指令来分别运行。

  1. 三角形网格二维浅水流动模拟
make test
  1. 三角形网格三维流动模拟
make test3d
  1. 混合网格二维浅水流动模拟
make test-quad
  1. 混合网格三维流动模拟
make test-quad-3d

注意:在每次模拟前,记得先使用make clobber取消编译;之后再使用上述指令进行编译+运行;运行结果生成在 /examples/cylinder/data中。
还需要注意的是使用make clobber会自动删除 /examples/cylinder/data 这个目录的文件夹(删除掉这个data文件夹),所以每次运行的结果记得先存好。

模拟结果

由上一节的内容可知,模型的时间步长dt=2.5e-3s,共运行1600步;ntout=80,故每个文件包含了 1+floor(Nsteps/Ntout) = 21个时间步的流场信息。
利用代码包/mfiles文件夹下的getcdata.m读取二进制文件,可以获得每一时间步及每个网格的流场数据,使用getcdata函数的示例参考SUNTANS模型学习(2)——用MATLAB工具提取Time Accuracy算例的计算数据。
以下展示二维浅水模拟结果中形成的卡门涡街(水平流速云图,T = 1600*dt):

  1. 三角形网格
    SUNTANS模型学习(3)——学习cylinder算例 (https://mushiming.com/)  第2张SUNTANS模型学习(3)——学习cylinder算例 (https://mushiming.com/)  第3张
  2. 混合网格
    SUNTANS模型学习(3)——学习cylinder算例 (https://mushiming.com/)  第4张
    SUNTANS模型学习(3)——学习cylinder算例 (https://mushiming.com/)  第5张
THE END

发表回复