CFD中文网

冠

各位大佬，我尝试在reactingTwoPhaseEulerFoam中添加空化模型，参考的是
www.tfd.chalmers.se/~hani/kurser/OS_CFD_2017/SuryaKaundinyaOruganti/Final_presentation_Surya_Oruganti.pdf
中的方式，但是在修改phaseChangeModel.C文件时产生了报错。
报错如下
0979155c-344d-4cdb-a9d3-449726064a30-图片.png
这部分我的代码与参考资料中给的一致：

namespace Foam
{
    defineTypeNameAndDebug(phaseChangeModel, 0);
    defineRunTimeSelectionTable(phaseChangeModel, components);
}

// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //

Foam::phaseChangeModel::phaseChangeModel
(
    const word& type,
    const dictionary& interfaceDict,
    const phaseModel& phase1,
    const phaseModel& phase2

)
:
    interfaceDict_(interfaceDict),
    phase1_(phase1),
    phase2_(phase2),
    phaseChangeModelCoeffs_(interfaceDict.optionalSubDict(type + "Coeffs")),
    pSat_("pSat", dimPressure, interfaceDict.lookup("pSat"))

{}


// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //

Foam::Pair<Foam::tmp<Foam::volScalarField>>
Foam::phaseChangeModel::vDotAlphal() const
{

    volScalarField alphalCoeff(1.0/phase1_.rho() - phase1_*(1.0/phase1_.rho() - 1.0/phase2_.rho()));
    Pair<tmp<volScalarField>> mDotAlphal = this->mDotAlphal();

    return Pair<tmp<volScalarField>>
    (
        alphalCoeff*mDotAlphal[0],
        alphalCoeff*mDotAlphal[1]
    );
}

Foam::Pair<Foam::tmp<Foam::volScalarField>>
Foam::phaseChangeModel::vDotP() const
{
    dimensionedScalar pCoeff(1.0/phase1_.rho() - 1.0/phase2_.rho());
    Pair<tmp<volScalarField>> mDotP = this->mDotP();

    return Pair<tmp<volScalarField>>(pCoeff*mDotP[0], pCoeff*mDotP[1]);
}

这种问题是什么原因导致的呢？我刚开始接触这种模型编译，不太懂

冠

各位大佬，最近想弄个OpenFOAM的冷凝模拟，想问下OpenFOAM自带的冷凝模型都有啥，在哪找呢

冠

我想询问一下，fluent中的mixture模型，在OpenFOAM中有相对应模型吗，在哪个位置呢？

冠

@李东岳谢谢老师的建议，我担心学习of的时间不太够

冠

@bestucan 谢谢老师，我还想问一下我uds方程里对流项没有密度，我直接把udf帮助文档里面对流项的代码中密度的地方改为1可以吧？
6dd6d91a-8e8a-46f6-b499-d4ce49fde849-68216d0572e0cd03f7d2bae473bb343.png

DEFINE_UDS_FLUX(No_Cxx_FLUX,f,t,i)                 //Cxx对流项    
{
    	cell_t c0, c1 = -1;
    	Thread *t0, *t1 = NULL;
    	real NV_VEC(psi_vec), NV_VEC(A), flux = 0.0;
    	c0 = F_C0(f,t);
    	t0 = F_C0_THREAD(f,t);
    	F_AREA(A, f, t);
    
   	if (BOUNDARY_FACE_THREAD_P(t)) 
    	{
       		real dens;
       
       		if (NNULLP(THREAD_STORAGE(t,SV_DENSITY)))
         		dens = F_R(f,t);  
       		else
         		dens = C_R(c0,t0); 
       		NV_DS(psi_vec, =, F_U(f,t), F_V(f,t), F_W(f,t), *, 1);
      		flux = NV_DOT(psi_vec, A); 
	}
    	else
    	{
     		c1 = F_C1(f,t);   
     		t1 = F_C1_THREAD(f,t);
     		NV_DS(psi_vec, =, C_U(c0,t0),C_V(c0,t0),C_W(c0,t0),*,1);
     		NV_DS(psi_vec, +=, C_U(c1,t1),C_V(c1,t1),C_W(c1,t1),*,1);
     		flux = NV_DOT(psi_vec, A)/2.0;
    	}
	return flux;
}

冠

不会是因为这个警告吧？
502eeef3-d4b5-40d8-bb39-e9fa15299eba-dc8a5ccdc03793d837070a5fdbcebea.png

