|
[Sponsors] |
dynamic contact angle udf returns no value to solver |
|
LinkBack | Thread Tools | Search this Thread | Display Modes |
November 1, 2012, 11:32 |
dynamic contact angle udf returns no value to solver
|
#1 |
New Member
Join Date: Aug 2011
Posts: 3
Rep Power: 15 |
I wrote dynamic contact angle for wall in VOF method but the function most of times returns no value to the solver . i even try to write simplest case of udf az follow bur the contact angle I see is still wrong here the complete code:
#include "udf.h" #include "sg_mphase.h" #include "sg_vof.h" #include "sg.h" #include "mem.h" #define VISCOSITY 0.001 #define SURF_TENS 0.0728 #define MYTRUE 1 #define MYFALSE 0 #define Hoff(x) acos(1-(2.0*tanh(5.16*(pow((x/(1+(1.31*pow(x, 0.99)))), 0.706))))) #define static_Con_Ang 30. #define index_source 3 /* This Code computes the normals of the VOF function*/ DEFINE_ADJUST(store_gradient, domain) { Thread *t; Thread **pt; cell_t c; int phase_domain_index = 1; /* Secondary Domain */ Domain *pDomain = DOMAIN_SUB_DOMAIN(domain,phase_domain_index); void calc_source(); Alloc_Storage_Vars(pDomain,SV_VOF_RG,SV_VOF_G,SV_N ULL); 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,domain,pt) if (FLUID_THREAD_P(t)) { Thread *ppt = pt[phase_domain_index]; begin_c_loop (c,t) { calc_source(c,t); } end_c_loop (c,t) } Free_Storage_Vars(pDomain,SV_VOF_RG,SV_VOF_G,SV_NU LL); } void calc_source(cell_t cell, Thread *thread) { real VOF_Val[3], VOF_Mag, source, VOF_Norm[3]; Thread *phaset; phaset= THREAD_SUB_THREAD(thread,1); if(C_VOF(cell,phaset)!=0.0 && N_TIME > 1) { /* The gradients of the VOF function are found in the x,y and z dir. */ if (NULLP(THREAD_STORAGE(phaset, SV_VOF_G))) { Message0("N_TIME = %d, ....show-grad:Gradient of VOF is not available \n ", N_TIME); Error("0"); } VOF_Val[0]=-C_VOF_G(cell,phaset)[0]; VOF_Val[1]=-C_VOF_G(cell,phaset)[1]; VOF_Val[2]=-C_VOF_G(cell,phaset)[2]; /* The magnitude of the VOF gradients is found so it can be normalized */ VOF_Mag=NV_MAG(VOF_Val); if(VOF_Mag!=0.0) { VOF_Mag=NV_MAG(VOF_Val); VOF_Norm[0]=VOF_Val[0]/VOF_Mag; VOF_Norm[1]=VOF_Val[1]/VOF_Mag; VOF_Norm[2]=0; } else { /* This is to avoid the divide by zero function*/ VOF_Norm[0]=0.0; VOF_Norm[1]=0.0; VOF_Norm[2]=0.0; } C_UDMI(cell,thread,0)=VOF_Norm[0]; C_UDMI(cell,thread,1)=VOF_Norm[1]; C_UDMI(cell,thread,2)=0; } source = 0.0; C_UDMI(cell, thread, index_source) = source; } DEFINE_SOURCE(VOF_Norms, cell, thread, dS, eqn) { real source; source = C_UDMI(cell, thread, index_source); return source; } /* This Define_profile code is designed to provide a dynamic contact angle for the VOF function*/ DEFINE_PROFILE(con_ang, t, pos) { /* First the various pointer variables are created*/ face_t f; cell_t c; real feta_d, vel_Val[3], cont_Line_Vel, VOF_Normal[3], cap_Num, static_Con_Rad, x_Bottom,h, x_Top, x_Bisect, hoff_Old, hoff_Cur, hoff_New, finish_Cond, inv_Hoff=0.0; int notConverged, itNum; Thread *t0,*pt; /* This code is designed to find the zero for the inverted hoffman function by finding the zero */ /* of the function at which the hoffman function results in the static contact angle */ /* First the variables are assigned*/ notConverged=MYTRUE; x_Bottom=0.001; x_Top=2.0; itNum=0; static_Con_Rad=((static_Con_Ang*M_PI)/180.); /* A while loop performs the bisection method, a simple but very stable zero finder */ while(notConverged) { /* The variables used in the bisection method are assigned and the hoffman */ /* functions are evaluated */ itNum++; hoff_Old=(Hoff(x_Bottom)- static_Con_Rad); hoff_Cur=(Hoff(x_Top)- static_Con_Rad); x_Bisect=(x_Bottom+x_Top)/2.0; hoff_New=(Hoff(x_Bisect)- static_Con_Rad); finish_Cond=fabs(1-(x_Bisect/x_Top)); /* The loop ends when the relative error is less than 1e-8 and the inverse */ /* hoffman value is stored for use later */ if(finish_Cond<0.00000001 || itNum>10000000) { inv_Hoff=x_Bisect; notConverged=MYFALSE; } /* Conditions for the bisection method */ if((hoff_Old*hoff_New)<0.0) { x_Top=x_Bisect; } if((hoff_Cur*hoff_New)<=0.0) { x_Bottom=x_Bisect; } } /* Now the main loop goes through all the faces in the boundary */ begin_f_loop(f,t) { /* The cell and phase threads are isolated */ c=F_C0(f,t); t0=THREAD_T0(t); pt= THREAD_SUB_THREAD(t0,1); /* The main formulation is only applied if the VOF is >0 */ if(C_VOF(c,pt)!=0.0) { /* The velocities are recorded in each direction */ vel_Val[0]=C_U(c,t0); vel_Val[1]=C_V(c,t0); vel_Val[2]=0; /* The VOF normals are brought in */ VOF_Normal[0]=C_UDMI(c,t0,0); VOF_Normal[1]=C_UDMI(c,t0,1); VOF_Normal[2]=0; /* The contact line vel. is calc from the dot product of VOF and Vel */ cont_Line_Vel=NV_DOT(vel_Val,VOF_Normal); /* The capillary number is found based on cont line vel. */ cap_Num=((VISCOSITY*cont_Line_Vel)/SURF_TENS); /* The dynamic contact angle is defined then stored in the profile */ C_UDMI(c,t0,4)=cont_Line_Vel; if ( cont_Line_Vel>0.000001) { h=48; } else if (cont_Line_Vel<-0.000001 ) { h=10; } else { h=30; } C_UDMI(c,t0,5)=h; feta_d=h; F_PROFILE(f,t,pos)=feta_d; } else { F_PROFILE(f,t,pos)=static_Con_Ang; } C_UDMI(c,t0,6)=F_PROFILE(f,t,pos); } end_f_loop(f,t) } and the simplest case code with just a constant value: #include "udf.h" DEFINE_PROFILE(CA_profile,t,i) { face_t f; begin_f_loop_all(f,t) { F_PROFILE(f,t,i)=30.0; } end_f_loop_all(f,t) } |
|
February 23, 2018, 10:45 |
|
#2 |
New Member
Join Date: Feb 2018
Posts: 6
Rep Power: 8 |
Hi shiraz, have you solved your problem? I have a similar problem as you did.
|
|
July 1, 2018, 04:33 |
|
#3 | |
New Member
Join Date: Jul 2018
Posts: 7
Rep Power: 8 |
Quote:
I have being stuck in this problem for several days! Nothing help could come out! If you have got the correct method for setting for coding the udf of dynamic contact angle, please give me a reply. Thank you very much! |
||
July 1, 2018, 17:01 |
|
#4 |
New Member
Join Date: Feb 2018
Posts: 6
Rep Power: 8 |
changing the unit of contact angle to rad in your code might help
|
|
July 3, 2018, 03:11 |
|
#5 |
New Member
Join Date: Jul 2018
Posts: 7
Rep Power: 8 |
It really works! But I don't know why this could happen. And FLUENT didn't illustrate it in help files. Maybe it's a bug? Anyway,
Thank you very much! Albert Lee |
|
July 3, 2018, 14:51 |
|
#6 |
New Member
Join Date: Feb 2018
Posts: 6
Rep Power: 8 |
yes, it's weird. And seems no one points it out clearly..... don't know why...
|
|
Thread Tools | Search this Thread |
Display Modes | |
|
|
Similar Threads | ||||
Thread | Thread Starter | Forum | Replies | Last Post |
dynamic contact angle implementation q | doubtsincfd | OpenFOAM Programming & Development | 7 | November 12, 2015 08:07 |
Slug Flow, interFoam, problems with Contact Angle | PrzemekPL | OpenFOAM Running, Solving & CFD | 13 | February 18, 2014 22:10 |
Dynamic Mesh solver | OF_NACA | OpenFOAM | 0 | June 21, 2011 03:03 |
Working directory via command line | Luiz | CFX | 4 | March 6, 2011 20:02 |
Dynamic contact angle with Fluent 6.3.26 | baechtel | FLUENT | 0 | May 21, 2009 23:55 |