CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > Software User Forums > ANSYS > FLUENT > Fluent UDF and Scheme Programming

Calculating flow variables only once in begin_f_loop()

Register Blogs Members List Search Today's Posts Mark Forums Read

Reply
 
LinkBack Thread Tools Search this Thread Display Modes
Old   September 24, 2019, 15:15
Default Calculating flow variables only once in begin_f_loop()
  #1
New Member
 
Join Date: Sep 2019
Posts: 24
Rep Power: 6
NonStopEagle is on a distinguished road
Hello everyone,
I have been trying to write a code for implementing maxwell first order slip boundary conitions at the diverging wall of a micro-nozzle. My code works well initially but eventually succumbs to sudden divergence. I believe i have figured out why this happens, the code is as follows :


#include "udf.h" /* must be at the beginning of every UDF */
/*#include "mem.h" /* must be at the beginning of this UDF */
//#include <math.h> /* must be at the beginning of this UDF */
/*#include "metric.h" /* must be at the beginning of this UDF */
#define k_B 1.3805e-23 /* Boltzman constant (j/K)*/
#define sigma 419e-12 /* molecular diameter*/
#define sigma_v 1 /* Tangential momentum accommodation coefficient*/
#define sigma_T 1 /*Thermal accommodation coefficient*/
#define pi 3.14159 /* pi number*/
#define sqrt_2 1.41421 /* sqrt(2) number */
#define len = 20e-5;
double relax = 0.1;
int counter = 0;
//#define relax 50e-3
DEFINE_PROFILE(wall_slip_velocity, t, index)
{
double T, P, density, MU, K, landa, a_m, du_da, dT_dx, u_f, DUY, DVY, DUX, DVX, u, v;
double A[ND_ND];
face_t f;
Thread* t0;
cell_t c0;
double DUR, DVR;


begin_f_loop(f, t)
{
t0 = THREAD_T0(t); /* adjacent cell thread to f */
c0 = F_C0(f, t); /* adjacent cell to f */
T = C_T(c0, t0); /* temperature of the cell */
P = C_P(c0, t0); /* peressure of the cell */
density = C_R(c0, t0); /* density of the cell*/
MU = C_MU_L(c0, t0); /* laminar viscosity of the cell */
u = C_U(c0, t0);
v = C_V(c0, t0);
/*DUR = C_U_RG(c0, t0)[1];
DVR = C_V_RG(c0, t0)[1];*/

K = C_K_L(c0, t0); /* thermal conductivity of the cell */
landa = (k_B * T) / (sqrt_2 * pi * sigma * sigma * P); /* mean free path line */

DUX = C_DUDX(c0, t0);
DVX = C_DVDX(c0, t0);
DUY = C_DUDY(c0, t0);
DVY = C_DVDY(c0, t0);
F_AREA(A, f, t);
a_m = NV_MAG(A); /* magnitude of the cell area vector */
du_da = NV_DOT(A, A) / a_m; /* n- component of the cell x-velocity reconstruction gradient(RG) vector */
dT_dx = 1; /* x-component of the cell temperature reconstruction gradient(RG) vector */
//F_PROFILE(f, t, index) = relax * (((2 - sigma_v) / sigma_v) * landa * (u / sqrt(u * u + v * v) * (((DUX + DVX) * (v / sqrt(u * u + v * v))) + ((DUY + DVY) * (u / sqrt(u * u + v * v)))) + 0.75 * (MU / (density * T)) * dT_dx);
//F_UDMI(f,t,1)= relax * (((2 - sigma_v) / sigma_v) * landa * (u / sqrt(u * u + v * v) * (((DUX + DVX) * (v / sqrt(u * u + v * v))) + ((DUY + DVY) * (u / sqrt(u * u + v * v)))) + 0.75 * (MU / (density * T)) * dT_dx);
F_UDMI(f, t, 0) = relax * (((2 - sigma_v) / sigma_v) * landa * (((DUX + DVX) * (v / sqrt(u * u + v * v))) - ((DUY + DVY) * (u / sqrt(u * u + v * v)))));
// F_PROFILE(f, t, index) = F_UDMI(f, t, 0);
//F_PROFILE(f, t, index) = F_UDSI(f, t, 0);
F_UDMI(f, t, 1) = landa;

}
end_f_loop(f, t)

}
DEFINE_PROFILE(temperature_jump, t, index)
{
face_t f;
cell_t c0;
Thread* t0;
double p, tem, g, lamda, num_den, temp, vis, temp1, x, i, mag, DT2;
double DT, prd, tc, tj, cp, a,u,v;
double QUANT;
g = 9.81;
a = 4.29866e-19;
begin_f_loop(f, t) // for each face: get face-id and thread
{
t0 = THREAD_T0(t); /* adjacent cell thread to f */
c0 = F_C0(f, t); /* adjacent cell to f */
p = C_P(c0, t0); // getting the pressure
tem = C_T(c0, t0); // getting the temperature
tc = C_K_L(c0, t0);
vis = C_MU_L(c0, t0);
cp = C_CP(c0, t0);
u = C_U(c0, t0);
v = C_V(c0, t0);
DT = C_T_G(c0, t0)[1];
DT2 = C_T_G(c0, t0)[0];
prd = cp * vis / (g * tc);

num_den = p / ((1.38064852e-23) * tem); // Calculation of the number density
temp1 = num_den * a;
lamda = 1 / (temp1);
QUANT = (((DT2) * (v / sqrt(u * u + v * v))) - ((DT) * (u / sqrt(u * u + v * v))));
tj = 380+relax *((((2 - sigma_v) / sigma_v) *1.14566*lamda * QUANT ));//0.15596 (380 + ((1.14566 * 1 * (lamda / prd) * ((0.9877 * DT)))) * 0.9877);
F_PROFILE(f, t, index) = tj;


}
end_f_loop(f, t)
}
DEFINE_EXECUTE_AT_END(define_relax)
{
if ((counter % 200) == 0)
{
relax = relax + 0.01;
printf("Relaxation : %f", relax);
}
counter = counter + 1;
}