冠

我想在fluent中使用黏弹性流体(FENE-P模型)，需要使用udf在动量方程中添加弹性项源项，并且使用uds添加6个标量方程（FENE-P模型的本构方程Cxx,Cxy,Cxz,Cyy,Cyz,Czz），与N-S方程耦合求解，编译成功后可以计算，但是uds的值变化很小（初始化时uds的值给的是1，计算结束后基本还是1，知识小数点后几位有变化），几乎体现不出弹性项的作用。求助各位大佬
（1）会不会是多个uds方程需要对i进行定义
在帮助手册找到这个（如果要求解多个标量，则可以在UDF中使用条件IF语句为每个i定义不同的通量函数。i=0与uds-0(要求解的第一个标量方程)相关联。）但是没找到具体怎么操作
FENE-P 模型的本构方程是
1db2aecf-bb28-4e18-83b2-b78012d300a0-8de57d7561f6f0ce3fa3533604fefd3.png
直角坐标下为
2c357f13-0b20-4832-8817-2d63da569bd0-4ad7186b8fd9973a10024ee26d56f59.png
49efe7aa-242a-40a1-99f0-77492f023026-e3f29d15bfaa922b3ad458e6135be3c.png
2f7eef5f-0dc2-481c-ae2e-fbf59aaf90f4-2675ed868f2111033ffc1d7af482a07.png
添加弹性项后的X方向上的动量方程
7da0ce9c-10e5-467d-bbc7-a5e298f58678-8985c92d1edb7e6c3b4912f0ccb09ae.png
添加弹性项后的y方向上的动量方程
3c27a4d7-a3de-4795-97bf-85f89b1ca171-7f1fc03c84026987030e60dd8546d30.png
添加弹性项后的z方向上的动量方程
3ba2dfa4-8db9-4ae7-ac16-e72c39f52553-4506fbbd7bf7138e779729d1836a718.png
具体的udf代码

#include "mem.h"
#define miup (0.00736)
#define lamda (2.47024)
#define L (40)


DEFINE_UDS_UNSTEADY(No_Cxx_unsteady,c,t,i,apu,su)        //cxx非稳态项
 {
 real physical_dt, vol, rho, phi_old;
 physical_dt = RP_Get_Real("physical-time-step");
 vol = C_VOLUME(c,t);
 
 *apu = -vol / physical_dt;/*implicit part*/
 phi_old = C_STORAGE_R(c,t,SV_UDSI_M1(0));
 *su = vol*phi_old/physical_dt;/*explicit part*/
 }


DEFINE_UDS_UNSTEADY(No_Cxy_unsteady,c,t,i,apu,su)         //cxy非稳态项
 {
 real physical_dt, vol, rho, phi_old;
 physical_dt = RP_Get_Real("physical-time-step");
 vol = C_VOLUME(c,t);
 
 *apu = -vol / physical_dt;/*implicit part*/
 phi_old = C_STORAGE_R(c,t,SV_UDSI_M1(1));
 *su = vol*phi_old/physical_dt;/*explicit part*/
 }

 DEFINE_UDS_UNSTEADY(No_Cxz_unsteady,c,t,i,apu,su)       //cxz非稳态项
 {
 real physical_dt, vol, rho, phi_old;
 physical_dt = RP_Get_Real("physical-time-step");
 vol = C_VOLUME(c,t);
 
 *apu = -vol / physical_dt;/*implicit part*/
 phi_old = C_STORAGE_R(c,t,SV_UDSI_M1(2));
 *su = vol*phi_old/physical_dt;/*explicit part*/
 }

 

 DEFINE_UDS_UNSTEADY(No_Cyy_unsteady,c,t,i,apu,su)       //cyy非稳态项
 {
 real physical_dt, vol, rho, phi_old;
 physical_dt = RP_Get_Real("physical-time-step");
 vol = C_VOLUME(c,t);
 
 *apu = -vol / physical_dt;/*implicit part*/
 phi_old = C_STORAGE_R(c,t,SV_UDSI_M1(3));
 *su = vol*phi_old/physical_dt;/*explicit part*/
 }

 DEFINE_UDS_UNSTEADY(No_Cyz_unsteady,c,t,i,apu,su)       //cyz非稳态项
 {
 real physical_dt, vol, rho, phi_old;
 physical_dt = RP_Get_Real("physical-time-step");
 vol = C_VOLUME(c,t);
 
 *apu = -vol / physical_dt;/*implicit part*/
 phi_old = C_STORAGE_R(c,t,SV_UDSI_M1(4));
 *su = vol*phi_old/physical_dt;/*explicit part*/
 }

 

 DEFINE_UDS_UNSTEADY(No_Czz_unsteady,c,t,i,apu,su)       //czz非稳态项
 {
     	real physical_dt, vol, rho, phi_old;
     	physical_dt = RP_Get_Real("physical-time-step");
     	vol = C_VOLUME(c,t);
 
     	*apu = -vol / physical_dt;/*implicit part*/
     	phi_old = C_STORAGE_R(c,t,SV_UDSI_M1(5));
     	*su = vol*phi_old/physical_dt;/*explicit part*/
 }

