bolang.rar_Stokes_二阶stokes_波浪udf_波浪udf模拟_波浪模拟

#include "udf.h" #include <math.h> #define PI 3.141592653 /*圆周率*/ #define D 5 /*水深*/ #define H 1 /*波高*/ #define T 4 /*周期*/ #define L 23 /*波长*/ #define FL 90 /*水槽长度*/ //x direction velocity// DEFINE_PROFILE(x_velocity, thread, position) { real x[ND_ND]; real t,y; real yu,k,w,Lu;/*波速、圆频率*/ face_t f; begin_f_loop(f,thread) { t = CURRENT_TIME; F_CENTROID(x,f,thread); y=x[1]; yu=y-D; k = 2*PI/L; w = 2*PI/T; Lu=D+(H/2)*cos(-1*w*t)+(PI*H*H/(8*L))*cosh(k*D)/(pow(sinh(k*D),3))*(cosh(2*k*D)+2)*cos(-2*w*t); if(y<=Lu) F_PROFILE(f,thread,position) =(PI*H/T)*(cosh(k*(yu+D))/sinh(k*D))*cos(-1*w*t)+0.75*(PI*PI*H*H/(T*L))*(cosh(2*k*(yu+D))/pow(sinh(k*D),4))*cos(-2*w*t); else F_PROFILE(f,thread,position) =0; } end_f_loop(f,thread) } //y direction velocity// DEFINE_PROFILE(y_velocity, thread, position) { real x[ND_ND]; real t,y; real yu,k,w,PP;/*波速、圆频率*/ face_t f; begin_f_loop(f,thread) { t = CURRENT_TIME; F_CENTROID(x,f,thread); y=x[1]; yu=y-D; k = 2*PI/L; w = 2*PI/T; PP=D+(H/2)*cos(-1*w*t)+(PI*H*H/(8*L))*cosh(k*D)/(pow(sinh(k*D),3))*(cosh(2*k*D)+2)*cos(-2*w*t); if(y<=PP) F_PROFILE(f,thread,position) =(PI*H/T)*(sinh(k*(yu+D))/sinh(k*D))*sin(-1*w*t)+0.75*(PI*PI*H*H/(T*L))*(sinh(2*k*(yu+D))/pow(sinh(k*D),4))*sin(-2*w*t); else F_PROFILE(f,thread,position) =0; } end_f_loop(f,thread) } //入口体积分数udf文件// DEFINE_PROFILE(inlet_voffactor, thread, position) { real x[ND_ND]; real t,y; real yu,k,w,Lu;/*波速、圆频率*/ face_t f; begin_f_loop(f,thread) { t = CURRENT_TIME; F_CENTROID(x,f,thread); y=x[1]; yu=y-D; k = 2*PI/L; w = 2*PI/T; Lu=D+(H/2)*cos(-1*w*t)+(PI*H*H/(8*L))*cosh(k*D)/(pow(sinh(k*D),3))*(cosh(2*k*D)+2)*cos(-2*w*t); if(y<=Lu) F_PROFILE(f,thread,position) =1; else F_PROFILE(f,thread,position) =0; } end_f_loop(f,thread) } //出口体积分数udf文件// DEFINE_PROFILE(outlet_voffactor, thread, position) { real x[ND_ND]; real y; face_t f; begin_f_loop(f,thread) { F_CENTROID(x,f,thread); y=x[1]; if(y<=D) F_PROFILE(f,thread,position) =1; else F_PROFILE(f,thread,position) =0; } end_f_loop(f,thread) } //x direction momentum// DEFINE_SOURCE(x_mom_source,cell,thread,dS,eqn) { real x[ND_ND]; real x_source; real t; real y; real xishu; t = CURRENT_TIME; y = x[1]; C_CENTROID(x,cell,thread); if(x[0]>=FL-L && x[0]<=FL) { if(y>=D-H && y<=D) { xishu =10*(x[0]-(FL-L))/L; x_source = -1*xishu*C_R(cell,thread)*C_U(cell,thread); dS[eqn] = -1*xishu; } } return x_source; } //y direction momentum// DEFINE_SOURCE(y_mom_source,cell,thread,dS,eqn) { real x[ND_ND]; real y_source; real t; real y; real xishu; t = CURRENT_TIME; y = x[1]; C_CENTROID(x,cell,thread); if(x[0]>=FL-L && x[0]<=FL) { if(y>=D-H && y<=D) { xishu = 10*(x[0]-(FL-L))/L; y_source = -1*xishu*C_R(cell,thread)*C_V(cell,thread); dS[eqn] = -1*xishu; } } return y_source; } //出口压力udf文件// DEFINE_PROFILE(outlet_pressure, thread, position) { real x[ND_ND]; real y; face_t f; begin_f_loop(f,thread) { F_CENTROID(x,f,thread); y=x[1]; if(y<=D) F_PROFILE(f,thread,position) =998.2*9.81*(D-y); else F_PROFILE(f,thread,position) =0; } end_f_loop(f,thread) }