NEW surface tension model for Fluent
 

Hello there!
I have made a new surface tension model. It gives no spurious currents within the volume of secondary phase, but there is always spurious velocity in one or two cells on interface. This velocity is increasing with time. What is it caused by?
I think I should make a change to pressure correction equation, but I am afraid that it is not possible in Fluent. right?

Please, any advice!

best regards
bohis

surf model
 

What is your model based on? Has it been published (just the model etc?) It would be interesting to see your surf. Ten. Force model.

just brief description
 

Hi, it is based on Height Functions and I reconstruct the length of interface in 2D cell. Ten. force is applied only where interface is present.
and what is your view on my mass imbalance in very few cells - spurious velocity?

cia!
Honza

Hi bohis.
Its been a long time now since you posted this thread, but can you share your surface tension model for fluent??? Actually I have been doing some simulations and parasitic currents at the interface have troubled me for some time. Your help would be of great value to me.
Regards,
Indrajit

surface tension
 

Hello!
U are right, it has been really a long time since I last did something with surface tension. Fluent calculates normals and curvatures from gradient of volume fraction and we can say that this is very inaccurate! I used so-called Height Functions to tackle this. The second thing is to formulate surface tension force as a volumetric force. I suggest U use Brackbills approach. Even with this strategy U should get better results. good luck!

Hi Bohis,
Even if I calculate the normals by height functions, how can I modify the surface tension expression in fluent?? Fluent uses Brackbill's model but I have found out that with this model parasitic currents increase for small drop radius and low speeds!!! These currents totally alter the shape of my drop giving unrealistic results!! Please suggest me something I could possibly do. Is there a way to modify fluent's surface tension expression???
Regards,
Indrajit

Hi,

Switch of Fluent surface tension model and write it youself (DEFINE_SOURCE(c,t)) I will do the same what FLUENT does, but U will use normals and curvatures computed using HEIGHT FUNCTIONS. Note that U can apply HF only on structured meshes...good luck, I cannot help you more..

Thanks Bohis,
I will see how that works out!! As also could you please tell me how do I change the smoothing value of VOF in fluent??? I read that after ('rpsetvar patch/vof? #t) I must set smoothing value to 2 to avoid spuriou currents.However I have not been able to find how it is done?
Regards,
Indrajit

Hi Bohis,
I wrote the UDF for height functions but it gives me segmentation violation error. I put the calculation of height functions and curvature in DEFINE_EXECUTE_AT_END so that they are not calculated every iteration ( DEFINE_SOURCE might do that!!) Please take a look at the following UDF and let me know how can I deal with segmentation violation error. I guess nesting of loops is causing the trouble.:
DEFINE_EXECUTE_AT_END(heightfunc)
{ cell_t cell1,c;
Thread *t,*pt,*pt1,*thread1;
Domain *d;
real xc[ND_ND],xc1[ND_ND],xc2[ND_ND];
d = Get_Domain(1);
real VOF,VOF1,sum=0,e=0,w=0,q=0,h1,h2,k;


/****************************************This part calulates the height functions***************************************** *************************/


thread_loop_c(thread1,d) /*loops over all cell threads in domain*/
{ printf("first change loop");
begin_c_loop_all (cell1,thread1)
{
pt= THREAD_SUB_THREAD(thread1,1);
if( C_VOF(cell1,pt)!=0 && C_VOF(cell1,pt)!=1)
{ sum=0;
C_CENTROID(xc,cell1,thread1);
VOF=C_UDMI(cell1,pt,0);
VOF1=C_UDMI(cell1,pt,1);
if(abs(VOF1)>= abs(VOF))
{ thread_loop_c(thread1,d) /*loops over all cell threads in domain*/
{
begin_c_loop_all (cell1,thread1)
{ pt1= THREAD_SUB_THREAD(thread1,1);
if( C_VOF(cell1,pt1)!=0 && C_VOF(cell1,pt1)!=1)
{ C_CENTROID(xc1,cell1,thread1);
if( xc1[0]==xc[0] && xc1[1]>=(xc[1]-(3*dely)) && xc1[1]<=(xc[1]+(3*dely)))
{ sum=sum+ (C_VOF(c,t)*dely);
}
}
}
end_c_loop_all (cell1,thread1)
}
}
C_UDMI(cell1,thread1,3)=sum;
}
}
end_c_loop_all (cell1,thread1)
}

/**********************************************This part caluclates the curvature***************************************** ******************************/

thread_loop_c(thread1,d) /*loops over all cell threads in domain*/
{ printf("second loop");
begin_c_loop_all (cell1,thread1)
{ pt1= THREAD_SUB_THREAD(thread1,1);
if( C_VOF(cell1,pt1)!=0 && C_VOF(cell1,pt1)!=1)
{ C_CENTROID(xc2,cell1,thread1);
thread_loop_c(t,d) /*loops over all cell threads in domain*/
{
begin_c_loop_all (c,t)
{ pt= THREAD_SUB_THREAD(t,1);
if( C_VOF(c,pt)!=0 && C_VOF(c,pt)!=1)
{ /*C_UDSI(c,t,0)=C_UDMI(c,t,0);
C_UDSI(c,t,1)=C_UDMI(c,t,1);
C_UDMI(c,t,6)= NV_MAG(C_UDSI_G(c,t,0))*C_UDMI(c,t,0)*0.0702*C_R(c ,t);
C_UDMI(c,t,7)=NV_MAG(C_UDSI_G(c,t,1))*C_UDMI(c,t,1 )*0.0702*C_R(c,t);*/
C_CENTROID(xc1,c,t);
if((xc2[0]+delx)==xc1[0])
{ e= C_UDMI(c,t,3);
}
if((xc2[0]-delx)==xc1[0])
{ w= C_UDMI(c,t,3);
}
if(xc2[0]==xc[1])
{ q= C_UDMI(c,t,3);
}
}
}
end_c_loop_all (c,t)
}
h1=(e-w)/delx;
h2=(e-(2*q)+w)/(delx*delx);
k=h2/pow((1+(h1*h1)),1.5);
C_UDMI(cell1,thread1,4) = 0.0702*k*(2e-6)*(h1)/volume;
C_UDMI(cell1,thread1,5)= 0.0702*k*(2e-6)*(-1)/volume;
}
}
end_c_loop_all (cell1,thread1)
}

}

Were you able to solve the segmentation error in your UDF?

You haven't specified VOF, VOF1, delx and dely anywhere in your code. Can you explain what they are?

Could you please tell me what C_UDMI(c,t,0) and C_UDMI(c,t,1) are in your UDF? They seem to have appeared from nowhere. I see no calculation done for these two terms.

Hi Indrajit,

Can you share the working correct form of this code. I will really appriciate for this.

Regards

Hi

A short hint for those who want do use the height function. The calculation for that is of order O(n^2). For large meshes this can take heaps of time. There is a function built in called Level-Set which does pretty much the same.

Hi IndrajitW,
Have you got the correct result with the UDF you have mentioned? I have been working on the wick-vapor interface for a long time but not with a good answer. Will you give me some hint?
Bill