DEFINE_UDS_FLUX(No_Cxx_FLUX,f,t,i)                 //Cxx对流项    
{
    	cell_t c0, c1 = -1;
    	Thread *t0, *t1 = NULL;
    	real NV_VEC(psi_vec), NV_VEC(A), flux = 0.0;
    	c0 = F_C0(f,t);
    	t0 = F_C0_THREAD(f,t);
    	F_AREA(A, f, t);
    
   	if (BOUNDARY_FACE_THREAD_P(t)) 
    	{
       		real dens;
       
       		if (NNULLP(THREAD_STORAGE(t,SV_DENSITY)))
         		dens = F_R(f,t);  
       		else
         		dens = C_R(c0,t0); 
       		NV_DS(psi_vec, =, F_U(f,t), F_V(f,t), F_W(f,t), *, 1);
      		flux = NV_DOT(psi_vec, A); 
	}
    	else
    	{
     		c1 = F_C1(f,t);   
     		t1 = F_C1_THREAD(f,t);
     		NV_DS(psi_vec, =, C_U(c0,t0),C_V(c0,t0),C_W(c0,t0),*,1);
     		NV_DS(psi_vec, +=, C_U(c1,t1),C_V(c1,t1),C_W(c1,t1),*,1);
     		flux = NV_DOT(psi_vec, A)/2.0;
    	}
	return flux;
} 



DEFINE_UDS_FLUX(No_Cxy_FLUX,f,t,i)                 //Cxy对流项    
{
    	cell_t c0, c1 = -1;
    	Thread *t0, *t1 = NULL;
    	real NV_VEC(psi_vec), NV_VEC(A), flux = 0.0;
    	c0 = F_C0(f,t);
    	t0 = F_C0_THREAD(f,t);
    	F_AREA(A, f, t);
    
   	if (BOUNDARY_FACE_THREAD_P(t)) 
    	{
       		real dens;
       
       		if (NNULLP(THREAD_STORAGE(t,SV_DENSITY)))
         		dens = F_R(f,t);  
       		else
         		dens = C_R(c0,t0); 
       		NV_DS(psi_vec, =, F_U(f,t), F_V(f,t), F_W(f,t), *, 1);
      		flux = NV_DOT(psi_vec, A); 
	}
    	else
    	{
     		c1 = F_C1(f,t);   
     		t1 = F_C1_THREAD(f,t);
     		NV_DS(psi_vec, =, C_U(c0,t0),C_V(c0,t0),C_W(c0,t0),*,1);
     		NV_DS(psi_vec, +=, C_U(c1,t1),C_V(c1,t1),C_W(c1,t1),*,1);
     		flux = NV_DOT(psi_vec, A)/2.0;
    	}
	return flux;
} 


DEFINE_UDS_FLUX(No_Cxz_FLUX,f,t,i)                 //Cxz对流项    
{
    	cell_t c0, c1 = -1;
    	Thread *t0, *t1 = NULL;
    	real NV_VEC(psi_vec), NV_VEC(A), flux = 0.0;
    	c0 = F_C0(f,t);
    	t0 = F_C0_THREAD(f,t);
    	F_AREA(A, f, t);
    
   	if (BOUNDARY_FACE_THREAD_P(t)) 
    	{
       		real dens;
       
       		if (NNULLP(THREAD_STORAGE(t,SV_DENSITY)))
         		dens = F_R(f,t);  
       		else
         		dens = C_R(c0,t0); 
       		NV_DS(psi_vec, =, F_U(f,t), F_V(f,t), F_W(f,t), *, 1);
      		flux = NV_DOT(psi_vec, A); 
	}
    	else
    	{
     		c1 = F_C1(f,t);   
     		t1 = F_C1_THREAD(f,t);
     		NV_DS(psi_vec, =, C_U(c0,t0),C_V(c0,t0),C_W(c0,t0),*,1);
     		NV_DS(psi_vec, +=, C_U(c1,t1),C_V(c1,t1),C_W(c1,t1),*,1);
     		flux = NV_DOT(psi_vec, A)/2.0;
    	}
	return flux;
} 