In the code the gradients are calculated everytime the begin_f_loop() runs, and every change of the gradient values brings some error in the resulting slip velocity. This erred slip velocity produces further deviations in the gradients, until finally the solution becomes unstable and diverges.
The solution to this, as i believe is to run the simulation in no-slip conditions till convergence and then implement slip such that the gradient and other flow variables are calculated only once initially at all the faces of the slip wall and then to use those values to calculate the slip velocity.
And i need your help to accomplish this.
Any help regarding this will be highly appreciated.
thanks in advance.
NonStopEagle is offline   Reply With Quote

Old   May 5, 2020, 05:29
Default
  #2
Member
 
Homayoon sohrabi
Join Date: May 2020
Posts: 56
Rep Power: 6
homay95 is on a distinguished road
Quote:
Originally Posted by NonStopEagle View Post
Hello everyone,
I have been trying to write a code for implementing maxwell first order slip boundary conitions at the diverging wall of a micro-nozzle. My code works well initially but eventually succumbs to sudden divergence. I believe i have figured out why this happens, the code is as follows :


#include "udf.h" /* must be at the beginning of every UDF */
/*#include "mem.h" /* must be at the beginning of this UDF */
//#include <math.h> /* must be at the beginning of this UDF */
/*#include "metric.h" /* must be at the beginning of this UDF */
#define k_B 1.3805e-23 /* Boltzman constant (j/K)*/
#define sigma 419e-12 /* molecular diameter*/
#define sigma_v 1 /* Tangential momentum accommodation coefficient*/
#define sigma_T 1 /*Thermal accommodation coefficient*/
#define pi 3.14159 /* pi number*/
#define sqrt_2 1.41421 /* sqrt(2) number */
#define len = 20e-5;
double relax = 0.1;
int counter = 0;
//#define relax 50e-3
DEFINE_PROFILE(wall_slip_velocity, t, index)
{
double T, P, density, MU, K, landa, a_m, du_da, dT_dx, u_f, DUY, DVY, DUX, DVX, u, v;
double A[ND_ND];
face_t f;
Thread* t0;
cell_t c0;
double DUR, DVR;


begin_f_loop(f, t)
{
t0 = THREAD_T0(t); /* adjacent cell thread to f */
c0 = F_C0(f, t); /* adjacent cell to f */
T = C_T(c0, t0); /* temperature of the cell */
P = C_P(c0, t0); /* peressure of the cell */
density = C_R(c0, t0); /* density of the cell*/
MU = C_MU_L(c0, t0); /* laminar viscosity of the cell */
u = C_U(c0, t0);
v = C_V(c0, t0);
/*DUR = C_U_RG(c0, t0)[1];
DVR = C_V_RG(c0, t0)[1];*/

K = C_K_L(c0, t0); /* thermal conductivity of the cell */
landa = (k_B * T) / (sqrt_2 * pi * sigma * sigma * P); /* mean free path line */

DUX = C_DUDX(c0, t0);
DVX = C_DVDX(c0, t0);
DUY = C_DUDY(c0, t0);
DVY = C_DVDY(c0, t0);
F_AREA(A, f, t);
a_m = NV_MAG(A); /* magnitude of the cell area vector */
du_da = NV_DOT(A, A) / a_m; /* n- component of the cell x-velocity reconstruction gradient(RG) vector */
dT_dx = 1; /* x-component of the cell temperature reconstruction gradient(RG) vector */
//F_PROFILE(f, t, index) = relax * (((2 - sigma_v) / sigma_v) * landa * (u / sqrt(u * u + v * v) * (((DUX + DVX) * (v / sqrt(u * u + v * v))) + ((DUY + DVY) * (u / sqrt(u * u + v * v)))) + 0.75 * (MU / (density * T)) * dT_dx);
//F_UDMI(f,t,1)= relax * (((2 - sigma_v) / sigma_v) * landa * (u / sqrt(u * u + v * v) * (((DUX + DVX) * (v / sqrt(u * u + v * v))) + ((DUY + DVY) * (u / sqrt(u * u + v * v)))) + 0.75 * (MU / (density * T)) * dT_dx);
F_UDMI(f, t, 0) = relax * (((2 - sigma_v) / sigma_v) * landa * (((DUX + DVX) * (v / sqrt(u * u + v * v))) - ((DUY + DVY) * (u / sqrt(u * u + v * v)))));
// F_PROFILE(f, t, index) = F_UDMI(f, t, 0);
//F_PROFILE(f, t, index) = F_UDSI(f, t, 0);
F_UDMI(f, t, 1) = landa;

}
end_f_loop(f, t)

}
DEFINE_PROFILE(temperature_jump, t, index)
{
face_t f;
cell_t c0;
Thread* t0;
double p, tem, g, lamda, num_den, temp, vis, temp1, x, i, mag, DT2;
double DT, prd, tc, tj, cp, a,u,v;
double QUANT;
g = 9.81;
a = 4.29866e-19;
begin_f_loop(f, t) // for each face: get face-id and thread
{
t0 = THREAD_T0(t); /* adjacent cell thread to f */
c0 = F_C0(f, t); /* adjacent cell to f */
p = C_P(c0, t0); // getting the pressure
tem = C_T(c0, t0); // getting the temperature
tc = C_K_L(c0, t0);
vis = C_MU_L(c0, t0);
cp = C_CP(c0, t0);
u = C_U(c0, t0);
v = C_V(c0, t0);
DT = C_T_G(c0, t0)[1];
DT2 = C_T_G(c0, t0)[0];
prd = cp * vis / (g * tc);

num_den = p / ((1.38064852e-23) * tem); // Calculation of the number density
temp1 = num_den * a;
lamda = 1 / (temp1);
QUANT = (((DT2) * (v / sqrt(u * u + v * v))) - ((DT) * (u / sqrt(u * u + v * v))));
tj = 380+relax *((((2 - sigma_v) / sigma_v) *1.14566*lamda * QUANT ));//0.15596 (380 + ((1.14566 * 1 * (lamda / prd) * ((0.9877 * DT)))) * 0.9877);
F_PROFILE(f, t, index) = tj;


}
end_f_loop(f, t)
}
DEFINE_EXECUTE_AT_END(define_relax)
{
if ((counter % 200) == 0)
{
relax = relax + 0.01;
printf("Relaxation : %f", relax);
}
counter = counter + 1;
}


