CFD中文网

benqing

我想把图片中的红色框的一部分输出，然后在paraview中做后处理，为了方便，我建立了一个同样的参数（黄色框），想用类似湍流模型里边的方法试一下，但是写不出来
我还尝试了另一种方法：在createfield.h里面定义了一个参数，这样会写出一个文件，但是写出的文件不会更新。
其中m是自定义的一个参数

benqing

@ct 没有，这个不该出现负数的吧

benqing

@ct 后来我直接用的k作为模化的湍动能，用UPrime2Mean的结果计算瞬时的resloved湍动能

benqing

@benqing 解决了，只算一个就可以，算流量比较方便
Q=A×V，A=(Ly-h)×Lz=(3h-h)×4.5h=9h×h，V=1，Q=9×0.028×0.028×1=0.007056，所以我的结果跟@心里的小漩涡的结果不太一样，后续看看计算结果

benqing

@心里的小漩涡小伙伴还在论坛吗？最近有兴趣用fluent跑跑这个case，但是这里面对周期性边界上mass flow rate和pressure Gradient设定不是很理解，这个数值是如何算出来的？

benqing

@李东岳在雷诺应力与URANS 中说：

我一直在找相关的sci，你这面有文章见过类似的操作么

URANS的话，没有见到过，这也只是我自己的想法。我用过PANS模型，有一些文献是将两者相加的，其实也没有理论的介绍，就是后处理的方法。
我找找看有没有URANS相关的文献！

benqing

@李东岳两者相加吧，我认为若URANS可以较好的解析流场的话，uPrime2Mean会占更大比重，但是不能忽略R。理论方面的话，可以参考PANS模型。

benqing

这个问题我解决了，宏C_WALL_DIST在高版本里面确实不能用了，从CFDonline上找到了一个可以用的代码，代码如下：

#include "udf.h"
#include "prox.h"

static cxboolean wall_dist_set = FALSE;//cxboolean布尔运算，与或非，对错

DEFINE_ON_DEMAND(set_wall_dist_udm0)
{
	#if !RP_HOST//表示非host节点的话，执行后面的代码，fluent并行的话，会启动一个 host节点和指定个数的node节点，host节点主要处理GUI和结果收集之类的任务node节点负责分区计算任务

	Domain *domain;
	Thread *t;
	cell_t c;

	if (! wall_dist_set)
	{
		domain = Get_Domain(1);//ROOT_DOMAIN_ID?

		Alloc_Storage_Vars(domain, SV_RTMP_0, SV_NULL);//Alloc_Storage_Vars 用于更改分配

		Calc_Cell_Wall_Distance_New(domain, SV_RTMP_0);

		thread_loop_c(t,domain)
		{
		begin_c_loop(c,t)
			{
				C_UDMI(c,t,0) = C_TMP0(c,t);
			}
		end_c_loop(c,t)
		}
		wall_dist_set = TRUE;
	}
	#endif /* !RP_HOST */
}

DEFINE_ON_DEMAND(reset_udm0)
{
	#if !RP_HOST

	Domain *domain;
	Thread *t;
	cell_t c;

	domain = Get_Domain(1);

	thread_loop_c(t,domain)
	{
	begin_c_loop(c,t)
		{
			C_UDMI(c,t,0) = 0.0;
		}
	end_c_loop(c,t)
	}

	#endif /* !RP_HOST */

	wall_dist_set = FALSE;
}

说实话这个代码我不能完全看懂，但是主要思想就是用DEFINE_ON_DEMAND定义两个命令reset_udm0和set_wall_dist_udm0，前者将壁面到最近壁面的距离reset为0，后者可以求出想要的结果，经过测试两个命令应该都是必须的，另外我放的代码和CFDonline上有一点不同，domain = Get_Domain(1);这一行在原代码里面是domain = Get_Domain(ROOT_DOMAIN_ID);注意修改，下图是我测试的一个水翼网格的结果，结果是对的，

再啰嗦几句代码使用方法，这个代码直接复制到一个.c文件中，在fluent中添加memory的数量为2，编译.c文件并load，在calculate之前，点Execute on Demand，并先后运行reset_udm0和set_wall_dist_udm0，然后就可以calculate了

benqing

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

benqing

你们的问题解决了吗，我遇到了同样的问题

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) 
	} 
}

benqing

@yfleung 这个具体没写过，找找定义模仿现有程序试试吧

benqing

@yfleung postChannel这个程序很好改，你可以自定义添加量，比如k和epsilon等等；求utau的话，需要壁面剪切应力(tau_wall)，u_tau=sqrt(tau_wall/密度)--其实在of中不需要除以密度，计算壁面剪切应力用wallShearStress，但是看代码的话，wallShearStress采用的速度是瞬时速度，其实这个工具适合计算定常，如果非定常的话，我修改过程序采用的UMean；这个时候得到的是面上的数据，你自己再做合适的后处理就行了

benqing

@yfleung 对于槽道流有一个非常好用的工具，叫postChannel，在channel395那个算例能用，这个工具可以做空间(槽道流中应该是流向和展向)的平均，比较容易用，具体你可以自己下看下，如果有问题可以再问。

benqing

@fangyuanaza 没有遇到过，你去网上找找原因吧，可以看看是不是版本原因，或者这个字典里其他的位置的错误，我用的2.3.0，感觉还是设定的问题

benqing

@fangyuanaza 我当时用face这个类型，你试试吧

benqing

@Mania 还有另外一种方法是安装两个版本的openfoam，然后后处理的时候用低版本的操作

benqing

@Mania 抱歉昨天忙了一天没有看到消息，没回复。
他的版本是比较低的版本，大概是2.3.0左右的，所以高版本下会出问题，你可以找下那个连接下有没有高版本的
或者另外一个帖子下有其他的链接，你可以试一下

benqing

@Mania 我的意思是直接在第一次计算的时候就用大括号，不知道会不会有影响。
有钝体的话，matlab作图确实会出现自动插值的图，可能会有问题。当然也是可以解决的，就是在有钝体的地方你要放上网格信息，并把值赋值，做成白色就可以。
你是想做展向平均吗，做展向平均的话，cfdonline上有别人写的程序

benqing

@Mania 不知道版本有没有影响，在2017版本里面公式内的变量需要用大括号{}，你可以看下tecplot help，下面这张图是17版本的help的提示

所以你换这个写的话，应该是{V34}=...，你可以试一下。
另外，其实你这个做三个面的平均的话，完全可以用matlab，非常简单，几行代码的事