DEFINE_UDS_FLUX(No_Cyy_FLUX,f,t,i)                 //Cyy对流项    
{
    	cell_t c0, c1 = -1;
    	Thread *t0, *t1 = NULL;
    	real NV_VEC(psi_vec), NV_VEC(A), flux = 0.0;
    	c0 = F_C0(f,t);
    	t0 = F_C0_THREAD(f,t);
    	F_AREA(A, f, t);
    
   	if (BOUNDARY_FACE_THREAD_P(t)) 
    	{
       		real dens;
       
       		if (NNULLP(THREAD_STORAGE(t,SV_DENSITY)))
         		dens = F_R(f,t);  
       		else
         		dens = C_R(c0,t0); 
       		NV_DS(psi_vec, =, F_U(f,t), F_V(f,t), F_W(f,t), *, 1);
      		flux = NV_DOT(psi_vec, A); 
	}
    	else
    	{
     		c1 = F_C1(f,t);   
     		t1 = F_C1_THREAD(f,t);
     		NV_DS(psi_vec, =, C_U(c0,t0),C_V(c0,t0),C_W(c0,t0),*,1);
     		NV_DS(psi_vec, +=, C_U(c1,t1),C_V(c1,t1),C_W(c1,t1),*,1);
     		flux = NV_DOT(psi_vec, A)/2.0;
    	}
	return flux;
} 

DEFINE_UDS_FLUX(No_Cyz_FLUX,f,t,i)                 //Cyz对流项    
{
    	cell_t c0, c1 = -1;
    	Thread *t0, *t1 = NULL;
    	real NV_VEC(psi_vec), NV_VEC(A), flux = 0.0;
    	c0 = F_C0(f,t);
    	t0 = F_C0_THREAD(f,t);
    	F_AREA(A, f, t);
    
   	if (BOUNDARY_FACE_THREAD_P(t)) 
    	{
       		real dens;
       
       		if (NNULLP(THREAD_STORAGE(t,SV_DENSITY)))
         		dens = F_R(f,t);  
       		else
         		dens = C_R(c0,t0); 
       		NV_DS(psi_vec, =, F_U(f,t), F_V(f,t), F_W(f,t), *, 1);
      		flux = NV_DOT(psi_vec, A); 
	}
    	else
    	{
     		c1 = F_C1(f,t);   
     		t1 = F_C1_THREAD(f,t);
     		NV_DS(psi_vec, =, C_U(c0,t0),C_V(c0,t0),C_W(c0,t0),*,1);
     		NV_DS(psi_vec, +=, C_U(c1,t1),C_V(c1,t1),C_W(c1,t1),*,1);
     		flux = NV_DOT(psi_vec, A)/2.0;
    	}
	return flux;
} 



DEFINE_UDS_FLUX(No_Czz_FLUX,f,t,i)                 //Czz对流项    
{
    	cell_t c0, c1 = -1;
    	Thread *t0, *t1 = NULL;
    	real NV_VEC(psi_vec), NV_VEC(A), flux = 0.0;
    	c0 = F_C0(f,t);
    	t0 = F_C0_THREAD(f,t);
    	F_AREA(A, f, t);
    
   	if (BOUNDARY_FACE_THREAD_P(t)) 
    	{
       		real dens;
       
       		if (NNULLP(THREAD_STORAGE(t,SV_DENSITY)))
         		dens = F_R(f,t);  
       		else
         		dens = C_R(c0,t0); 
       		NV_DS(psi_vec, =, F_U(f,t), F_V(f,t), F_W(f,t), *, 1);
      		flux = NV_DOT(psi_vec, A); 
	}
    	else
    	{
     		c1 = F_C1(f,t);   
     		t1 = F_C1_THREAD(f,t);
     		NV_DS(psi_vec, =, C_U(c0,t0),C_V(c0,t0),C_W(c0,t0),*,1);
     		NV_DS(psi_vec, +=, C_U(c1,t1),C_V(c1,t1),C_W(c1,t1),*,1);
     		flux = NV_DOT(psi_vec, A)/2.0;
    	}
	return flux;
} 