In the code the gradients are calculated everytime the begin_f_loop() runs, and every change of the gradient values brings some error in the resulting slip velocity. This erred slip velocity produces further deviations in the gradients, until finally the solution becomes unstable and diverges.
The solution to this, as i believe is to run the simulation in no-slip conditions till convergence and then implement slip such that the gradient and other flow variables are calculated only once initially at all the faces of the slip wall and then to use those values to calculate the slip velocity.
And i need your help to accomplish this.
Any help regarding this will be highly appreciated.
thanks in advance.
hello! do you know how i could write udf cod for these functions below:
1. V= L(du/dy)
2. T)f = T)s + L(dT/dy)
T)f is fluid temperature and T)s is solid temperature and L is slip length
i would appreciate for any help. thanks
homay95 is offline   Reply With Quote

Old   May 5, 2020, 06:35
Default Slip Velocity
  #3
Senior Member
 
vinerm's Avatar
 
Vinerm
Join Date: Jun 2009
Location: Nederland
Posts: 2,946
Blog Entries: 1
Rep Power: 36
vinerm will become famous soon enough
You can use the code that you quoted in your last as well as current post. However, there is one common mistake in both of these codes; both use velocity to define slip boundary. That brings inconsistency and requires extra care from user end. The approach is to use flux profile, i.e., shear stress for velocity slip and heat flux for temperature slip. But if you are not well versed with UDFs, then, I'd recommend you use the code that you quoted in your last post. You can modify it as per your equations.
__________________
Regards,
Vinerm

