基于数值模拟的复杂地形风场风资源评估方法

《基于数值模拟的复杂地形风场风资源评估方法》由会员上传分享，免费在线阅读，更多相关内容在教育资源-天天文库。

万方数据第30卷第3期空气动力学学报V01．30，No．32012年6月ACTAAERODYNAMICASINICAJun．，2012文章编号：0258—1825(2012)03一0415-06基于数值模拟的复杂地形风场风资源评估方法梁思超1，张晓东1，康顺1’2，康雅兰3，赵永锋4(1．华北电力大学电站设备状态监测与控制教育部重点实验室，北京1022062．西安现代控制技术研究所，陕西西安710065；3．中国电力工程顾问集团华北电力设计院工程有限公司，北京100120；4．中国福霖风能工程有限公司，北京100034)摘要：与采用线性模型的WAsP软件相比，CFD方法具有可以逼真模拟复杂地形三维风场大气流动的优势，因此成为了复杂地形风场风资源评估的发展方向。为了把CFD软件包FINE／TURBO和带有壁面函数的}￡湍流模型用于复杂地形地貌风电场的风资源评估和微观选址中去．本文以南澳岛风场作为研究对象，在中性大气条件假设下，对风场区域以30。风向为间隔进行数值模拟。研究采用不同入口边界条件对数值计算结果的影响。比对模拟结果发现，风场中的风加速与湍流强度基本上不随来流条件发生变化。基于这一规律，结合全年测风数据评估风场的风资源分布状况，建立了一套风资源评估的数值模拟及后处理方法，为工程实际提供了有力参考。关键词：湍流模型；微观选址；数值模拟；风资源；大气边界层中图分类号：TK89文献标识码：A0引言风电场风资源评估和微观选址的好坏将直接关系到风电场发电量的多少。经过大量工程实践，基于线性模型的WAsP软件被认为不适用于复杂地形风场。计算流体力学(CFD)作为一种很有前途的手段，在这一领域正不断发展。湍流数值模拟方法主要有直接模拟(DNS)、大涡模拟(LES)、分离涡模拟(DES)以及雷诺时均模拟(RANS)。从精度上看，在实际山地Bolund山测风实验的盲评中，排名靠前的均为两方程的雷诺时均数值模拟方法，优于大涡模拟和风洞实验结果。综合计算精度，运算时间，计算机资源，雷诺时均法还没有受到其他任何模拟方法的挑战[1]。带有壁面函数的k-e湍流模型，由于同时考虑了风速，湍动能以及地表粗糙度，适于风工程的研究。为把这一用于工业流动问题的模型用于风工程研究，前人对湍流模型参数进行了修正[2_3]。2009年在NUMECA软件平台下，提出其他相关参数对模拟的影响，作出进一步修正n]。使用该修正模型的Ask—ervein山数值模拟结果很好地与实验值吻合。为了将这一模拟方法用于实际风场的风资源评估和微观选址中去，以实际南澳岛风场为研究对象。结合测绘数据和卫星数据等间隔生成不同风向下的计算域，进行全周数值模拟。1南澳岛简介为了简化问题，降低边界条件设置难度，选取岛屿地形南澳岛风场作为研究对象。南澳岛濒临台湾海峡，处于南海北部，距离汕头市莱芜岛约10km，面积106．45km2，岛上山地面积占93．6％。岛上平均风速达8．44m／s，年有效发电小时数达7215h，年平均有效风能密度为678W／m2，有东北和西南两个主风向，风力资源居世界最佳之列。图1南澳岛Fig．1Nanaoisland*收稿日期：201l一01—08；修订日期：201l-02—26基金项目：国家国际科技合作计划资助(2010DFA64600)作者简介：梁思超(1984一)，男，辽宁沈阳人，博士研究生，主要从事：风电场微观选址的数值研究 万方数据416空气动力学学报第30卷瓣嘲i一■22·24口20-22口18-20口16·t8口14．16■12-】4■10一12●8-10●6-8■4．6●2—4■0-2图2固定测风塔位置图36号塔风向玫瑰图Fig．2FixedmastspositionsFig．3Windroseof6#南澳风电场资源勘测期间，设立了若干固定测风塔和移动测风塔。其中，l、2、3、4、6、7号测风塔为固定测风塔，测有全年风资源数据，位置如图2所示。从6号测风塔70m高度处风向玫瑰图(图3)中可以看到，南澳岛的主风向为40。左右。2数值方法采用NUMECA的FINE／TURB0软件包，该软件采用时间相关法求解湍流Navier—Stokes方程、中心节点的有限体积离散、显式Runge—Kutta法，全多重网格初场处理以及多重网格迭代加速收敛。选用的湍流模型为带有壁面函数的肛￡湍流模型。2．1k-e两方程湍流模型标准肛￡湍流模型三维标准输送方程为：云(川+砉(咖，)一蠹『(户+等h，ak，7，+Pt—IDe(1)易(一+丢(J0EUj)一蠹[(产+等)妾卜c，Pk--C吐P警@’k和￡为湍动能及耗散率，湍动能生成项R为：P。一一_1)面警(3)式(3)中，U，和U：分别为j方向的平均风速和脉动速度。湍流粘度为：肚一f0(：。L(4)湍流模型参数参数C。、C。、cE。、靠、瞑在模拟大气边界层时设置为C。=0．03，C．=1．21，C。2—1．92，诉一1．0，匪=1．3【⋯。为排除软件对流场中卢。／P的影响，设置MUCI。IP=l×109。数值迭代时，由于湍动能是从入口不断向出口发展的，迭代速度很慢，不利于计算的收敛。将参数KEGRID的值设置为大于全多重网格加速重数的最大值。此时，湍动能k会在最粗网格上进行迭代，为较细网格的迭代提供好的初场，加速收敛，防止发散。计算中KEGRID取4。2．2壁面函数在FINE／TURBO软件中，不考虑粘性底层的存在，粗糙度‰为粗糙地表上风速为零点到粗糙地表顶端的距离。壁面函数的表达式如下：兰：土1n兰型+B(5)“此Z“其中，“’为摩擦速度，与壁面剪切应力r删有关，“’一~／‰||／J0。卡门常数*一般取0．4l。z。为表面粗糙度，d为零偏移量，定义“(‰+d)一0。当d≠0时，相当于把壁面移到距真实壁面距离d的地方。对于粗糙壁。A，>70时，常数B取零。其中，女?一A；_y+／y，∥’=yu+／v。计算中取d=0，B=0。3计算域与边界条件设置实际风工程中，由于经费有限，只测绘风场内有限的地形区域。而在CFD中，为了考虑周围地形对风场空气流动的影响，需要更大区域的地形数据来生成计算域和网格。针对这一问题，提出结合网上可免费获取的SRTM卫星地形数据，补全测绘地形数据的方案。卫星数据绝对水平和高程精度分别为20m和16m，仅作为地形补充，精度足够。南澳岛风场二期的测绘地形数据与卫星数据如图4、图5所示。地形处理时考虑了WGS84与北京54坐标的转换、二维插值、两地形间的平滑¨]。图4测绘地形图图5SRTM卫星地形图Fig．4MappingterraindataFig．5SRTMterraindata以30。为间隔，分别进行12个风向的数值模拟。图6为0。、30。、60。风向下的计算域和网格。由于风向正交，如0。、90。、180。、270。，这四个风向采用相同的计算域和网格。于是，模拟中仅需要三个计算域．根据风向设置相应的边界条件(即风速的大小和方 万方数据第3期粱思超等：基于数值模拟的复杂地形风场风资源评估方法向)即可。另做40。主风向计算域，这样就有4个不同的计算域。一◆_t图6风向角为0。、30。、60。的网格和计算域Fig．6Gridsanddomainsof0。。30。and60。degree计算域尺寸为18000m×18000m×4000m。水平方向网格单元尺寸为90m×90m。计算域z、y、z三个方向的网格数为200、200、76。垂直方向网格分布分为三层结构：底层、中层和顶层(图7)。底层网格共高lOm，共有3个网格单元，按等延展比分布，第一层网格高度2．3m。中层内每个网格高度相等，均为5m，共12个。顶层网格也按等延展比分布，随着高度的增加网格单元尺度增大。由于中层和底层内每个网格单元具有固定的高度，在后处理时，可以在一定高度上生成一个与地面起伏趋势相同的网格面，便于分析模拟结果。以40。风向计算为例(图8)，设置边界条件时，来流方向边界为人口边界，对面一侧为出口边界，两侧为镜像边界，计算域顶端为入口边界，地面为粗糙固壁。将计算域划分为若干块，便于在不同固壁区域设置不同的地表粗糙度。设置海面粗糙度为0．002m，陆地粗糙度为0．3m。=j=：：三：j中层‘：i：：：．=：：=二顶层基量jj_!‘’三i；厘挥受i蒌誊熏i垂黧蠹i蠢i盔旗和底层zLY图7垂直方向网格分布图8边界条件设置Fig．7VerticalgriddistributionFig．8Boundaryconditions入1：1边界上速度、湍动能、以及耗散率由以下模型确定‘“。t一芫c，一丢，2㈣￡一生(8)膣其中，心取0．4，C。一0．03。h。为大气边界层高度，h。一“’／6f，厂为科里奥利参数，，一2nsina(1AI)，其中，n为地球自转的角速度，A为风场所处纬度位置。认为大气边界层顶端为湍动能消失的地方，当高度超过大气边界层高度时，速度不再变化，湍动能和耗散率为零。计算中，n一7．292×10一rad，A取23．5。，‰取0．002m。以50m为特征高度，以风速2、4、6、8、10、12、14、16、18、20m／s作为特征值，根据式(6)分别计算“’，生成来流条件。4模拟结果分析4．1不同来流条件对计算结果的影响为了研究不同来流风速对于模拟结果的影响，引入无量纲量风加速AS，定义式为：△s一坠粤掣(9)U，。八nJUk(^)为距地表高h处的当地风速，Uref(h)为参考风速。取50m高处的风速进行比较。各向同性时，湍流强度JT与湍动能k和风速U的关系为：七一号u2IT2(10)以o。和90。风向计算为例，分别选取6个位置(如图9)。其中，A点和B点位于背风坡；C点和D点位于山顶；E点和F点位于迎风坡。图10和图11分别为各个位置50m高度上的风加速和湍流强度与入口特征风速的关系。虽然各点位置不同，流动状况也不同．但所处位置上的风加速和湍流强度基本上不随来流发生变化，其余风向模拟结果均存在这一规律。这对于后处理带来了极大的方便，每个方向仅作一个风况的特征计算即可。图90。和90。风向下特征位置Fig．9Characteristicpositionsof0。and90。 万方数据418空气动力学学报第30卷1．n0j宅＆∞．0．-1．1A0j’口＆”·0．1仉．A．B．C，DE。F9t．A．B·C，DEvFEj7{·．·一。。-·04；：：；：：；：；；皇5j-。‘‘‘‘‘‘‘‘‘‘3J5{11‘．⋯⋯⋯．{‘········’100々々々0々々00¨o—w夺sipl}／1(m-矿s")夯扩一inrflowls广peeldr／(m矿．s)inflowspeedl(m：snflowspee／I-．s-i面图100。风向下特征位置风加速和湍流强度Fig．10Speedupandintensityof0。01．A。B．C，DE。F3．01-A．B·C，DE，F5{赫．～--·-·．．o{：zI?zz：z：z巴撼5j⋯⋯⋯·o’孙；；：；；；；；；I嗡i蠡赢纛≯以56i骊酽图1190。风向下特征位置风加速和湍流强度Fig．11Speedupandintensityof90。4．2计算结果比较由于以上规律的存在，可以认为：同风向下、不同风速时，风场不同位置两处风速比例关系不变。因此，在一个风向下，当知道某一位置的风速U．时，可采用公式(11)计算另一位置风速U2。其中，U，，和U妲为该风向下，对应位置模拟值。rU，一=娶尘(11)U，l为验证数值模拟的精度，进行各测风塔间的风速互相推算比较。除6号塔外，其余各塔1、2、3、4、7均有50m高度处的数据。以7号塔风向为风向基准，在风向和风速稳定的时间段内选取3个小时平均数据，见表l。由式(11)，以各测风塔风速计算其它测风塔位置处风速，并比较差异，得到表2～表4。表1样本风速Table1Samplewindspeeddata200211280710031．85．16．07．06．77．72002122703：0060．87．18．79．211．510，12002092320：0090，510．010．77．69．78．5表230。风向计算结果比较Table230。resultscomparisons表360。风向计算结果比较Table360。resultscomparisons邀!：：：：18．7％29．8％39．2％17．5％29．5％23．2％33．5％9．7％342．5％30．2％13．4％17．5％464．5％50．3％15．S％35．7％721．2％10．7％14．9％26．3％表490。风向计算结果比较Table490。resultscomparisons从三个风向的互推差异来看(表2、表3、表4)，数值模拟结果存在较大误差。产生误差的原因可能来自两方面。一是数值模拟网格的太稀，不能很好地表现地面的起伏状况，简单的壁面函数不足以描述复杂的粗糙地表。二是由于安装，损坏，附近障碍物影响等原因造成的测量数据不准。但无论怎样，这一方法在风资源评估当中是值得借鉴的。4．3主风向下的流场分析及微观选址由于40。为该风场的主风向，研究该风向下的流场分布状况是十分必要的。图12中标出的7个位置，综合空气流动状况，风能密度和湍流强度，适于安装风力机。最下面的区域由于面积较大，地形平坦，便于施工，应优先考虑。图12适于安装风机的位置及风速、风能密度、湍流强度分布Fig．12Suitablepositionsandwindspeed，windenergydensity，turbulenceintensitydistribution蝴嘣咧蝴m吼乱眠蛾慨跏坜撕＆船舱咪聊螂慨乱孔地置嘣慨班螂孤％■趴喘嘣蜥掰住“¨玑蛳瑞洲慨己硒阻吼班螂慨咧m拍儿&％％7O82蚰"∞∞％％％％6438勰加柚拍嘣吲咪班蚰能u乏 万方数据第3期梁思超等：基于数值模拟的复杂地形风场风资源评估方法419风场区域内有4个位置(图13)处于背风处，流动状况复杂，风能密度小，流动不平稳。尤其是最下那个区域，处于高山的背风坡，十分陡峭，又由于左侧山地的阻挡作用，一股气流斜向流入，加剧了该处的流动分离。这些位置的湍流强度明显大于周围，因此，要避免在这些位置安装风力机。图14年平均风能密度云图Fig．14Annualaveragewindenergydensitycontour6结论在相同计算域下，地表附近风加速和湍流强度模拟值基本上不随入口参考风速的变化而发生变化。因此，当知道风场某一位置的风速和湍流强度时，可根据模拟值计算出另一位置的风速和湍流强度。基于这一思想，建立起的数值模拟和后处理方法对于实际风场的风资源评估具有积极的指导作用。图13不适于安装风机的位置及风速、风能密度、湍流强度分布参考文献Fig．13Unsuitablepositionsandwindspeed，windenergydensity，burbulenceintensitydistribution5全周模拟结果的后处理实际工程当中，为了简化计算，多采用威布尔分布来拟合测风数据。描述风速的概率分布。但很多情况下，风速概率分布不满足威布尔分布。这时，可采用处理时间序列的方式来评估风资源。即，每一个测风数据均进行累计计算。处理时，选择一位于高处、受地形影响较小的测风塔作为基准测风塔，用来判定该时刻的风况属于哪个风向的数值计算。在风场区域内布置二维网格面，利用式(11)计算节点i采用第j个测风塔的数据点得到的风速【，。。由于测风塔往往有n个，可采用反距离加权法计算该结点的风速U。，即：”T”一U，=∑鲁／∑万1(12)i一1uqi=1¨l】其中d。为节点i到测风塔j的距离。用每一时刻测风数据计算得到风能密度，累计相加后平均，即得到每个网格点上的年平均风能密度。得到风场区域年平均风能密度云图，如图14。HANJALICK，KENJERESS．Somedevelopmentsinturbulencemodelingforwindandenvironmentalengi—neering[-J'1．JournalofWindEngineeringandIndus—trialAerodynamics，2008，96(10一11)：1537—1570．BECHMANNA，S()RENSENNN，joHNSENS，eta1．HybridRANS／LESmethodforhighReynoldsnum—bers，appliedtOatmosphericflowovercomplexterrain[A]．TheScienceofMakingTorquefromWind[C]．Lyngby，DENMARK，2007．No．012054．YANGY，GUM，CHENSQ，eta1．Newinflowboundaryconditionsformodellingtheneutralequilibri—umatmosphericboundarylayerincomputationalwindengineering[J]．J．ofWindEngineeringandIndus—trialAerodynamics，2009，97：88—95．梁思超，张晓东，康顺．复杂地形风场扰流数值模拟方法[J]．工程热物理学报．2011，32(6)：945—948．梁思超．复杂地形风场数值处理及扰流模拟方法研ED]．[硕士论文]．北京：华北电力大学，2010．ZHANGXD．CFDsimulationofneutralABLflowsrR]．Ris口-R-1688(EN)，2009． 万方数据420空气动力学学报第30卷AssessmentmethodofwindresourcesincomplexterrainbasedonnumericalsimulationLIANGSi-cha01，ZHANGXiao—don91，KANGShunl”，KANGYa-lan3，ZHA0Yong-fen94(1．KeyLab．ofStateSupervisionandControlofPowerPlantEquipment，NorthChinaElectricPorwerUniversity．Beij’ing102206，China}2．Xi’口nModernControlTechnologyResearchInstitute，xi'an710065，China；3．NorthChinaPowerEngineeringCo．，LtdofChinaPowerEngineeringConsultingGroup，Beijing100120．China，4．ChinaFulinWindPotwrEngineeringCo．·Led。Beijing100034，China)Abstract：CFDmethodhastheadvantageofwellsimulatingtheflowoverwindfarmofcomplexterrain，comparedwiththemethodofusingthesoftwareWAsPwithlinearmodel，thusbecomingthedevelopmentdirectionofassessmentofwindresourcesincomplexterrain．InordertousetheCFDpackageFINE／TURBOandstandardk-eturbulencemodelwithwalIfunctioninthewindresourcesassessmentandmicrositingofwindfarminthecomplexterrain，inthispaper，IslandNanaohadbeensimulatedinwinddirectionsatinter-valsof30degrees，undertheassumptionofneutralcondition，tOinvestigatetheinfluenceofinletconditionontheresults．Thesimulationresultsshowthatthewindspeed—upandintensitydon’tchangeessentiallyunderdifferentinletconditions．Basedonthislaw，withwholeyearwinddata，anassessmentmethodwithpostprocessinghasbeenbuilt，providingastrongreferenceforengineeringpractice．Keywords：turbulencemodel；micrositing；numericalsimulation；windresource；atmosphericboundarylayer●●h●_·◆¨I-●IlH■4·◆__◆悼_●●州⋯●●忡。●-◆。I’◆_。●}q’●_’’-忡-I’'_’●_-◆¨●_“●一I-◆_·◆．_●·“◆_’¨■■_}●一●·◆川H■q’◆悻-_-●_-●●●q“◆●。◆●◆I“．一n删～一¨¨劓～-●__敬告读者《空气动力学学报》互联网期刊平台已正式开通使用。该平台集作者投稿、专家审稿、审稿进度查询和刊物网上订阅等功能于一体。作者可登录进行在线投稿和状态查询。审稿专家可登录在线审稿。网址：http：／／kqdlxxb．cars．org．cn。欢迎各位专家、作者和读者访问!今后，本刊将不再接受纸质稿件和电子邮件投稿。作者在线投稿时，须先进行注册登记，除提交文章电子版(PDF格式)并签署版权转让协议外，还须上传保密审查单扫描件。详情请参考网页上的《投稿指南》。欢迎大家踊跃投稿!专家在线审稿时，可点击个人邮箱审稿邀请函直接进入审稿页面，也可在地址栏输入网址，然后在主页面的左上角登录区按注册信息登录后审稿。如有其它疑问，请与编辑部联系。电话(传真)：0816—2463138。感谢大家的关注与支持!《空气动力学学报》编辑部2012年6月20日 基于数值模拟的复杂地形风场风资源评估方法作者：梁思超，张晓东，康顺，康雅兰，赵永锋，LIANGSi-chao，ZHANGXiao-dong，KANGShun，KANGYa-lan，ZHAOYong-feng作者单位：梁思超,张晓东,LIANGSi-chao,ZHANGXiao-dong(华北电力大学电站设备状态监测与控制教育部重点实验室,北京,102206)，康顺,KANGShun(华北电力大学电站设备状态监测与控制教育部重点实验室,北京102206;西安现代控制技术研究所,陕西西安710065)，康雅兰,KANGYa-lan(中国电力工程顾问集团华北电力设计院工程有限公司,北京,100120)，赵永锋,ZHAOYong-feng(中国福霖风能工程有限公司,北京,100034)刊名：空气动力学学报英文刊名：ActaAerodynamicaSinica年，卷(期)：2012,30(3)参考文献(6条)1.HANJALICK;KENJERESSSomedevelopmentsinturbulencemodelingforwindandenvironmentalengineering2008(10-11)2.BECHMANNA;SORENSENNN;JOHNSENSHybridRANS/LESmethodforhighReynoldsnumbers,appliedtoatmosphericflowovercomplexterrain[外文期刊]20073.YANGY;GUM;CHENSQNewinflowboundaryconditionsformodellingtheneutralequilibriumatmosphericboundarylayerincomputationalwindengineering20094.梁思超;张晓东;康顺复杂地形风场扰流数值模拟方法2011(06)5.梁思超复杂地形风场数值处理及扰流模拟方法研20106.ZHANGXDCFDsimulationofneutralABLflows[Risφ-R-1688(EN)]2009本文链接：http://d.g.wanfangdata.com.cn/Periodical_kqdlxxb201203023.aspx