UDF divergence of a vector
 

i am working with the eulerian model and I would like to add a new term in the dispersed phase momentum equation. I can add easy terms using the macro DEFINE_SOURCE but i need the divergence of a vector. I made some easy probes to calculate the divergence of vector but i didn't suceed. Is there any macro to calculate the divergence of a vector? Is there anyone who has made a udf to calculate the divergence? I created some udf that use the divergence theorema to calculate the volume of a cell and the result is good, but when i tried to probe the udf in the steady state continuity equation the udf failed. If theres is someone interested i would send him my udfs. Sorry I am a spanish girl and i explain the best i could.

Did you ever get this? I am working on the same thing. I just submitted a support question to FLUENT, I will post the result.

Casey

Yeah, I did. I cant get the simulation to converge at first. I opened up a service request with fluent, but they said the UDF is running fine and that the problem must be somewhere else. I am not happy with the fact that the simulation does not converge for the first few time steps, so I am going to keep looking into it, but below is the divergence part of the UDF I wrote.
Casey

/* Include statements */

#include "udf.h"

/************************************************** *********************
UDF for specifying frictional pressure
************************************************** **********************/

DEFINE_PROPERTY(pf_cfifteen, cp, tp)
{
real np = 5.;
real nn = 2.;
real esmax = 0.63;
real esmin = 0.5;
real bp = 0.6237;
real sm_val = 1.0e-12;
real P[ND_ND];
real delt, pfr, nc, fric, Smag, vfrac, dudx, dudy, dvdx, dvdy, r, radv, DIV;

/* Retrieve volume fraction and magnitude of the deviatoric Smag tensor */
/*x direction is axial direction, y direction is radial so dvdy = d(radialvelocity)/dradial */

vfrac = C_VOF(cp,tp);
dudx = C_U_G(cp,tp)[0];
dudy = C_U_G(cp,tp)[1];
dvdx = C_V_G(cp,tp)[0];
dvdy = C_V_G(cp,tp)[1];
radv = C_V(cp,tp);
C_CENTROID(P,cp,tp);
r = P[1];
DIV = dvdy + radv/r + dudx;
Smag = sqrt(pow((dvdy-(1./3.)*DIV),2.)+2.*pow((0.5*(dudy-dvdx)),2.)+pow((dudx-(1./3.)*DIV),2.));

__________________________________________________ ____

I had to leave out the rest of the UDF because it pertains to my research and I have not published yet. Below here it just uses the DIV in "if" statements, so any user that wants to use the above DIV part of the code you just need the folowing:
__________________________________________________ __________
dudx = C_U_G(cp,tp)[0];
dudy = C_U_G(cp,tp)[1];
dvdx = C_V_G(cp,tp)[0];
dvdy = C_V_G(cp,tp)[1];
radv = C_V(cp,tp);
C_CENTROID(P,cp,tp);
r = P[1];
DIV = dvdy + radv/r + dudx;
Smag = sqrt(pow((dvdy-(1./3.)*DIV),2.)+2.*pow((0.5*(dudy-dvdx)),2.)+pow((dudx-(1./3.)*DIV),2.));

__________________________________________________ _________

Hope this helps everyone!

Casey

Hi Casey,
if the case is 2D than we should divide the DIV by 2 right?
thanks in advance!