UDF，非线性雷诺应力，USDI？

benqing

想在fluent中用UDF实现，非线性雷诺应力，能从CFDonline上找到一个相关的代码，CFDonline链接，我在fluent中能够编译成功，但是在开始计算的时候出现问题，报错：

==============================================================================

Node 0: Process 13836: Received signal SIGSEGV.

==============================================================================
MPI Application rank 0 exited before MPI_Finalize() with status 2
 The fl process could not be started.

测试发现是C_UDSI这部分的问题，在我的理解中UDS是用于建立标量方程的，这个里面采用UDSI计算了UU等雷诺应力分量，有考虑在UDS对话框中添加UDS，但是对这些量并没有输运方程，不知道应该如何处理，希望有人帮忙解决，多谢！
代码贴在下面

/* turbulence model: nonlinear k-epsilon turbulence model by Ehrhard (1999) */ /* UDF (User Defined Function) Routine for FLUENT */ 
/* see also the paper: */ /* Implementation and comparison of different turbulence models for three dimensional wall jets with Fluent */ 
/* FLUENT CFD Forum 2005 Bad Nauheim, Germany */ /* Heschl Ch*., Sanz W.**, Klanatsky P.* */ 
/* * FachhochschulstudiengçÏge Burgenland GmbH, SteinamangerstraÅ½ße 21, A-7423 Pinkafeld */ 
/* ** TU Graz, Institut of Thermal Turbomachinery and Machine Dynamics, Inffeldgasse 25, A-8010 Graz */ 
/* Contact: christian.heschl@fh-burgenland.at, wolfgang.sanz@tugraz.at, peter.klanatsky@fh-burgenland.at */

/* Die Autoren ?ernhemen keine Verantwortung oder Haftung f? den Inhalt der unten angef?rten UDF-Routinen */

#include "udf.h" 
#include "math.h"

/* Turbulence model constants */ 
const real C_1=-0.2; 		const real C_2=0.4; 	const real C_5=0.0; 
const real CE_1=1.44; 		const real CE_2=1.92; 
const real SIG_K=1.0; 		const real SIG_E=1.3; 
const real C_MU_ke=0.09; 	const real kappa=0.42; 	const real E_log=9.793;

/* User-defined scalars */ 
enum { UU, VV, WW, UV, UW, VW };

DEFINE_SOURCE(u_source, c, t, dS, eqn) 
{ 
	real source; 
	dS[eqn]=0.0; 
	source= - C_R(c,t) * ( C_UDSI_G(c,t,UU)[0] + C_UDSI_G(c,t,UV)[1] + C_UDSI_G(c,t,UW)[2] );//C_R(c,t) = rho density 

	return source; 
}

DEFINE_SOURCE(v_source, c, t, dS, eqn) 
{ 	real source;

	dS[eqn]= 0.0; source= - C_R(c,t) * ( C_UDSI_G(c,t,UV)[0] + C_UDSI_G(c,t,VV)[1] + C_UDSI_G(c,t,VW)[2] );

	return source; 
}

DEFINE_SOURCE(w_source, c, t, dS, eqn) 
{	 real source;

	dS[eqn]= 0.0; source= - C_R(c,t) * ( C_UDSI_G(c,t,UW)[0] + C_UDSI_G(c,t,VW)[1] + C_UDSI_G(c,t,WW)[2] );

	return source; 
}

DEFINE_TURBULENT_VISCOSITY(user_mu_t,c,t) 
{
	real S,W,c_mu,T,mu_t; real S11, S12, S13, S21, S22, S23, S31, S32, S33; real W11, W12, W13, W21, W22, W23, W31, W32, W33;

	S11 = 0.5*( C_DUDX(c,t)+C_DUDX(c,t) ); S12 = 0.5*( C_DUDY(c,t)+C_DVDX(c,t) ); S13 = 0.5*( C_DUDZ(c,t)+C_DWDX(c,t) ); 
	S21 = 0.5*( C_DVDX(c,t)+C_DUDY(c,t) ); S22 = 0.5*( C_DVDY(c,t)+C_DVDY(c,t) ); S23 = 0.5*( C_DVDZ(c,t)+C_DWDY(c,t) ); 
	S31 = 0.5*( C_DWDX(c,t)+C_DUDZ(c,t) ); S32 = 0.5*( C_DWDY(c,t)+C_DVDZ(c,t) ); S33 = 0.5*( C_DWDZ(c,t)+C_DWDZ(c,t) ); 
	W11 = 0.5*( C_DUDX(c,t)-C_DUDX(c,t) ); W12 = 0.5*( C_DUDY(c,t)-C_DVDX(c,t) ); W13 = 0.5*( C_DUDZ(c,t)-C_DWDX(c,t) ); 
	W21 = 0.5*( C_DVDX(c,t)-C_DUDY(c,t) ); W22 = 0.5*( C_DVDY(c,t)-C_DVDY(c,t) ); W23 = 0.5*( C_DVDZ(c,t)-C_DWDY(c,t) ); 
	W31 = 0.5*( C_DWDX(c,t)-C_DUDZ(c,t) ); W32 = 0.5*( C_DWDY(c,t)-C_DVDZ(c,t) ); W33 = 0.5*( C_DWDZ(c,t)-C_DWDZ(c,t) );

	T=C_K(c,t)/C_D(c,t);//C_K(c,t) = turbulent kinitic energy(TKE); C_D(c,t)  = epsilon

	S=T*sqrt( 2*(S11*S11+S12*S12+S13*S13+S21*S21+S22*S22+S23*S23 +S31*S31+S32*S32+S33*S33) ); 
	W=T*sqrt( 2*(W11*W11+W12*W12+W13*W13+W21*W21+W22*W22+W23*W23 +W31*W31+W32*W32+W33*W33) );

	c_mu=MIN( 1./(0.9*pow(S,1.4)+0.4*pow(W,1.4)+3.5) , 0.15 );

	mu_t=c_mu*T*C_K(c,t);

	return mu_t; 
}

