采用大涡模型进行风工程计算的方法1

风工程计算指导（第一版）-内部资料ADINA技术部2011年11月1前言风工程计算中流场尺度大，建筑风工程的流场计算域可达百米甚至数百米，模型单元数量多，收敛过程需要的迭代步数多。由此导致的计算量大以及实际存在的收敛困难是必须要克服的问题。在ADINA中，全耦合流场的控制体积或者有限元计算方法尽管精度非常高，但由于风工程模型网格数量多（10~1000万），流场自由度全耦合方法会耗费大量时间。我们意见一般采用FCBI-C算法，即SimpleC算法，在空间离散和时间积分上选择二阶精度算法，还包括其它一些关键建模环节的合理处理，可实现高效，高精度的风工程计算。以下我们用一个方块障碍的风洞试验模拟来说明主要问题。示例模型说明方块障碍物尺寸为1*1*1米，风洞高度2米。计算模型初步尺寸设定为，上游取3米，下流6米，侧面3米。网格划分尺度参加附带模型。初步计算的网格单元数量为57408。分析一：常规瞬态风场计算（入口速度不变）第一步采用K-E湍流模型,获得平均速度下的稳态流场风工程计算指导（第一版）-内部资料ADINA技术部2011年11月2流态。只计算一步，输出map文件给后续计算使用，作为后续计算的初始流场状态。主要参数设置如下：***********************************************************SOLUTIONPARAMETERS-Step1***********************************************************总体参数设置：分析类型STEADY-STATE，不可压缩，采用K-E湍流模型，输出map文*件MAP-OUTP=YES*MAP-FORM=YES；采用FCBI-C算法。MASTERANALYSIS=STEADY-STATEMODEX=EXECUTETSTART=0.00000000000000,IDOF=1TURBULEN=KEHYDRO=YESSTREAM=YESTRACTB=DEFAULT,IRINT=DEFAULTAUTOMATI=NOSOLVER=AMG1COMPRESS=NOFSINTERA=NO,NMASS=0MASSCOUP=NOMAP-OUTP=YESMAP-FORM=YESNONDIMEN=YES,MAXSOLME=0MTOTM=2RECL=3000ALE=NOTHERMAL-=NO,UPWINDIN=CONTROL-VOLUMEMESHUPDA=CURRENTMESHADAP=NO,COUPLING=ITERATIVEPOROUS-C=NOCELL-BCD=YESVOF=NOFCBI=CENTER,TURB-ITE=COUPLEDEM-MODEL=NO*无量纲化：注意长度单位为流场核心区域单元平均边长，速度为入口平均速度。NONDIMENSIONINPUT=YESOUTPUT=NOL=.05,XREF=.0YREF=.0ZREF=.0,V=20.0D=1.0CP=1.0,T=1.0TREF=.0*方程残差为0.0001，方程类型为质量方程；参数残差0.0001，参数类型为所有参数；速度*压力场耦合算法SIMPLEC；空间离散精度一阶First；最大迭代次数999。OUTER-ITERATEPSER=1.E-4ERTYPE=MASS,EPSVR=1.E-4VRTYPE=ALLMAXITE=999MINITE=5,PVCOUP=SIMPLECORDER=FIRSTFTSTYPE=NOFTS=1.E+21,VPT-MAXI=1VP-MAXIT=1KE-MAXIT=1VIDPRT=DEFAULTITSMAX=0,MS-SOLVE=DEFAULT*松弛因子：速度0.35,压力0.85，湍流动能0.4,耗散率0.4。OVAR-CONTROLRELAX-PRI-SSVELOCITY=0.35,PRESSURE=0.85TEMPERAT=0.95,TURB-K=0.4TURB-E=0.4*****************************************************************INNERITERATIONSCALES-内迭代过程设置****************************************************************最少迭代次数：速度4，压力25，湍流5。IVAR-CONTROLITE-MINVELOCITY=4PRESSURE=25TEMPERAT=2TURB-K=5*方程残差：速度0.01，压力0.01，湍流动能0.0001。风工程计算指导（第一版）-内部资料ADINA技术部2011年11月3IVAR-CONTROLRED-EQNVELOCITY=0.01,PRESSURE=0.01TEMPERAT=0.1,TURB-K=1.E-4*变量残差：速度0.01，压力0.01，湍流动能0.0001。IVAR-CONTROLRED-VARVELOCITY=0.01,PRESSURE=0.01TEMPERAT=0.1,TURB-K=1.E-4**************************************************************MATERIALPPROPERTIES**************************************************************注意入口施加了0.005的扰动分量。MATERIALTURBULENT-KENAME=1XMU=1.E-5,C1=1.44C2=1.92C3=0.8,CMU=0.09SIGMAK=1.0,SIGMAE=1.3SIGMAT=0.9DW=70.0,VONK=0.4CP=0.0,XKCON=0.0BETA=0.0QB=0.0,RHO=1.0TREF=0.0,GRAV-X=0.0GRAV-Y=0.0,GRAV-Z=0.0SIGMA=0.0TYPE=STANDARD,KAPPA=1.0E+20CV=0.0MDESCRIP='NONE'*注意模型其它设置1入口速度剖面，ADINA中用空间函数1~4定义。2入口施加湍流载荷，一般可输入平均速度和湍流强度让ADINA自动计算K，E数值；3设置输出TXT格式的map文件，见总体控制中参数。4对称边界条件采用Wall模拟，wall设置为sliding类型。风工程计算指导（第一版）-内部资料ADINA技术部2011年11月4第一步命令流文件original-ke.in。第二步采用大涡模型进行风场模拟，要求时间积分和空间离散采用二阶精度格式。主要设置如下：***********************************************************SOLUTIONPARAMETERS***********************************************************总体参数设置：分析类型TRANSIENT，其它参数不变。MASTERANALYSIS=TRANSIENTMODEX=EXECUTETSTART=0.00000000000000,IDOF=1TURBULEN=NOHYDRO=YESSTREAM=YESTRACTB=DEFAULT,IRINT=DEFAULTAUTOMATI=NOSOLVER=AMG1COMPRESS=NOFSINTERA=NO,NMASS=0MASSCOUP=NOMAP-OUTP=YESMAP-FORM=YESNONDIMEN=YES,MAXSOLME=0MTOTM=2RECL=3000ALE=NOTHERMAL-=NO,UPWINDIN=CONTROL-VOLUMEMESHUPDA=CURRENTMESHADAP=NO,COUPLING=ITERATIVEPOROUS-C=NOCELL-BCD=YESVOF=NOFCBI=CENTER,TURB-ITE=COUPLEDEM-MODEL=NO*TIMESTEPNAME=DEFAULT@CLEAR5000.2@*时间积分采用二阶精度格式。ANALYSISTRANSIENTALPHA=0.707106781186547METHOD=COMPOSI*NONDIMENSIONINPUT=YESOUTPUT=NOL=.05,XREF=.0YREF=.0ZREF=.0,V=10.0D=1.0CP=1.0,T=1.0TREF=.0*空间离散阶次second。其它参数不变。OUTER-ITERATEPSER=1.E-5ERTYPE=MASS,EPSVR=1.E-4VRTYPE=ALLMAXITE=180MINITE=5,PVCOUP=SIMPLECORDER=SECONDFTSTYPE=NOFTS=1.E+21,VPT-MAXI=1VP-MAXIT=1KE-MAXIT=1VIDPRT=DEFAULTITSMAX=0,MS-SOLVE=DEFAULT*参数不变。OVAR-CONTROLRELAX-PRI-SSVELOCITY=0.35,风工程计算指导（第一版）-内部资料ADINA技术部2011年11月5PRESSURE=0.85TEMPERAT=0.95,TURB-K=0.4TURB-E=0.4****************************************************************INNERITERATIONSCALES-内迭代设置-不变***************************************************************IVAR-CONTROLITE-MINVELOCITY=4PRESSURE=25TEMPERAT=2TURB-K=5*IVAR-CONTROLRED-EQNVELOCITY=0.01,PRESSURE=0.01TEMPERAT=0.1,TURB-K=1.E-4*IVAR-CONTROLRED-VARVELOCITY=0.01,PRESSURE=0.01TEMPERAT=0.1,TURB-K=1.E-4**************************************************************MATERIALPPROPERTIES**************************************************************风场瞬态计算采用大涡模型MATERIALTURBULENT-LESNAME=1XMU=1.E-5,SGS=SMAGRANDOM=YESINTENSITY=0.005INTERVAL=5,CP=0.0XKCON=0.026BETA=0.0,QB=0.0XKD=0.1TPR=1.0,RHO=1.0TREF=0.0,GRAV-X=0.0GRAV-Y=0.0,GRAV-Z=0.0SIGMA=0.0,KAPPA=1.0E+20CV=0.0*模型其它设置：1删除入口施加湍流载荷；2读入上一步输出的map文件到当前数据库中。INITIAL-MAPPFILENAME='step1.map'EXTERNAL=ALL@CLEAR'X-VELOCITY''Y-VELOCITY''Z-VELOCITY''PRESSURE'@*3由于计算步数多，结果文件可能巨大，因此控制只输出节点结果。PORTHOLESAVEDEFA=NOFORMATTE=NOINPUT-DA=1NODESAVE=YES风工程计算指导（第一版）-内部资料ADINA技术部2011年11月6此模型命令流为original-les.in。第三步（可选）一般大涡模型的发展需要过程，因此需要更好的结果可能需要模拟更多的物理时间才能够获得，尤其当尺度比较大的情况。观察第二步获得流场的流动发展情况，即关注点的压力或者速度的平均值和周期性是否完好，如果完好则不必执行第三步的计算。