DEFINE_SOURCE(No_Cxx_source, c, t, dS, eqn)             //Cxx源项
{  
	real source;
    real fR;
    fR=(L*L-3)/(L*L-C_UDSI(c,t,0)-C_UDSI(c,t,3)-C_UDSI(c,t,5));
	source=2*(C_UDSI(c,t,0)*C_DUDX(c,t)+C_UDSI(c,t,1)*C_DUDY(c,t)+C_UDSI(c,t,2)*C_DUDZ(c,t))-(fR*C_UDSI(c,t,0)-1)/lamda;
	dS[eqn]=0.0;

	return source;
}

DEFINE_SOURCE(No_Cxy_source, c, t, dS, eqn)             //Cxy源项
{  
	real source;
    real fR;
    fR=(L*L-3)/(L*L-C_UDSI(c,t,0)-C_UDSI(c,t,3)-C_UDSI(c,t,5));
	source=(C_UDSI(c,t,0)*C_DVDX(c,t)+C_UDSI(c,t,1)*C_DVDY(c,t)+C_UDSI(c,t,2)*C_DVDZ(c,t)+\
    C_UDSI(c,t,1)*C_DUDX(c,t)+C_UDSI(c,t,3)*C_DUDY(c,t)+C_UDSI(c,t,4)*C_DUDZ(c,t))-\
    (fR*C_UDSI(c,t,1))/lamda;
	dS[eqn]=0.0;

	return source;
}


DEFINE_SOURCE(No_Cxz_source, c, t, dS, eqn)             //Cxz源项
{  
	real source;
    real fR;
    fR=(L*L-3)/(L*L-C_UDSI(c,t,0)-C_UDSI(c,t,3)-C_UDSI(c,t,5));
	source=(C_UDSI(c,t,0)*C_DWDX(c,t)+C_UDSI(c,t,1)*C_DWDY(c,t)+C_UDSI(c,t,2)*C_DWDZ(c,t)+\
    C_UDSI(c,t,2)*C_DUDX(c,t)+C_UDSI(c,t,4)*C_DUDY(c,t)+C_UDSI(c,t,5)*C_DUDZ(c,t))-\
    (fR*C_UDSI(c,t,2))/lamda;
	dS[eqn]=0.0;

	return source;
}




DEFINE_SOURCE(No_Cyy_source, c, t, dS, eqn)             //Cyy源项
{  
	real source;
    real fR;
    fR=(L*L-3)/(L*L-C_UDSI(c,t,0)-C_UDSI(c,t,3)-C_UDSI(c,t,5));
	source=2*(C_UDSI(c,t,1)*C_DVDX(c,t)+C_UDSI(c,t,3)*C_DVDY(c,t)+C_UDSI(c,t,4)*C_DVDZ(c,t))-(fR*C_UDSI(c,t,3)-1)/lamda;
	dS[eqn]=0.0;

	return source;
}


DEFINE_SOURCE(No_Cyz_source, c, t, dS, eqn)             //Cyz源项
{  
	real source;
    real fR;
    fR=(L*L-3)/(L*L-C_UDSI(c,t,0)-C_UDSI(c,t,3)-C_UDSI(c,t,5));
	source=(C_UDSI(c,t,1)*C_DWDX(c,t)+C_UDSI(c,t,3)*C_DWDY(c,t)+C_UDSI(c,t,4)*C_DWDZ(c,t)+\
    C_UDSI(c,t,2)*C_DVDX(c,t)+C_UDSI(c,t,4)*C_DVDY(c,t)+C_UDSI(c,t,5)*C_DVDZ(c,t))-\
    (fR*C_UDSI(c,t,4))/lamda;
	dS[eqn]=0.0;

	return source;
}







DEFINE_SOURCE(No_Czz_source, c, t, dS, eqn)             //Czz源项
{  
	real source;
    real fR;
    fR=(L*L-3)/(L*L-C_UDSI(c,t,0)-C_UDSI(c,t,3)-C_UDSI(c,t,5));
	source=2*(C_UDSI(c,t,2)*C_DWDX(c,t)+C_UDSI(c,t,4)*C_DWDY(c,t)+C_UDSI(c,t,5)*C_DWDZ(c,t))-(fR*C_UDSI(c,t,5)-1)/lamda;
	dS[eqn]=0.0;

	return source;
}


