[arunraj]

Dear ALL,

How to find the index of facets lie in a boundary. I need to use it to call the face index of shadow face. I cannot THREAD_SHADOW as my boundary is velocity inlet condition. Please tell a method to workaround this issue.

begin_f_loop(f,t)
{
F_CENTROID(x,f,t);
xx=x[0];
c0 = F_C0(f,t);
t0 = THREAD_T0(t);
f_shadow = f_wick_vapor[facet_index];
Message("f_shadow is %d\n",f_shadow);
c0_shadow = F_C0(f_shadow,t_shadow);
t0_shadow = THREAD_T0(t_shadow);

[AlexanderZ]

as far as I know, velocity inlet may not have shadow walls, it is external boundary condition. It means, you may check its face ID in fluent GUI

[arunraj]

Thank you Alexander for your kind reply. I need to use it in UDF. How can I use? The author I am following has called the face index in his work. He used array to get the ID.
using ff_shadow=ff_VAPOR_WICK [iii];

However, I am not able to successfully implement this into mine. Please let me know the work around for this problem.

I have pasted only part of his UDF.

Numerical Analysis Of Phase Change, Heat Transfer And Fluid Flow Within Miniature Heat Pipes.

declaration:

face_t ff_VAPOR_WICK [1000];

DEFINE_INIT(SHADOW_READING, domain)
{
zone_ID = MARZ_VAPOR_WICK_ID;
tt= Lookup_Thread(domain,zone_ID);
zone_ID = MARZ_WICK_VAPOR_ID;
tt_shadow = Lookup_Thread(domain,zone_ID);
iii=-1;
begin_f_loop(ff, tt)
{
iii=iii+1;
F_CENTROID(x,ff,tt);
xx = x[0];
jjj=-1;
begin_f_loop(ff_shadow, tt_shadow)
{
jjj=jjj+1;
F_CENTROID(x,ff_shadow, tt_shadow);
xx_shadow = x[0];
if (fabs((xx_shadow-xx)/TOTAL_LENGHT) < 1.0E-6)
{
ff_VAPOR_WICK [iii]=ff_shadow;
iii_VAPOR_WICK[jjj]=iii;
}
}
end_f_loop(ff_shadow, tt_shadow)
}
end_f_loop(ff, tt)
}

DEFINE_ADJUST(WICK_VAPOR_INTERFACE, domain)
{
zone_ID = MARZ_VAPOR_WICK_ID;
tt = Lookup_Thread(domain,zone_ID);
zone_ID = MARZ_WICK_VAPOR_ID;
tt_shadow = Lookup_Thread(domain,zone_ID);
iii=-1;
begin_f_loop(ff, tt)
{
iii=iii+1;
F_CENTROID(x,ff,tt);
xx = x[0];
ff_shadow=ff_VAPOR_WICK [iii];
c0 = F_C0(ff,tt);
t0 = F_C0_THREAD(ff,tt);
temp_Cell=C_UDSI(c0,t0,0);
temp_face=F_UDSI(ff,tt,0);
k_Cell=C_K_L(c0,t0);
ro_Cell=C_R(c0,t0);
cp_Cell=C_CP(c0,t0);
v_face=F_V(ff,tt);
P_Cell=C_P(c0,t0)-PRESSURE_VAPOR_0;
BOUNDARY_FACE_GEOMETRY(ff,tt,A,ds,es,A_by_Es,dr0);
c0_shadow = F_C0(ff_shadow,tt_shadow);
t0_shadow = F_C0_THREAD(ff_shadow,tt_shadow);
temp_Cell_shadow=C_UDSI(c0_shadow,t0_shadow,0);
temp_face_shadow=F_UDSI(ff_shadow,tt_shadow,0);
k_Cell_shadow=C_K_L(c0_shadow,t0_shadow);
ro_Cell_shadow=C_R(c0_shadow,t0_shadow);
cp_Cell_shadow=C_CP(c0_shadow,t0_shadow);
v_face_shadow=F_V(ff_shadow,tt_shadow);
BOUNDARY_FACE_GEOMETRY(ff_shadow,
tt_shadow,A,ds_shadow,ess,A_by_Ess,dr0s);

[AlexanderZ]

where this code comes from?
few variables are not defined at all. On the other hand, looks like this is just a part of code.

UDF uses MARZ_VAPOR_WICK_ID and MARZ_WICK_VAPOR_ID to define which face it is, this information comes from Fluent graphic interface, you may check ID there.

Frankly speaking, I don't understand the logic of this code, why the loop inside loop was used.
If you need to know the cell which is adjusted to interface (face) you wanna use c0/t0 from one side and c1/t1 from the other side, ID of that face comes from fluent GUI


AND ONCE AGAIN, *from my experience* , if you have inlet boundary condition -> it is not interface, so you don't have any shadow, or whatever. Inlet is external boundary condition

make a sketch of what you exactly want to do

[arunraj]

Thank you for your kind reply. We are trying to solve the heat pipe simulation using single phase approximation. There are three zones, wall (solid), wick (porous), and vapor (Fluid). The interface between wall and wick is a coupled wall. But the interface between the wick and vapor is considered as two-sided velocity inlet condition. This is why we use the term 'SHADOW'. It is just to represent the other side. Here f and t are the face and thread of Vapor-wick interface which is vapor side of the interface.

Even the term 'interface' can also give wrong notion. It is not interface rather it is a two-sided velocity inlet. Also, even using the term velocity inlet can be wrong as it is applied in the internal boundary.

The real trick is this simulation is to model the mass transfer through vapor-wick interface and wick-vapor interface. The mass transfer is computed based on interface temperature, which is computed from energy balance. To obtain the temperature, I did face loop in vapor-wick interface but need to access the other side of face (i.e. wick-vapor interface). However, Fluent UDF has no commonly used macros for accessing the facets index on the other side.

Could you please suggest a method to access the face on other side for two-sided cases. I believe it can be done with the help "facet index". However, I may be wrong as I have only limited experience.

Thank you for showing interest in my work.

Now coming to your questions.

Q1. where this code comes from?
few variables are not defined at all. On the other hand, looks like this is just a part of code.

The code is available in the below thesis. Please have a look with the help of link.

Numerical Analysis Of Phase Change, Heat Transfer And Fluid Flow Within Miniature Heat Pipes. Mehdi Famouri, University of South Carolina.

https://scholarcommons.sc.edu/cgi/vi...37&context=etd

Q2.
UDF uses MARZ_VAPOR_WICK_ID and MARZ_WICK_VAPOR_ID to define which face it is, this information comes from Fluent graphic interface, you may check ID there.

Yes, you are correct. The IDs can be obtained using fluent graphic interface.

Q3.
Frankly speaking, I don't understand the logic of this code, why the loop inside loop was used.

Loop inside loop was used for matching the data on either side of the interface.

Q3.
If you need to know the cell which is adjusted to interface (face) you wanna use c0/t0 from one side and c1/t1 from the other side, ID of that face comes from fluent GUI

c1, t1 will not work for two-side velocity inlet. It works if I consider it as an interior. I want to apply velocity condition using DEFINE_PROFILE. I can not do it in interior.

Q4.
AND ONCE AGAIN, *from my experience* , if you have inlet boundary condition -> it is not interface, so you don't have any shadow, or whatever. Inlet is external boundary condition

make a sketch of what you exactly want to do

Yes, you are right. inlet boundary condition is not interface. However, we are just using the term to say the boundary condition we are applying at the interface. It has two-sided velocity inlet. So we are using term 'SHADOW'

[arunraj]

I am sorry Alexander. I am afraid to paste the entire code as it not my thesis. However, you can access using the link.

[AlexanderZ]

I've never heard about 2-sided inlet before.
How did you make it? does it mean, there there is no link between inlet and inlet-shadow (in your case), no heat/mass transfer ? So you can put space between WICK and VAPOR regions, and it will be same?

[arunraj]

Thank you Alexander for your time and effort. I really appreciate.

I've never heard about 2-sided inlet before.

That's why I like this problem very much. The only reason two-sided inlet is used is the velocity is entering/exiting radially which can be specified only velocity boundary condition. If I considered it as 2-sided wall, I can only apply mass source to the cells near the boundary. However, the direction of mass transfer is lost. To my knowledge, it is not possible to specify the direction of mass using DEFINE_SOURCE.

How did you make it? does it mean, there there is no link between inlet and inlet-shadow (in your case), no heat/mass transfer ?
Yes. there is no heat/mass transfer. However, it is taken care by the boundary condition. We use energy balance to provide temperature at the inlet and shadow and velocity is computed based on kinetic theory.

So you can put space between WICK and VAPOR regions, and it will be same?
I am not sure about this. But it is not physical to have space in between wick and vapor region.

Now, can you please suggest how to get the face of shadow velocity inlet. The thesis I am following has used the concept of array to get the face index. He is using ff_WICK_VAPOR [iii]=ff_shadow; to get the index of face. I believe is WICK_VAPOR is name of the boundary. But I couldn't find any other reference other than this to understand how to get face index for shadow velocity inlet.

[AlexanderZ]

the reason, why I'm asking about space between inlet and shadow-inlet is I wanna imagine, how your mesh looks like.
space is not physical, but it will not effect anyhow, if heat transfer and mass transfer are defined by you (not by coupled boundary condition of fluent, which is default for wall/shadow on boundary between fluid and solid, or which is true in interface)

as far as I understand the code from thesis works like this:
you need ID of WICK and ID of VAPOR faces, comes from fluent GUI, so you know it

you make a first loop over all finite faces in WICK face,
inside each of finite face you have other face loop over VAPOR finite faces

checking distance between WICK finite face and each VAPOR finite face you find which one is me closest, actually where the distance is less than 1.0E-6

if (fabs((xx_shadow-xx)/TOTAL_LENGHT) < 1.0E-6)

if it is true you write in the array the VAPOR face and index

ff_VAPOR_WICK [iii]=ff_shadow;
iii_VAPOR_WICK[jjj]=iii;

that is how you link WICK face with VAPOR face

in DEFINE_ADJUST macro you are taken parameters from cells adjust to each finite face of WICK and apply these parameters to cells adjusted to respective VAPOR faces (using the link, defined before)

I hope that is what you've asked for

[arunraj]

Thank you Alexander for you kind reply. If his method is right, why do I always get the SEGMENTATION FAULT, when I execute the below UDF. I have followed similar approach. I used wick_vapor and vapor_wick as the name of two-sided velocity inlet. However, I used in DEFINE_ON_DEMAND as my method is slightly different. Mine is not working in DEFINE_ADJUST too.
Can you please point out my mistakes as I have been trying so hard on this problem?

#include "udf.h"
#include "math.h"
#include "sg_udms.h"
#include "sg.h"
#include "stdio.h"
#include "mem.h"
#include "dpm.h"
#include "surf.h"

face_t f_wick_vapor[750];
real rho_vapor;

DEFINE_INIT(vapor_wick_shadow, d)
{
face_t f, f_shadow;
face_t f_vapor_wick[750];
real x[ND_ND];
real xx, xx_shadow;
int iii, jjj;
real total_length = 0.5;
real dist;

Thread *t = Lookup_Thread(d, 10); //* ID 10 vapor_wick_interface*//
Thread *t_shadow = Lookup_Thread(d, 11); //* ID 11 wick_vapor_interface*//
iii=-1;
begin_f_loop(f,t)
{
iii=iii+1;
F_CENTROID(x,f,t);
xx=x[0];
jjj=-1;
begin_f_loop(f_shadow,t_shadow)
{
jjj=jjj+1;
F_CENTROID(x,f_shadow,t_shadow);
xx_shadow=x[0];

dist = fabs((xx_shadow-xx)/total_length);

if (dist <= 1.0E-6)
{
f_vapor_wick[iii]=f_shadow;
}
}
end_f_loop(f_shadow,t_shadow)
}
end_f_loop(f,t)
}

DEFINE_INIT(wick_vapor_shadow, d)
{
face_t f, f_shadow;
real x[ND_ND];
real xx, xx_shadow;
int iii, jjj;
real total_length = 0.5;
real dist;

Thread *t = Lookup_Thread(d, 11); //* ID 11 wick_vapor_interface*//
Thread *t_shadow = Lookup_Thread(d, 10); //* ID 10 vapor_wick_interface*//
iii=-1;
begin_f_loop(f,t)
{
iii=iii+1;
F_CENTROID(x,f,t);
xx=x[0];
jjj=-1;
begin_f_loop(f_shadow,t_shadow)
{
jjj=jjj+1;
F_CENTROID(x,f_shadow,t_shadow);
xx_shadow=x[0];
dist = fabs((xx_shadow-xx)/total_length);

if (dist <= 1.0E-6)
{
f_wick_vapor[iii]=f_shadow;
}
}
end_f_loop(f_shadow,t_shadow)
}
end_f_loop(f,t)
}

DEFINE_ON_DEMAND(density_interface)
{
Domain *d;
d = Get_Domain(1);
face_t f, f_shadow;

cell_t c, c0, c0_shadow;
Thread *t0;
Thread *t0_shadow;

real ps, temp;
real R = 8314.0;
real M = 22.98977;
real po = 1596.0;
real to = 818.0;
real RR;

real x[ND_ND], A[ND_ND], dr0[ND_ND], es[ND_ND], dr0s[ND_ND], ess[ND_ND];
real ds, dss, A_by_Es, A_by_Ess, k_ds, k_dss;
real temp_cell, temp_face, temp_cell_shadow, temp_face_shadow, press_cell;
real k_vapor, k_liquid, k_eff_wick, den_vapor, den_liquid;
real k_wick = 66.73824;
real porosity = 0.33;
real hfg = 4217920.0;
real flux_vapor, flux_wick;

real v_interface_vapor[750], v_interface_wick[750];
real temp_interface[750];
real xx;
int iii;
int iii_vapor_wick[750];

Thread *td = Lookup_Thread(d, 4); /* ID 4 is vapor doamin */
Thread *t = Lookup_Thread(d, 10); //* ID 10 vapor_wick_interface*//
Thread *t_shadow = Lookup_Thread(d, 11); //* ID 11 wick_vapor_interface*//

begin_c_loop(c,td)
{
ps = C_P(c,td) + op_pres;
temp = C_UDSI(c,td,0);
RR=R/M;
C_UDMI(c,td,0) = ps/(RR*temp);
}
end_c_loop(c,td)

iii=-1;
begin_f_loop(f,t)
{
iii=iii+1;
F_CENTROID(x,f,t);
xx=x[0];
c0 = F_C0(f,t);
t0 = THREAD_T0(t);
f_shadow = f_wick_vapor[iii];
c0_shadow = F_C0(f_shadow,t_shadow);
t0_shadow = THREAD_T0(t_shadow);
k_vapor = C_K_L(c0,t0);
den_vapor = C_R(c0,t0);
temp_cell = C_UDSI(c0,t0,0);
temp_face = F_UDSI(f,t,0);
press_cell = C_P(c0,t0);
BOUNDARY_FACE_GEOMETRY(f,t,A,ds,es,A_by_Es,dr0);
k_ds = k_vapor/ds;
flux_vapor = k_ds*(temp_cell-temp_face);
k_liquid = C_K_L(c0_shadow,t0_shadow);
k_eff_wick = (porosity*k_liquid)+((1-porosity)*k_wick);
den_liquid = C_R(c0_shadow,t0_shadow);
temp_cell_shadow = C_UDSI(c0_shadow,t0_shadow,0);
temp_face_shadow = F_UDSI(f_shadow,t_shadow,0);
BOUNDARY_FACE_GEOMETRY(f_shadow,t_shadow,A,dss,ess ,A_by_Ess,dr0s);
k_dss = k_eff_wick/dss;
flux_wick = k_dss*(temp_cell_shadow-temp_face_shadow);
v_interface_vapor[iii] = (flux_vapor-flux_wick)/(hfg*C_UDMI(c,td,0));
v_interface_wick[iii] = v_interface_vapor[iii]*(C_UDMI(c,td,0)/den_liquid);
RR=R/M;
temp_interface[iii] = 1/((1/to)-((RR/hfg)*log(press_cell/po)));
}
end_f_loop(f,t)
}

[arunraj]

Is it possible to get the face index based on the below UDF.

#include "udf.h"
#include "surf.h"
#include "cxsurf.h"
#include "cxiface.h"

#define SID 14

DEFINE_ON_DEMAND(line_integrals)
{
Surface *s;
int i, j, k, index;
float v[2][3];
float face_val, line_int, line_ave, curve_length, tot_length;

/* s is the surface having data of surface id SID
s->np no of points.
s->nf no of faces.
s->points points to array of np points.
s->ip pointer to interpolation vector of np points.
s->facets facet list of the surface. length of array is nfl.
s->cells pointer to cells of each point.
s->nfl facet list length.
*/
s = SurfaceList+SID;

tot_length = 0.;
line_int = 0.;

for (i = 0, k = 0; k < s->np; i +=s->facets[i]+1, k++)
{
if (s->facets[i] != 2) continue; /* makes sure only line is there */
/*facet value for line surfaces will be the cell value in
which the facet lies*/
face_val = SV_Cell_Field(SVFIELD_CURRENT, s->cells+k);
for (j = 0; j < s->facets[i]; j++)
{
index = s->facets[i+j+1];
/* sets point coord x y z*/
SET_POINT_COORD(v[j][0], s->points[index], X_DIM);
SET_POINT_COORD(v[j][1], s->points[index], Y_DIM);
SET_POINT_COORD(v[j][2], s->points[index], Z_DIM);
}
/* Length calculation */
curve_length = sqrt(SQR(v[0][0]-v[1][0])+SQR(v[0][1]-v[1][1])+SQR(v[0][2]-v[1][2]));
tot_length += curve_length;
line_int += curve_length * face_val;
}

line_ave = line_int/tot_length;
Message("Line integral of currently displayed variable on line of ID %d : %6f\n", SID, line_int);
Message("Line length average of currently displayed variable on line of ID %d : %6f\n", SID, line_ave);
}

[AlexanderZ]

first of all, put d = Get_Domain(1); after all variable definitions, now it is on top (it should not work actually)

you are using user-defined scalars, are the defined?
temp = C_UDSI(c,td,0);

are you sure it is not zero in the whole domain? cause you are using this: C_UDMI(c,td,0) = ps/(RR*temp); which most likely is dividing by zero

always compile the code


I have no idea about second code, cause never've seen these commands before

[arunraj]

Hi Alexander, I modified the code as per your instruction. But SEGMENTATION FAULT error exist when i execute DEFINE_ON_DEMAND. I did not complete the entire set up. I even tried by replacing C_UDSI with C_T and turned ON energy, but still the same error exists. I have attached the steps below for UDSI case. I define 4 memory locations as well.

#include "udf.h"
#include "math.h"
#include "sg_udms.h"
#include "sg.h"
#include "stdio.h"
#include "mem.h"
#include "dpm.h"
#include "surf.h"

face_t f_wick_vapor[750];
real rho_vapor;

DEFINE_INIT(vapor_wick_shadow, d)
{
face_t f, f_shadow;
face_t f_vapor_wick[750];
real x[ND_ND];
real xx, xx_shadow;
int iii, jjj;
real total_length = 0.5;
real dist;

Thread *t = Lookup_Thread(d, 10); //* ID 10 vapor_wick_interface*//
Thread *t_shadow = Lookup_Thread(d, 11); //* ID 11 wick_vapor_interface*//
iii=-1;
begin_f_loop(f,t)
{
iii=iii+1;
F_CENTROID(x,f,t);
xx=x[0];
jjj=-1;
begin_f_loop(f_shadow,t_shadow)
{
jjj=jjj+1;
F_CENTROID(x,f_shadow,t_shadow);
xx_shadow=x[0];

dist = fabs((xx_shadow-xx)/total_length);

if (dist <= 1.0E-6)
{
f_vapor_wick[iii]=f_shadow;
}
}
end_f_loop(f_shadow,t_shadow)
}
end_f_loop(f,t)
}

DEFINE_INIT(wick_vapor_shadow, d)
{
face_t f, f_shadow;
real x[ND_ND];
real xx, xx_shadow;
int iii, jjj;
real total_length = 0.5;
real dist;

Thread *t = Lookup_Thread(d, 11); //* ID 11 wick_vapor_interface*//
Thread *t_shadow = Lookup_Thread(d, 10); //* ID 10 vapor_wick_interface*//
iii=-1;
begin_f_loop(f,t)
{
iii=iii+1;
F_CENTROID(x,f,t);
xx=x[0];
jjj=-1;
begin_f_loop(f_shadow,t_shadow)
{
jjj=jjj+1;
F_CENTROID(x,f_shadow,t_shadow);
xx_shadow=x[0];
dist = fabs((xx_shadow-xx)/total_length);

if (dist <= 1.0E-6)
{
f_wick_vapor[iii]=f_shadow;
}
}
end_f_loop(f_shadow,t_shadow)
}
end_f_loop(f,t)
}

DEFINE_ON_DEMAND(density_interface)
{
face_t f, f_shadow;
cell_t c, c0, c0_shadow;
Thread *t0;
Thread *t0_shadow;
real ps, temp;
real R = 8314.0;
real M = 22.98977;
real po = 1596.0;
real to = 818.0;
real RR;
real x[ND_ND], A[ND_ND], dr0[ND_ND], es[ND_ND], dr0s[ND_ND], ess[ND_ND];
real ds, dss, A_by_Es, A_by_Ess, k_ds, k_dss;
real temp_cell, temp_face, temp_cell_shadow, temp_face_shadow, press_cell;
real k_vapor, k_liquid, k_eff_wick, den_vapor, den_liquid;
real k_wick = 66.73824;
real porosity = 0.33;
real hfg = 4217920.0;
real flux_vapor, flux_wick;
real v_interface_vapor[750], v_interface_wick[750];
real temp_interface[750];
real xx;
int iii;
int iii_vapor_wick[750];
Domain *d;
d = Get_Domain(1);
Thread *td = Lookup_Thread(d, 4); /* ID 4 is vapor doamin */
Thread *t = Lookup_Thread(d, 10); //* ID 10 vapor_wick_interface*//
Thread *t_shadow = Lookup_Thread(d, 11); //* ID 11 wick_vapor_interface*//

begin_c_loop(c,td)
{
ps = C_P(c,td) + op_pres;
temp = C_UDSI(c,td,0);
RR=R/M;
C_UDMI(c,td,0) = ps/(RR*temp);
}
end_c_loop(c,td)

iii=-1;
begin_f_loop(f,t)
{
iii=iii+1;
F_CENTROID(x,f,t);
xx=x[0];
c0 = F_C0(f,t);
t0 = THREAD_T0(t);
f_shadow = f_wick_vapor[iii];
c0_shadow = F_C0(f_shadow,t_shadow);
t0_shadow = THREAD_T0(t_shadow);
k_vapor = C_K_L(c0,t0);
den_vapor = C_R(c0,t0);
temp_cell = C_UDSI(c0,t0,0);
temp_face = F_UDSI(f,t,0);
press_cell = C_P(c0,t0);
BOUNDARY_FACE_GEOMETRY(f,t,A,ds,es,A_by_Es,dr0);
k_ds = k_vapor/ds;
flux_vapor = k_ds*(temp_cell-temp_face);
k_liquid = C_K_L(c0_shadow,t0_shadow);
k_eff_wick = (porosity*k_liquid)+((1-porosity)*k_wick);
den_liquid = C_R(c0_shadow,t0_shadow);
temp_cell_shadow = C_UDSI(c0_shadow,t0_shadow,0);
temp_face_shadow = F_UDSI(f_shadow,t_shadow,0);
BOUNDARY_FACE_GEOMETRY(f_shadow,t_shadow,A,dss,ess ,A_by_Ess,dr0s);
k_dss = k_eff_wick/dss;
flux_wick = k_dss*(temp_cell_shadow-temp_face_shadow);
v_interface_vapor[iii] = (flux_vapor-flux_wick)/(hfg*C_UDMI(c,td,0));
v_interface_wick[iii] = v_interface_vapor[iii]*(C_UDMI(c,td,0)/den_liquid);
RR=R/M;
temp_interface[iii] = 1/((1/to)-((RR/hfg)*log(press_cell/po)));
}
end_f_loop(f,t)
}

[AlexanderZ]

did you initialize the case before execute on demand?

[arunraj]

I initialized, did patch and then execute DEFINE_ON_DEMAND, still the same error.

[AlexanderZ]

so it means the problem is in DEFINE_ON_DEMAND(density_interface)

v_interface_vapor[iii] = (flux_vapor-flux_wick)/(hfg*C_UDMI(c,td,0));
v_interface_wick[iii] = v_interface_vapor[iii]*(C_UDMI(c,td,0)/den_liquid);

C_UDMI(c,td,0) is not defined, cause td is not defined you are inside begin_f_loop(f,t) loop


and Message0 macro in your code to define which line leads to an error
most likely this line

 v_interface_vapor[iii] = (flux_vapor-flux_wick)/(hfg*C_UDMI(c,td,0));