PM to be used if and only if you do not want something to be shared publicly. PM is considered to be of the least priority.
vinerm is offline   Reply With Quote

Old   May 6, 2020, 04:31
Default
  #4
Member
 
Homayoon sohrabi
Join Date: May 2020
Posts: 56
Rep Power: 6
homay95 is on a distinguished road
First of all thank you so much for your response, and yes I just started a project on michrochannels in ANSYS-FLUENT and I'm not professional in writing UDF codes. so you are saying that I can use this code as velocity slip and temperature jump for a solid-liquid interface?
homay95 is offline   Reply With Quote

Old   May 6, 2020, 04:34
Default Code
  #5
Senior Member
 
vinerm's Avatar
 
Vinerm
Join Date: Jun 2009
Location: Nederland
Posts: 2,946
Blog Entries: 1
Rep Power: 36
vinerm will become famous soon enough
The code that you quoted in your first post. But do note that you will have to modify the equation being used since you want to use a different equation.
__________________
Regards,
Vinerm

PM to be used if and only if you do not want something to be shared publicly. PM is considered to be of the least priority.
vinerm is offline   Reply With Quote

Old   May 6, 2020, 04:39
Default
  #6
Member
 
Homayoon sohrabi
Join Date: May 2020
Posts: 56
Rep Power: 6
homay95 is on a distinguished road
I understand but as I'm not so into programming, could you please be more specific? I mean which part of the code should be changed according to my equations?
Have you written a UDF for this (Maxwell's boundary condition) before?
homay95 is offline   Reply With Quote

Old   May 6, 2020, 04:50
Default Udf
  #7
Senior Member
 
vinerm's Avatar
 
Vinerm
Join Date: Jun 2009
Location: Nederland
Posts: 2,946
Blog Entries: 1
Rep Power: 36
vinerm will become famous soon enough
The UDF given at

UDF for slip boundary condition

is rather well documented. Final value applied is whatever is given on rhs of F_PROFILE. If you have C or C++ before, then modification would be easy. Even if you have not, then you can remove some of the lines within begin_f_loop{}end_f_loop but keep F_PROFILE line. Change the equation on its rhs with your equation. The gradient you need for your equation is already given in that code as strain rate magnitude.
__________________
Regards,
Vinerm

PM to be used if and only if you do not want something to be shared publicly. PM is considered to be of the least priority.
vinerm is offline   Reply With Quote

Old   May 6, 2020, 05:02
Default
  #8
Member
 
Homayoon sohrabi
Join Date: May 2020
Posts: 56
Rep Power: 6
homay95 is on a distinguished road
okay thank you very much for your help. and there is just one more question: do you know how i could interpret more than one UDF to ANSYS-FLUENT?
homay95 is offline   Reply With Quote

Old   May 6, 2020, 05:06
Default Interpretation
  #9
Senior Member
 
vinerm's Avatar
 
Vinerm
Join Date: Jun 2009
Location: Nederland
Posts: 2,946
Blog Entries: 1
Rep Power: 36
vinerm will become famous soon enough
You can interpret only one file, however, no body stops you from putting any number of functions within the file. So, just put everything in single file.
__________________
Regards,
Vinerm

PM to be used if and only if you do not want something to be shared publicly. PM is considered to be of the least priority.
vinerm is offline   Reply With Quote

Old   May 6, 2020, 05:18
Default
  #10
Member
 
Homayoon sohrabi
Join Date: May 2020
Posts: 56
Rep Power: 6
homay95 is on a distinguished road
okay. now I get it. and I'm really sorry for taking your time but as I said I'm starter in ANSYS-FLUENT, so I have so much questions. where should I apply this UDF for interface? I mean for the wall or for the shadow-wall between solid and liquid?
homay95 is offline   Reply With Quote

Old   May 6, 2020, 05:23
Default Wall and Shadow
  #11
Senior Member
 
vinerm's Avatar
 
Vinerm
Join Date: Jun 2009
Location: Nederland
Posts: 2,946
Blog Entries: 1
Rep Power: 36
vinerm will become famous soon enough
Open the wall boundary and it will show you the cell zone it is associated with. You have to assign to the wall or its shadow, whichever is associated with the fluid zone.
__________________
Regards,
Vinerm

PM to be used if and only if you do not want something to be shared publicly. PM is considered to be of the least priority.
vinerm is offline   Reply With Quote

Old   May 6, 2020, 05:42
Default
  #12
Member
 
Homayoon sohrabi
Join Date: May 2020
Posts: 56
Rep Power: 6
homay95 is on a distinguished road
So I should apply the UDF twice, once in the "Momentum-specified shear" as slip velocity and once in "Thermal-Temperature" as temperature jump, and after that should I leave the solid adjacent wall as "Coupled", is it right?
homay95 is offline   Reply With Quote

Old   May 6, 2020, 05:48
Default Wall-Shadow Boundary
  #13
Senior Member
 
vinerm's Avatar
 
Vinerm
Join Date: Jun 2009
Location: Nederland
Posts: 2,946
Blog Entries: 1
Rep Power: 36
vinerm will become famous soon enough
Once you change it from coupled to temperature, it changes for both, i.e., for wall and its shadow.
__________________
Regards,
Vinerm

PM to be used if and only if you do not want something to be shared publicly. PM is considered to be of the least priority.
vinerm is offline   Reply With Quote

Old   May 6, 2020, 06:10
Default
  #14
Member
 
Homayoon sohrabi
Join Date: May 2020
Posts: 56
Rep Power: 6
homay95 is on a distinguished road
Sorry but i didn't get it, after interpreting the UDF which includes both slip velocity and temperature jump in one file, where should I specifically assign this UDF file in the fluid zone. I mean for Momentum and Thermal seperately.
homay95 is offline   Reply With Quote

Old   May 6, 2020, 06:40
Default Functions
  #15
Senior Member
 
vinerm's Avatar
 
Vinerm
Join Date: Jun 2009
Location: Nederland
Posts: 2,946
Blog Entries: 1
Rep Power: 36
vinerm will become famous soon enough
Though you have single file, there are two functions, each with different name. You will see both names under all boundaries. Just select the name you want to apply.
__________________
Regards,
Vinerm

PM to be used if and only if you do not want something to be shared publicly. PM is considered to be of the least priority.
vinerm is offline   Reply With Quote

Old   May 6, 2020, 06:59
Default
  #16
Member
 
Homayoon sohrabi
Join Date: May 2020
Posts: 56
Rep Power: 6
homay95 is on a distinguished road
okay then. thank you very very much for your help. I really owe you and sorry again for taking your time
best regards,
Homayoon
homay95 is offline   Reply With Quote

Old   June 6, 2020, 08:36
Default
  #17
Member
 
Homayoon sohrabi
Join Date: May 2020
Posts: 56
Rep Power: 6
homay95 is on a distinguished road
hello again. i've written a udf code for velocity slip boundary condition as follows:
#include "udf.h"
double dudy ,vel;
double current_vel;
DEFINE_PROFILE(slipvelocity, thread, position)
{
face_t f;
cell_t c;
Thread *tc;
double b = 0.000001;

begin_f_loop(f, thread)
{
c = F_C0(f,thread);
tc = THREAD_T0(thread);
dudy = C_U_G(c,tc)[1];
vel = b*dudy;
F_PROFILE(f, thread, position) = vel;
}
end_f_loop(f, thread)
}

is it correct? or does it need any modification?
homay95 is offline   Reply With Quote

Old   June 6, 2020, 10:37
Default Code
  #18
Senior Member
 
vinerm's Avatar
 
Vinerm
Join Date: Jun 2009
Location: Nederland
Posts: 2,946
Blog Entries: 1
Rep Power: 36
vinerm will become famous soon enough
The code is syntactically correct, except that you have global variables. Define all variables within DEFINE_ macro and not outside it.
__________________
Regards,
Vinerm

PM to be used if and only if you do not want something to be shared publicly. PM is considered to be of the least priority.
vinerm is offline   Reply With Quote

Old   June 6, 2020, 10:45
Default
  #19
Member
 
Homayoon sohrabi
Join Date: May 2020
Posts: 56
Rep Power: 6
homay95 is on a distinguished road
okay. thank you very much for your help. one more question is should i hook it as shear stress or moving wall velocity component in boundary conditions?
homay95 is offline   Reply With Quote

Old   June 6, 2020, 11:05
Default Boundary Condition
  #20
Senior Member
 
vinerm's Avatar
 
Vinerm
Join Date: Jun 2009
Location: Nederland
Posts: 2,946
Blog Entries: 1
Rep Power: 36
vinerm will become famous soon enough
Depends on what the UDF returns. Better approach is to write UDF for shear stress, however, your code returns velocity and should be used for velocity.
__________________
Regards,
Vinerm

PM to be used if and only if you do not want something to be shared publicly. PM is considered to be of the least priority.
vinerm is offline   Reply With Quote

Reply

Tags
fluent - udf - parallel, slip b.c., slip flow, udf and programming

Thread Tools Search this Thread
Search this Thread:

Advanced Search
Display Modes

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off
Trackbacks are Off
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
UDF program for calculating the mass flow outlet with a pulsed output flow Benito89 Fluent UDF and Scheme Programming 9 August 21, 2019 23:57
Match Pressure Inlet/Outlet Boundary Condition Mass Flow Rate MSchneid Fluent UDF and Scheme Programming 3 February 23, 2019 06:00
mass flow in is not equal to mass flow out saii CFX 12 March 19, 2018 05:21
Calculating mass flow rate at multiphase flows Kuslo187 OpenFOAM Post-Processing 1 August 21, 2015 18:11
How to update polyPatchbs localPoints liu OpenFOAM Running, Solving & CFD 6 December 30, 2005 17:27


All times are GMT -4. The time now is 10:59.