DEFINE_SOURCE(x_momentum, c, t, dS, eqn)             //x方向动量源项
{  
	real source;
    real cxx = C_UDSI(c,t,0);
    real cxy = C_UDSI(c,t,1);
    real cxz = C_UDSI(c,t,2);
    real cyx = C_UDSI(c,t,1);
    real cyy = C_UDSI(c,t,3);
    real cyz = C_UDSI(c,t,4);
    real czx = C_UDSI(c,t,2);
    real czy = C_UDSI(c,t,4);
    real czz = C_UDSI(c,t,5);
    real NV_VEC(gcxx);
    real NV_VEC(gcxy);
    real NV_VEC(gcxz);
    real NV_VEC(gcyx);
    real NV_VEC(gcyy);
    real NV_VEC(gcyz);
    real NV_VEC(gczx);
    real NV_VEC(gczy);
    real NV_VEC(gczz);
    NV_DS(gcxx, =,C_UDSI_G(c,t,0)[0],C_UDSI_G(c,t,0)[1],C_UDSI_G(c,t,0)[2],*,1);
    NV_DS(gcxy, =,C_UDSI_G(c,t,1)[0],C_UDSI_G(c,t,1)[1],C_UDSI_G(c,t,1)[2],*,1);
    NV_DS(gcxz, =,C_UDSI_G(c,t,2)[0],C_UDSI_G(c,t,2)[1],C_UDSI_G(c,t,2)[2],*,1);
    NV_DS(gcyx, =,C_UDSI_G(c,t,1)[0],C_UDSI_G(c,t,1)[1],C_UDSI_G(c,t,1)[2],*,1);
    NV_DS(gcyy, =,C_UDSI_G(c,t,3)[0],C_UDSI_G(c,t,3)[1],C_UDSI_G(c,t,3)[2],*,1);
    NV_DS(gcyz, =,C_UDSI_G(c,t,4)[0],C_UDSI_G(c,t,4)[1],C_UDSI_G(c,t,4)[2],*,1);
    NV_DS(gczx, =,C_UDSI_G(c,t,2)[0],C_UDSI_G(c,t,2)[1],C_UDSI_G(c,t,2)[2],*,1);
    NV_DS(gczy, =,C_UDSI_G(c,t,4)[0],C_UDSI_G(c,t,4)[1],C_UDSI_G(c,t,4)[2],*,1);
    NV_DS(gczz, =,C_UDSI_G(c,t,5)[0],C_UDSI_G(c,t,5)[1],C_UDSI_G(c,t,5)[2],*,1);
    source = miup*(L*L-3)/lamda/pow(L*L-cxx-cyy-czz,2)*((gcxx[0]+gcxy[1]+gcxz[2])*(L*L-cxx-cyy-czz)+\
    cxx*(gcxx[0]+gcyy[0]+gczz[0])+cxy*(gcxx[1]+gcyy[1]+gczz[1])+cxz*(gcxx[2]+gcyy[2]+gczz[2]));
	return source;
}