DEFINE_ADJUST(rsm_adjust,domain) 
{

	Thread *t; cell_t c; real T; real mu_t, c_mu; 
	real S11, S12, S13, S21, S22, S23, S31, S32, S33; 
	real W11, W12, W13, W21, W22, W23, W31, W32, W33; 
	real S, W, Sij, Wij; 
	real C_3, C_4, C_6, C_7; 
	real P11, P22, P33, Pkk; 
	real tau_w, u_tauw, u_star, y_star, Gk, u_mag;

/* Set the turbulent viscosity */ 
	thread_loop_c(t,domain) 
	if (FLUID_THREAD_P(t)) 
	{ 
		begin_c_loop(c,t) 
			{
			S11 = 0.5*( C_DUDX(c,t)+C_DUDX(c,t) ); S12 = 0.5*( C_DUDY(c,t)+C_DVDX(c,t) ); S13 = 0.5*( C_DUDZ(c,t)+C_DWDX(c,t) ); 
			S21 = 0.5*( C_DVDX(c,t)+C_DUDY(c,t) ); S22 = 0.5*( C_DVDY(c,t)+C_DVDY(c,t) ); S23 = 0.5*( C_DVDZ(c,t)+C_DWDY(c,t) );
			S31 = 0.5*( C_DWDX(c,t)+C_DUDZ(c,t) ); S32 = 0.5*( C_DWDY(c,t)+C_DVDZ(c,t) ); S33 = 0.5*( C_DWDZ(c,t)+C_DWDZ(c,t) ); 
			W11 = 0.5*( C_DUDX(c,t)-C_DUDX(c,t) ); W12 = 0.5*( C_DUDY(c,t)-C_DVDX(c,t) ); W13 = 0.5*( C_DUDZ(c,t)-C_DWDX(c,t) ); 
			W21 = 0.5*( C_DVDX(c,t)-C_DUDY(c,t) ); W22 = 0.5*( C_DVDY(c,t)-C_DVDY(c,t) ); W23 = 0.5*( C_DVDZ(c,t)-C_DWDY(c,t) ); 
			W31 = 0.5*( C_DWDX(c,t)-C_DUDZ(c,t) ); W32 = 0.5*( C_DWDY(c,t)-C_DVDZ(c,t) ); W33 = 0.5*( C_DWDZ(c,t)-C_DWDZ(c,t) );

			T=C_K(c,t)/C_D(c,t);

			S=T*sqrt( 2*(S11*S11+S12*S12+S13*S13+S21*S21+S22*S22+S23*S23 +S31*S31+S32*S32+S33*S33) ); 
			W=T*sqrt( 2*(W11*W11+W12*W12+W13*W13+W21*W21+W22*W22+W23*W23 +W31*W31+W32*W32+W33*W33) );

			Sij=S11*S11+S12*S12+S13*S13+S21*S21+S22*S22+S23*S23+S31*S31+S32*S32+S33*S33; 
			Wij=W11*W11+W12*W12+W13*W13+W21*W21+W22*W22+W23*W23+W31*W31+W32*W32+W33*W33;

			c_mu=MIN( 1./(0.9*pow(S,1.4)+0.4*pow(W,1.4)+3.5) , 0.15 );

			mu_t=c_mu*T*C_K(c,t);

			C_3=2.0-exp(-SQR(S-W)); 
			C_4=-32.0*SQR(c_mu); 
			C_6=-16.0*SQR(c_mu); 
			C_7=16.0*SQR(c_mu);

			C_UDSI(c,t,UU)= C_1*mu_t*T* (S11*S11+S12*S12+S13*S13-1./3.*Sij) + 2.*C_2*mu_t*T*(W12*S12+W13*S13) + C_3*mu_t*T*(W12*W12+W13*W13-1./3.*Wij) - 2.*C_4*mu_t*T*T * (W12*(S11*S12+S12*S22+S13*S23)+W13*(S11*S13+S12*S23+S13*S33)) + C_6*mu_t*T*T* (S11*Sij) + C_7*mu_t*T*T* (S11*Wij);//User-Defined Scalar(UDS) cell variables

			C_UDSI(c,t,VV)= C_1*mu_t*T* (S12*S12+S22*S22+S23*S23-1./3.*Sij) + 2.*C_2*mu_t*T*(W23*S23-W12*S12) + C_3*mu_t*T*(W12*W12+W23*W23-1./3.*Wij) + 2.*C_4*mu_t*T*T * (W12*(S11*S12+S12*S22+S13*S23)-W23*(S12*S13+S22*S23+S23*S33)) + C_6*mu_t*T*T* (S22*Sij) + C_7*mu_t*T*T* (S22*Wij);

			C_UDSI(c,t,WW)= C_1*mu_t*T* (S13*S13+S23*S23+S33*S33-1./3.*Sij) - 2.*C_2*mu_t*T*(W13*S13+W23*S23) + C_3*mu_t*T*(W13*W13+W23*W23-1./3.*Wij) + 2.*C_4*mu_t*T*T * (W13*(S11*S13+S12*S23+S13*S33)+W23*(S12*S13+S22*S23+S23*S33)) + C_6*mu_t*T*T* (S33*Sij) + C_7*mu_t*T*T* (S33*Wij);

			C_UDSI(c,t,UV)= C_1*mu_t*T* (S11*S12+S12*S22+S13*S23) + C_2*mu_t*T*(W12*(S22-S11)+W13*S23+W23*S13) + C_3*mu_t*T*(W13*W23) + C_4*mu_t*T*T* (W12*(S11*S11+S13*S13-S22*S22-S23*S23)-W13*(S12*S13+S22*S23+S23*S33)-W23*(S11*S13+S12*S23+S13*S33)) + C_6*mu_t*T*T* (S12*Sij) + C_7*mu_t*T*T* (S12*Wij);

			C_UDSI(c,t,UW)= C_1*mu_t*T* (S11*S13+S12*S23+S13*S33) + C_2*mu_t*T*(W13*(S33-S11)+W12*S23-W23*S12) -C_3*mu_t*T*(W12*W23) + C_4*mu_t*T*T* (W13*(S11*S11+S12*S12-S23*S23-S33*S33)-W12*(S11*S13+S22*S23+S23*S33)+W23*(S11*S12+S12*S22 +S13*S23)) + C_6*mu_t*T*T* (S13*Sij) + C_7*mu_t*T*T* (S13*Wij);

			C_UDSI(c,t,VW)= C_1*mu_t*T* (S12*S13+S22*S23+S23*S33) + C_2*mu_t*T*(W23*(S33-S22)-W12*S13-W13*S12) + C_3*mu_t*T*(W12*W13) + C_4*mu_t*T*T* (W23*(S12*S12+S22*S22-S13*S13-S33*S33)+W12*(S11*S13+S12*S23+S13*S33)+W13*(S11*S12+S12*S22+S13*S23)) + C_6*mu_t*T*T* (S23*Sij) + C_7*mu_t*T*T* (S23*Wij);

			C_UDMI(c,t,0)= -2.*mu_t*S11 + 2./3.*C_K(c,t) + C_UDSI(c,t,UU); 
			C_UDMI(c,t,1)= -2.*mu_t*S22 + 2./3.*C_K(c,t) + C_UDSI(c,t,VV); 
			C_UDMI(c,t,2)= -2.*mu_t*S33 + 2./3.*C_K(c,t) + C_UDSI(c,t,WW); 
			C_UDMI(c,t,3)= -2.*mu_t*S12 + C_UDSI(c,t,UV); 
			C_UDMI(c,t,4)= -2.*mu_t*S13 + C_UDSI(c,t,UW); 
			C_UDMI(c,t,5)= -2.*mu_t*S23 + C_UDSI(c,t,VW);

			} 
		end_c_loop(c,t) 
	} 
}

bestucan

先用单线程试试。另外，开fluent是从x64 commander开的么，还是直接双击

benqing

@bestucan 谢谢，这个问题解决了，需要设置UDS的个数，不论有没有标量方程

CFD中文网

UDF，非线性雷诺应力，USDI？