Fluent-UDF造波源程序集

法一：边界造波法程序一：inlet.c#includeudf.h/*#includestdio.h(输入输出)*/#includestdlib.h#includemath.h#defineHH0.5/*waverhight*//*不同波形需要修改的波形参数*/#defineLL0.95/*波长*/#defineg9.81#definepi3.1415925#defined0.5/*waterdeepth*/DEFINE_PROFILE(x_velocity,thread,index){realkk=2.0*pi/LL;realc=sqrt(g*tanh(kk*d)/kk);/*参见一般性公式的推导*/realTT=LL/c;realx[ND_ND];/*thiswillholdthepositionvector*/realy=0;reals=0;realct=0;/*相位角*/face_tf;realt=CURRENT_TIME;realu=0;t=RP_Get_Real(flow-time);begin_f_loop(f,thread)/*loopsoverallfacesinthethreadpassedintheDEFINEmacroargument*/{F_CENTROID(x,f,thread);y=x[1];s=y+d;ct=kk*(x[0]-c*t);/*参见一般性公式的推导*/if(y0.5*HH*sin(ct))/*水面以下，其中，0.5*HH*sin(ct)为波面方程*/u=pi*HH*cosh(kk*s)*sin(ct)/(TT*sinh(kk*d));/*x方向速度分量公式，参见一般性公式的推导*/elseu=0.0;/*水面以上流体单位速度矢量的x方向分量*/F_PROFILE(f,thread,index)=u;}end_f_loop(f,thread)}DEFINE_PROFILE(y_velocity,thread,index){realkk=2.0*pi/LL;realc=sqrt(g*tanh(kk*d)/kk);realTT=LL/c;realx[ND_ND];/*thiswillholdthepositionvector*/realy=0;reals=0;realct=0;face_tf;realt=CURRENT_TIME;realv=0;t=RP_Get_Real(flow-time);begin_f_loop(f,thread){F_CENTROID(x,f,thread);y=x[1];s=y+d;ct=kk*(x[0]-c*t);if(y0.5*HH*sin(ct))v=pi*HH*sinh(kk*s)*cos(ct)/(TT*sinh(kk*d));elsev=0.0;F_PROFILE(f,thread,index)=v;}end_f_loop(f,thread)}DEFINE_PROFILE(voffactor,thread,index){realkk=2.0*pi/LL;realc=sqrt(g*tanh(kk*d)/kk);realTT=LL/c;realx[ND_ND];/*thiswillholdthepositionvector*/realy=0;reals=0;realct=0;face_tf;realt=CURRENT_TIME;t=RP_Get_Real(flow-time);begin_f_loop(f,thread){F_CENTROID(x,f,thread);y=x[1];s=y+d;ct=kk*(x[0]-c*t);if(y0.5*HH*sin(ct))F_PROFILE(f,thread,index)=1.0;elseF_PROFILE(f,thread,index)=0.0;}end_f_loop(f,thread)程序二：wave.c（此程序同程序一大致相同）#includeudf.h/*mustbeatthebeginningofeveryUDFyouwritenew9case*/realAA=0.5;/*waveramplitude*/realLL=0.95;/*不同波形需要修改的波形参数*/realTT=0.78;realpi=3.1415926realkk=2.0*pi/TT;realww=2.0*pi/LL;realh=0.5;/*waterdeepth*/realux=1.0;/*此为何变量？*/DEFINE_PROFILE(x_velocity,thread,index){realx[ND_ND];/*thiswillholdthepositionvector*/realy;face_tf;realt=CURRENT_TIME;realu;begin_f_loop(f,thread)/*loopsoverallfacesinthethreadpassedintheDEFINEmacroargument*/{F_CENTROID(x,f,thread);y=x[1];if(y(AA*cos(kk*x[0]-ww*t)))/*(-a*sin(w*t)+2.0为入口的波面随时间的变化*//*u=a*w*(exp(k*(y+h))+exp(-k*(y+h)))*-sin(w*t)/(exp(k*h)-exp(-k*h))*/u=9.8*AA*kk/ww*cosh(kk*(y+h))*cos(kk*x[0]-ww*t)/cosh(kk*h);elseu=0.0;/*水面以上流体单位速度矢量的x方向分量*/F_PROFILE(f,thread,index)=u;}end_f_loop(f,thread)}DEFINE_PROFILE(y_velocity,thread,index){realx[ND_ND];/*thiswillholdthepositionvector*/realy;face_tf;realt=CURRENT_TIME;realv;t=RP_Get_Real(flow-time);begin_f_loop(f,thread)/*loopsoverallfacesinthethreadpassedintheDEFINEmacroargument*/{F_CENTROID(x,f,thread);y=x[1];if(y(AA*cos(kk*x[0]-ww*t)))v=9.8*AA*kk/ww*sinh(kk*(y+h))*sin(kk*x[0]-ww*t)/cosh(kk*h);elsev=0.0;F_PROFILE(f,thread,index)=v;}end_f_loop(f,thread)}DEFINE_PROFILE(voffactor,thread,index){realx[ND_ND];/*thiswillholdthepositionvector*/realy;face_tf;realt=CURRENT_TIME;begin_f_loop(f,thread)/*loopsoverallfacesinthethreadpassedintheDEFINEmacroargument*/{F_CENTROID(x,f,thread);y=x[1];if(y(AA*cos(kk*x[0]-ww*t)))F_PROFILE(f,thread,index)=1.0;elseF_PROFILE(f,thread,index)=0.0;}end_f_loop(f,thread)}法二：推波板造波法（动网格）程序三：pushboard.c#includestdio.h#includeudf.h#defineT5.8/*周期*/#defineS1.04/*冲程*/#definePI3.1415926DEFINE_CG_MOTION(pushboard,dt,cg_vel,cg_omega,time,dtime){realu=0;realww=0;ww=2*PI/T;if(time=2*T)u=ww*S*time*cos(ww*time)/(4*T);elseu=ww*S*cos(ww*time)/2;cg_vel[0]=u;}法三：摇板造波法程序四：11.c（三维造波）#includeudf.h#definek0.77#definew2.75#defineA0.059/*振幅*/#definel1.6/*波长*/DEFINE_CG_MOTION(yaoban,dt,vel,omega,time,dtime){vel[0]=0.0;/*x方向的速度*/vel[1]=0.0;vel[2]=0.0;omega[0]=0.0;omega[1]=0.0;omega[2]=atan(A/l)*w*cos(w*time);/*z方向的角频率*/}程序五：b.c#includeudf.hstaticrealG=9.81;staticrealVLA_M=100;staticrealVLA_F=0.0;staticrealVLA_V=0.0;DEFINE_CG_MOTION(vla,dt,vel,omega,time,dtime){Thread*t;Domain*d;realdv,CG[ND_ND],force[2],moment[2];if(!Data_Valid_P())return;NV_S(vel,=,0);NV_S(omega,=,0);d=THREAD_DOMAIN(DT_THREAD((Dynamic_Thread*)dt));t=DT_THREAD(dt);NV_S(CG,=,0.0);Compute_Force_And_Moment(d,t,CG,force,moment,FALSE);dv=dtime*(force[1]-G*VLA_M)/VLA_M;VLA_V+=dv;vel[1]=VLA_V;VLA_F=force[1];Message(time=%fF(y)=%fV(y)=%f\n,time,force[1],VLA_V);}