DEFINE_SOURCE(y_momentum, c, t, dS, eqn)             //y方向动量源项
{  
	real source;
    real cxx = C_UDSI(c,t,0);
    real cxy = C_UDSI(c,t,1);
    real cxz = C_UDSI(c,t,2);
    real cyx = C_UDSI(c,t,1);
    real cyy = C_UDSI(c,t,3);
    real cyz = C_UDSI(c,t,4);
    real czx = C_UDSI(c,t,2);
    real czy = C_UDSI(c,t,4);
    real czz = C_UDSI(c,t,5);
    real NV_VEC(gcxx);
    real NV_VEC(gcxy);
    real NV_VEC(gcxz);
    real NV_VEC(gcyx);
    real NV_VEC(gcyy);
    real NV_VEC(gcyz);
    real NV_VEC(gczx);
    real NV_VEC(gczy);
    real NV_VEC(gczz);
    NV_DS(gcxx, =,C_UDSI_G(c,t,0)[0],C_UDSI_G(c,t,0)[1],C_UDSI_G(c,t,0)[2],*,1);
    NV_DS(gcxy, =,C_UDSI_G(c,t,1)[0],C_UDSI_G(c,t,1)[1],C_UDSI_G(c,t,1)[2],*,1);
    NV_DS(gcxz, =,C_UDSI_G(c,t,2)[0],C_UDSI_G(c,t,2)[1],C_UDSI_G(c,t,2)[2],*,1);
    NV_DS(gcyx, =,C_UDSI_G(c,t,1)[0],C_UDSI_G(c,t,1)[1],C_UDSI_G(c,t,1)[2],*,1);
    NV_DS(gcyy, =,C_UDSI_G(c,t,3)[0],C_UDSI_G(c,t,3)[1],C_UDSI_G(c,t,3)[2],*,1);
    NV_DS(gcyz, =,C_UDSI_G(c,t,4)[0],C_UDSI_G(c,t,4)[1],C_UDSI_G(c,t,4)[2],*,1);
    NV_DS(gczx, =,C_UDSI_G(c,t,2)[0],C_UDSI_G(c,t,2)[1],C_UDSI_G(c,t,2)[2],*,1);
    NV_DS(gczy, =,C_UDSI_G(c,t,4)[0],C_UDSI_G(c,t,4)[1],C_UDSI_G(c,t,4)[2],*,1);
    NV_DS(gczz, =,C_UDSI_G(c,t,5)[0],C_UDSI_G(c,t,5)[1],C_UDSI_G(c,t,5)[2],*,1);
    source = miup*(L*L-3)/lamda/pow(L*L-cxx-cyy-czz,2)*((gcyx[0]+gcyy[1]+gcyz[2])*(L*L-cxx-cyy-czz)+\
    cyx*(gcxx[0]+gcyy[0]+gczz[0])+cyy*(gcxx[1]+gcyy[1]+gczz[1])+cyz*(gcxx[2]+gcyy[2]+gczz[2]));
	return source;
}

DEFINE_SOURCE(z_momentum, c, t, dS, eqn)             //z方向动量源项
{  
	real source;
    real cxx = C_UDSI(c,t,0);
    real cxy = C_UDSI(c,t,1);
    real cxz = C_UDSI(c,t,2);
    real cyx = C_UDSI(c,t,1);
    real cyy = C_UDSI(c,t,3);
    real cyz = C_UDSI(c,t,4);
    real czx = C_UDSI(c,t,2);
    real czy = C_UDSI(c,t,4);
    real czz = C_UDSI(c,t,5);
    real NV_VEC(gcxx);
    real NV_VEC(gcxy);
    real NV_VEC(gcxz);
    real NV_VEC(gcyx);
    real NV_VEC(gcyy);
    real NV_VEC(gcyz);
    real NV_VEC(gczx);
    real NV_VEC(gczy);
    real NV_VEC(gczz);
    NV_DS(gcxx, =,C_UDSI_G(c,t,0)[0],C_UDSI_G(c,t,0)[1],C_UDSI_G(c,t,0)[2],*,1);
    NV_DS(gcxy, =,C_UDSI_G(c,t,1)[0],C_UDSI_G(c,t,1)[1],C_UDSI_G(c,t,1)[2],*,1);
    NV_DS(gcxz, =,C_UDSI_G(c,t,2)[0],C_UDSI_G(c,t,2)[1],C_UDSI_G(c,t,2)[2],*,1);
    NV_DS(gcyx, =,C_UDSI_G(c,t,1)[0],C_UDSI_G(c,t,1)[1],C_UDSI_G(c,t,1)[2],*,1);
    NV_DS(gcyy, =,C_UDSI_G(c,t,3)[0],C_UDSI_G(c,t,3)[1],C_UDSI_G(c,t,3)[2],*,1);
    NV_DS(gcyz, =,C_UDSI_G(c,t,4)[0],C_UDSI_G(c,t,4)[1],C_UDSI_G(c,t,4)[2],*,1);
    NV_DS(gczx, =,C_UDSI_G(c,t,2)[0],C_UDSI_G(c,t,2)[1],C_UDSI_G(c,t,2)[2],*,1);
    NV_DS(gczy, =,C_UDSI_G(c,t,4)[0],C_UDSI_G(c,t,4)[1],C_UDSI_G(c,t,4)[2],*,1);
    NV_DS(gczz, =,C_UDSI_G(c,t,5)[0],C_UDSI_G(c,t,5)[1],C_UDSI_G(c,t,5)[2],*,1);
    source = miup*(L*L-3)/lamda/pow(L*L-cxx-cyy-czz,2)*((gczx[0]+gczy[1]+gczz[2])*(L*L-cxx-cyy-czz)+\
    czx*(gcxx[0]+gcyy[0]+gczz[0])+czy*(gcxx[1]+gcyy[1]+gczz[1])+czz*(gcxx[2]+gcyy[2]+gczz[2]));
	return source;
}
code_text

