fluent udf DEFINE_ADJUST求两相界面力能load但不能计算

guohf

我用vof法计算气液两相流，想先找出两相界面，并求界面电密度，用DEFINE_ADJUST和DEFINE_ON_DEMAND都试过，能顺利load进去，但一hook这个界面电密度文件fluent就提示致命错误。下面是我的部分程序，请高手指导一下，先谢了。
为方便，附上附件，最后一个udf DEFINE_ADJUST(calm_rhos, domain)出错。

28 秒前 上传
点击文件名下载附件




#include "stdio.h"
#include "udf.h"
#include "sg.h"
#include "sg_mphase.h"
#include "math.h"
#include "flow.h"
#include "mem.h"
#include "metric.h"
#include "id.h"

#define EPS_1 4.25   /* 相对介电常数*/
#define EPS_0 8.8541878e-12  /* electrical permittivity 真空介电常数c^2/Nm^2*/


/* Define which user-defined scalars to use. */
enum                 /*枚举类型名*/
{
phi,    /* 电势*/  
elrho,      /* 电密度*/
N_REQUIRED_UDSI    /*存取用户自定义标量输运方程的数目*/
};

enum                 /*枚举类型名*/
{
Ex,
Ey,
Ez,
vfx,
vfy,
vfz,
rhos,  /* 界面电密度*/
N_REQUIRED_UDM    /*存取用户自定义标量输运方程的数目*/
};
DEFINE_ADJUST(store_VOF_gradient, d)  /*存储vof体积分数的梯度*/
{
Thread *t;
Thread *ppt;
Thread **pt;
cell_t c;
int phase_domain_index=1;
Domain *pDomain = DOMAIN_SUB_DOMAIN(d,phase_domain_index);
Alloc_Storage_Vars(pDomain,SV_VOF_RG,SV_VOF_G,SV_NULL);
Scalar_Reconstruction(pDomain, SV_VOF,-1,SV_VOF_RG,NULL);
Scalar_Derivatives(pDomain,SV_VOF,-1,SV_VOF_G,SV_VOF_RG, Vof_Deriv_Accumulate);
mp_thread_loop_c(t,d,pt)
{
   if (FLUID_THREAD_P(t))
   {
    ppt = pt[phase_domain_index];
    begin_c_loop (c,t)
    {
     C_UDMI(c,t,vfx) = C_VOF_G(c,ppt)[0];
     C_UDMI(c,t,vfy) = C_VOF_G(c,ppt)[1];
     C_UDMI(c,t,vfz) = C_VOF_G(c,ppt)[2];
    }
    end_c_loop (c,t)
   }
}
Free_Storage_Vars(pDomain,SV_VOF_RG,SV_VOF_G,SV_NULL);
}

DEFINE_ADJUST(Calc_E_xyz, domain)  /*存储电场强度*/
{
Thread *t;
cell_t c;
face_t f;
Thread *t0;
cell_t c0;
real dr0[ND_ND],dr1[ND_ND], xf[ND_ND], dy;

if (n_uds < N_REQUIRED_UDSI)
   Internal_Error("not enough user-defined scalars allocated");

if(n_udm < N_REQUIRED_UDM)
   Error("Not enough udm allocated\n");

/* Do nothing if gradient isn't allocated yet. */
if (!Data_Valid_P()) return;
thread_loop_c(t,domain)
{
   if (FLUID_THREAD_P(t))
    if (NULL != THREAD_STORAGE(t,SV_UDS_I(phi)) && NULL != T_STORAGE_R_NV (t,SV_UDSI_G(phi)))
    {
     begin_c_loop(c,t)       /*体循环*/
     {
      C_UDMI(c,t,Ex) = -1.*C_UDSI_G(c,t,phi)[0];    /*Ex*/
      C_UDMI(c,t,Ey) = -1.*C_UDSI_G(c,t,phi)[1];
      C_UDMI(c,t,Ez) = -1.*C_UDSI_G(c,t,phi)[2];
     }
     end_c_loop(c,t)
      /* Message("x electric strength: %g\n", C_UDMI(c,t,Ex)); */
    }
}
}

DEFINE_ADJUST(calm_rhos, domain)   /* 界面电密度,出问题的地方，请帮忙指正*/
{
Thread *t;
Thread *pri_th,*sec_th;
cell_t c;
real xc[ND_ND], EEp[3], EEs[3], nn[3], EE[3], rr;
/*  Domain *domain = Get_Domain(1); Get domain pointer */
int ID = 1;
t =Lookup_Thread(domain, ID); /* mixture-level thread pointer */
pri_th = THREAD_SUB_THREAD(t,0);
sec_th = THREAD_SUB_THREAD(t,1);

thread_loop_c(t,domain)
{
   begin_c_loop(c,t)       /*体循环*/
   {
    C_CENTROID(xc,c,t);
    if (C_VOF(c,sec_th)>0. && C_VOF(c,sec_th)<1.)
    {
     EEs[0]=C_UDMI(c,sec_th,Ex);/* 次相电场强度 */
     EEs[1]=C_UDMI(c,sec_th,Ey);
     EEs[2]=C_UDMI(c,sec_th,Ez);
     EEp[0]=C_UDMI(c,pri_th,Ex);/* 基相电场强度 */
     EEp[1]=C_UDMI(c,pri_th,Ey);
     EEp[2]=C_UDMI(c,pri_th,Ez);
     nn[0]=C_UDMI(c,sec_th,vfx);/* vof梯度、法向矢量 */
     nn[1]=C_UDMI(c,sec_th,vfy);
     nn[2]=C_UDMI(c,sec_th,vfz);
     NV_V(EEs,=,EEp);
     rr=EPS_0*EPS_1;
     NV_VS_VS(EE, =,EEp,*,EPS_0, -, EEs,*,rr);
     C_UDMI(c,t,rhos) = NV_DOT(nn,EE);
    }
    else
    {
     C_UDMI(c,t,rhos) =0;
    }
   }
   end_c_loop(c,t)
}
}

guohf

求指点下了，谢谢 :L

CFD项目工作室

可以发来我们看下  应该可以处理的

guohf

CFD项目工作室 发表于 2017-3-2 09:51
可以发来我们看下  应该可以处理的

非常感谢，我把二维的case和data，UDF发上来，请帮我看下吧，已困扰了好几个月了，谢谢了。