冠

各位大佬，我想在OpenFOAMv2206中植入我在OpenFOAM9中修改的模型。但是在编译过程中出现了错误：

9cb8d1c3-cfce-407c-a09c-cdecbf702dd3-图片.png

我的代码在OF9中是可以运行的，但是在OpenFOAMv2206中就报错了。看其他的代码感觉OpenFOAM9和OpenFOAMv2206没有太大的区别，这个错误是为什么呢？

Foam::autoPtr<Foam::viscoelasticLogLaw> Foam::viscoelasticLogLaw::New
(
    const word& name,
    const volVectorField& U,
    const surfaceScalarField& phi,
    const dictionary& dict
)
{
    const word typeName(dict.get<word>("type"));

    Info<< "Selecting viscoelasticLog model " << typeName << endl;

    dictionaryConstructorTable::iterator cstrIter =
        dictionaryConstructorTablePtr_->find(typeName);

    if (cstrIter == dictionaryConstructorTablePtr_->end())
    {
        FatalErrorIn
        (
            "viscoelasticLogLaw::New(const word& name, const volVectorField&, "
            "const surfaceScalarField&)"
        )   << "Unknown viscoelasticLogLaw type " << typeName
            << endl << endl
            << "Valid viscoelasticLogLaw types are :" << endl
            << dictionaryConstructorTablePtr_->sortedToc()
            << exit(FatalError);
    }

    return autoPtr<viscoelasticLogLaw>(cstrIter()(name, U, phi, dict));
}

冠

@vbcwl 老哥我这么理解对不对，就是速度设置成U+，也就是U/Ut，那这个按照自带的算例Ret=395的话，是不是就是U=Ubar/ut=0.1335/0.079呢

冠

各位大佬，我想复现一下Ret=395的一个槽道流，但是有一个疑惑的点是平均速度应该如何设置呢？Ret和传统的Re之间有换算关系么

冠

各位大佬，在DNS槽道流case中想求取yPlus进行数据分析，但OpenFOAM里边的yPlus工具只适用于LES和Rans，有适用于DNS的yPlus工具吗

冠

@guangzhou 这个不影响吧，我用rheoTool的目的是他的粘弹性流体库。他的求解器对我来说限制有点多，所以我就直接把他当库文件用自己新写了一个DNS求解器。如果觉得他的不适合的话就引用库自己改一个就好了啊

冠

@guangzhou 这个求解器在高维森贝格数下会产生稳定性问题，所以一般采用对数重构法去解决这个问题。论文可以参考张红娜老师在17年的一篇论文《Numerical Simulation of Heat Transfer Process of Viscoelastic Fluid Flow at High Weissenberg Number by Log-Conformation Reformulation》，或者可以看看rheoTool工具箱相关介绍

冠

破案了，我用OF7就不太行，改成OF8就可以了，可能是版本问题或者是我安装OF7的时候有啥问题导致的，改个版本就好了

冠

@李东岳东岳老师好，这个扰动已经添加了，就是他一开始自带的，我也用perturbU添加了扰动，这两个都是在生成湍流后又变成了这个样子的

冠

各位大佬，我在icoFoam中添加了fvOptions源项，


        fvVectorMatrix UEqn
        (
            fvm::ddt(U)
          + fvm::div(phi, U)
          - fvm::laplacian(nu, U)
           ==
            fvOptions(U)
        );
UEqn.relax();

之后使用meanVelocityForce，测试了channel395，但发现云图最终跑成了这个样子：
654f20ee-a065-4c3a-8673-17ffc2069392-图片.png
想问一下各位大佬，这个是我哪里做的有问题呢？

冠

@hhh 解决了，需要现在ICEM中设置周期性边界，之后就可以了
我参照的是https://blog.csdn.net/muyangzixue/article/details/108631943

冠

@冠竹就是我是否可以通过一个简单的开方，得到脉动速度u'v'w'呢？

冠

@coolhhh 大佬，我还有个问题，这个fieldAvergae中的prime2Mean求取的是雷诺应力，而雷诺应力定义为流体微元表面上脉动动量输运的平均值，那我应该怎样去用雷诺应力去计算脉动速度呢

冠

@coolhhh 感谢大佬

CFD中文网

冠竹

帖子