[fonk]

Good Morning,

I wrote a UDF to compute the mass transfer occurring in a multi-phase flow. I was able to successfully compile, load and hook the function however, when doing the solution initialization (hybrid initialization), I get the normal message saying the initialization is done but then fluent keeps computing for unknown reasons(fluent processes using 100% of my CPU).

I managed to identify the problems as the commands PRF_GRSUM1 and PRF_GISUM1 using comments.

Basically I compiled:

variables declarations

cell loop that compute the sum of the thermal conductivities as well as cell number within the node's interior cells that presents only one phase.

Global sum of the thermal conductivities and cell number over all nodes.

(lot's of code set as comments)

return(m_to_v) (previously declared as equal to 0)

As a result, if PRF_GRSUM1 is not put as a comment fluent shut down after initialization(it is too quick for me to see if there is any error message before that). If PRF_GRSUM1 is put as comment and PRF_GISUM1 is not, fluent keeps computing after the "Hybrid Initialization is done" message. If I put all global sum as comment and "uncomment" the rest of the code, initialization ends normally and I can launch the calculation but it results in this error:

Node 9: Process 11636: Received signal SIGSEGV

MPI Application rank 5 exited before MPI_Finalize() with status 2


Then Fluent shutdown.

So I have 2 issues:

• infinite computation/shutting down with global sum
• SIGSEGV error

My code for the UDF
(You will see some questions in the code as comments. I used it to remind me that I wasn't sure about the code when I wrote it and that I need to check/look into it. You can disregard them or answer them as it pleases you.)
mass_transfer_2.c <- easier to read

#include "udf.h"

DEFINE_MASS_TRANSFER(mass_transfer, cell, mix_thread, from_index, from_species_index, to_index, to_species_index)

{ 
real summ_lambda = 0.0;
real summ_lambda_tot=0.0;
real mean_lambda_unsaturated = 0.0;
real summ_q=0.0;
real summ_q_cis=0.0;
real m_to_v=0.0;
int cell_compt=0;
int cell_compt_tot=0;
real L_heat = 1158.1; /* Heat of vaporization (kJ/kg) for saturation properties such as T = 300K (27C) and P = 1.0611 MPa */
real Vol=1.0;
real Vol_cis=1.0;


#if !RP_HOST /*careful about which variables are also used in the host*/
		
	real q_c_i_in_int_2=0.0;
	real q_c_i_in_int_2_cis=0.0;
	real NV_VEC(q_c_i_in_int_1);
	real NV_VEC(q_c_i_in_int_1_cis);
	real NV_VEC(gradT);
	real NV_VEC(gradT_cis);
	real NV_VEC(area_vec);
	real NV_VEC(area_vec_cis);
	real threshold_two_phase=0.000001;
	real threshold_one_phase=0.999999;
	real NV_VEC(A_cs);
	real NV_VEC(A_cis);
	real NV_VEC(pos);
	real NV_VEC(pos_cis);
	real NV_VEC(ex);
	real NV_VEC(ey);
	real NV_VEC(ez);
	real c0_centroid[ND_ND];
	real c0_cis_centroid[ND_ND];
	real c1_centroid[ND_ND];
	real c1_cis_centroid[ND_ND];
	face_t f;
	face_t f_cis;
	cell_t c0, c0_cis;
	cell_t c1, c1_cis;
	Thread *gas,*liq, *ft, *ft_cis;
	int n;
	int m;
	liq = THREAD_SUB_THREAD(mix_thread, from_index); 	/* From phase interaction in Fluent */ 
	gas = THREAD_SUB_THREAD(mix_thread, to_index);		/* From phase interaction in Fluent */ 
	/* Compute the mean thermal conductivity of the gas phase (unsaturated pahse during the boiling process) */
	/*Do I need to use the thread mix_thread or the gas thread ?*/
        begin_c_loop_int(cell,gas) /* Is it performing a loop on every cell in the gas phase? Does it also include cells in which the two phases exists? I need gas phase only*/
	{
	if (C_VOF(cell,gas) > threshold_one_phase)        /*because I need gas only*/
	{
	summ_lambda += C_K_L(cell,gas);  
	cell_compt+=1;
	}
        }
	end_c_loop_int(cell, gas)
#endif
	
#if RP_NODE /*Nothing in serial as this is used to sum the data of all compute node*/	
	summ_lambda_tot = PRF_GRSUM1(summ_lambda);
	cell_compt_tot = PRF_GISUM1(cell_compt);
#endif

node_to_host_real_1(summ_lambda_tot); 	/* Pass the information from node 0 to host*/
node_to_host_int_1(cell_compt_tot);		/* Pass the information from node 0 to host*/
	
#if !RP_NODE /* Calculate the mean thermal conductivity if host or serial*/
	mean_lambda_unsaturated = summ_lambda_tot/cell_compt_tot;
	Message("mean thermal conductivity of the unsaturated phase (gas) is %f (W/m.K) \n", mean_lambda_unsaturated);
#endif

host_to_node_real_1(mean_lambda_unsaturated); /*Give the result of the global mean to the compute nodes*/


#if !RP_HOST


	begin_c_loop_int(cell, mix_thread) /*loop needs to be on internal cells only so that it is possible to make calculation with properties from regular exteriors cells  */
	/* Using a loop on internal + regular exterior cells will results in doing calcualtion twice in case the interface is at a boundary. */ 
	{
	if (C_VOF(cell,gas) > threshold_two_phase && C_VOF(cell,gas) < threshold_one_phase)        /*If it is an interface (gas present but also liquid) then...*/
	{
	/*...performs calculation of the imported heat */

	Vol += C_VOLUME(cell, mix_thread);
	c_face_loop(cell, mix_thread, n)        /* loop over faces of cell, do n change from 0 to number of face cell automatically? */

	{/*Here it is not needed to know if it is principal face or not */

	f = C_FACE(cell, mix_thread, n);        /* get face n of cell */
	ft = C_FACE_THREAD(cell, mix_thread, n);      /* get face thread of face f */
	c0 = F_C0(f, ft);        /* get c0, cell at one side of the face f */
	c1 = F_C1(f, ft);        /* get c1, cell at the other side of face f */
	F_AREA(A_cs, f, ft);       /* get the vector area of the face (magnitude = surface area),direction and sign correspond to the normal vector that point out of the domain */
	NV_VS(area_vec, = , A_cs , / , -NV_MAG(A_cs));	       /* compute the normal vector of the face using the vector area */
	C_CENTROID(c0_centroid, c0, mix_thread);	       /* get the centroid of cell c0 */
	C_CENTROID(c1_centroid, c1, mix_thread);        /* get the centroid of cell c1 */
	NV_VV(pos, =, c1_centroid, -, c0_centroid);       /* get the position vector for gradient calculation */
	NV_VS(gradT, =, pos, *, ((C_T(c1, mix_thread) - C_T(c0, mix_thread))/(NV_MAG2(pos))));	/* gradient calculation */
	NV_VS(q_c_i_in_int_1, =, gradT, *, (-1)*mean_lambda_unsaturated*NV_MAG(A_cs));
	q_c_i_in_int_2 = NV_DOT(q_c_i_in_int_1,area_vec);
					
	if (q_c_i_in_int_2<0)
	{
	q_c_i_in_int_2 = 0;
	}
	summ_q += q_c_i_in_int_2;/* summ of the heat towards the gas phase */
					

	/*summ_q += (NV_DOT(NV_VS(q_c_i_in_int_1, =, gradT, *, (-1)*mean_lambda_unsaturated*NV_MAG(A_cs)) , area_vec),0);*/
	if(summ_q<0)
	{
	Message("WARNING imported heat is negative %f\n",summ_q);
	}


	if (cell == c0) /*if cell c0 is the cell with the interface then... */
	{
	if (C_VOF(c1,liq) > threshold_one_phase) /* ...look if c1 is liquid phase (saturated phase) ... */
	{
	Vol_cis += C_VOLUME(c1, mix_thread); /*... and make calculation of imported heat if it is the case*/
	c_face_loop(c1, mix_thread, m)
	{
	f_cis = C_FACE(c1, mix_thread, m); 
	ft_cis = C_FACE_THREAD(c1, mix_thread, m); 
	F_AREA(A_cis,f_cis,ft_cis);
	c0_cis = F_C0(f_cis,ft_cis); 
	c1_cis = F_C1(f_cis,ft_cis); 
	NV_VS(area_vec_cis, = , A_cis , / , -NV_MAG(A_cis)); 
	C_CENTROID(c0_cis_centroid, c0_cis,mix_thread);	
	C_CENTROID(c1_cis_centroid, c1_cis,mix_thread); 
	NV_VV(pos_cis, =, c1_cis_centroid, -, c0_cis_centroid); 
	NV_VS(gradT_cis, =, pos_cis, *, ((C_T(c1_cis, mix_thread) - C_T(c0_cis,mix_thread))/NV_MAG2(pos_cis))); 
	NV_VS(q_c_i_in_int_1_cis, =, gradT_cis, *, (-1)*mean_lambda_unsaturated*NV_MAG(A_cis));
	q_c_i_in_int_2_cis = NV_DOT(q_c_i_in_int_1_cis,area_vec_cis);
							
	if (q_c_i_in_int_2_cis < 0)
	{
	q_c_i_in_int_2_cis = 0;
	}

	summ_q_cis += q_c_i_in_int_2_cis;
	/*summ_q_cis += (NV_DOT(NV_VS(q_c_i_in_int_1_cis, = , gradT_cis, *, -mean_lambda_unsaturated*NV_MAG(A_cis)), -area_vec_cis),0);*/
	if (summ_q_cis<0)
	{
	Message("WARNING imported heat is negative %f\n",summ_q_cis);
	}
        }
	}
	}
	else
	{
	if (C_VOF(c0,liq) > threshold_one_phase)
	{
	Vol_cis += C_VOLUME(c0, mix_thread);
	c_face_loop(c0, mix_thread, m)
	{
	f_cis = C_FACE(c0, mix_thread, m);
	ft_cis = C_FACE_THREAD(c0, mix_thread, m);
	F_AREA(A_cis, f_cis, ft_cis);
	c0_cis = F_C0(f_cis,ft_cis); 
	c1_cis = F_C1(f_cis,ft_cis); 
	NV_VS(area_vec_cis, = , A_cis , / , NV_MAG(A_cis)); 
	C_CENTROID(c0_cis_centroid, c0_cis,mix_thread);	
	C_CENTROID(c1_cis_centroid, c1_cis,mix_thread); 
	NV_VV(pos_cis, =, c1_cis_centroid, -, c0_cis_centroid);
	NV_VS(gradT_cis, =, pos_cis, *, ((C_T(c1_cis, mix_thread) - C_T(c0_cis,mix_thread))/NV_MAG2(pos_cis)));
	NV_VS(q_c_i_in_int_1_cis, =, gradT_cis, *, (-1)*mean_lambda_unsaturated*NV_MAG(A_cis));
	q_c_i_in_int_2_cis = NV_DOT(q_c_i_in_int_1_cis,area_vec_cis);
								
	if (q_c_i_in_int_2_cis < 0)
	{
	q_c_i_in_int_2_cis = 0;
	}

	summ_q_cis += q_c_i_in_int_2_cis;
	/*summ_q_cis += (NV_DOT(NV_VS(q_c_i_in_int_1_cis, = , gradT_cis, *, -mean_lambda_unsaturated*NV_MAG(A_cis)), -area_vec_cis),0);*/	
	if (summ_q_cis<0)
	{
	Message("WARNING imported heat is negative %f\n ",summ_q_cis);
	}
        }	
	}
	}
	}
	}
	}
	end_c_loop_int(cell, mix_thread)
#endif

#if RP_NODE /*in case of parallel calulations, global sum of the imported heat over all nodes*/
	summ_q=PRF_GRSUM1(summ_q);
	summ_q_cis=PRF_GRSUM1(summ_q_cis);
	Vol = PRF_GRSUM1(Vol);
	Vol_cis = PRF_GRSUM1(Vol_cis);
#endif

	/* Send data to host for calculation */

node_to_host_real_2(summ_q,summ_q_cis);

#if !RP_NODE
	m_to_v =  summ_q/(L_heat*Vol) + summ_q_cis/(L_heat*Vol_cis);  
	Message("calculated mass transfer flux is %f (kg/(m^3.s)) \n",m_to_v);
	return(m_to_v);
#endif
}

I use Ansys 19.0
The system is one flat plate of heat exchanger with liquid water on one side and liquid ammonia to be evaporated on the other side.

I use implicit VOF as multiphase model with sharp interface. 3 Phases: vapor ammonia, liquid ammonia and water.

Viscous model: k-epsilon, Realizable, Enhanced Wall Treatment

Pressure-velocity scheme coupled + coupled with Volume Fractions
Spatial discretization: Gradient: :Least Squares Cell based
Pressure: PRESTO!
Momentum, Turbulent Kinetic Energy, Turbulent dissipation rate and Energy: First order upwind (I tried QUICK and second order upwind but it resulted in a floating point exception).

The model works without mass transfer. It also works with the Evaporation-condensation model (500 iterations was not enough to meet convergence criteria but no error occurred. I did not try more than 500 iterations as I just wanted to check for errors.)

Mesh picture is attached.
Capture.jpg
Nodes: 93636
Elements: 87120

Your ideas would be appreciated.
Thank you for your help.

[vinerm]

There are a few issues with the code. I still haven't gone through whole of it but the ones I identified are as follows.

1. One of the comment is not closed with */ in the beginning. Since /* */ can span multiple lines, you do not need to use it on each line. Start it at the beginning and close it after all the initial comments have been written.

2. You need to transfer cell_compt to node before using it since those are defined before !RP_HOST, i.e., in the host. Actually, the whole code does not have anything that requires host. So, you can put whole of the code within !RP_HOST.

3. Use volume based averaging instead of arithmatic. So, prefer C_K_L * C_VOLUME and then divide by total volume. Not an error but improvement.

4. Your values for cut-off fractions are very strict. Use 0.9 and 0.1 instead of 0.999999 and 0.000001. Numerically, even 1 is represented as 0.999999. And even with 0.9 and 0.1, liquid and gas will be present only near the interface. So, you will get correct results.

5. You are determining face indices in cell loop and then fetching same cell information from each face. This is making it go into a very big loop. To fetch gradients, simply use C_T_G. I could not understand all the equations that you are using but temperature gradient can be fetched using C_T_G.

I have not checked if the phase indices are correctly being used or not. Those could also be reason for infinite loop. Another addition you can do is to call PRF_GSYNC() command before using global sums and reductions. This ensure that summing and reduction is done only after all nodes have completed their work.

[fonk]

Good morning,

Thank you for your answer Vinerm,

1)Regarding the remark about the multiple /**/, I fixed it

2)

You need to transfer cell_compt to node before using it since those are defined before !RP_HOST, i.e., in the host.

From my understanding of the Fluent customization manual, commands that are not wrap in any #if (!)RP_NODE/HOST/PARALLEL are computed by both host and nodes. cf 7.3.5 Parallel UDF Example:

#include "udf.h"
 
 DEFINE_ADJUST(face_av,domain)
 {
 /* Variables used by serial, host, node versions */
 int surface_thread_id=0;
 real total_area=0.0;
 real total_force=0.0;
 
 /* "Parallelized" Sections */
 #if !RP_HOST  /* Compile this section for computing processes only (serial
          and node) since these variables are not available
          on the host */
     Thread* thread;
     face_t face;
     real area[ND_ND];
 #endif /* !RP_HOST */
...

Please rectify me if I am wrong. I did try to put the definition within #if !RP_HOST but the problem remained.

Actually, the whole code does not have anything that requires host. So, you can put whole of the code within !RP_HOST.

I use global reduction macro to and then pass the result to the host for further computing (within the host). one of the part of the code I am referring to:

	
#if RP_NODE /*Nothing in serial as this is used to sum the data of all compute node*/	
summ_lambda_tot = PRF_GRSUM1(summ_lambda);
vol_compt_tot = PRF_GRSUM1(vol_compt);
#endif
/* Pass the information from node 0 to host*/
node_to_host_real_2(summ_lambda_tot,vol_compt_tot); 	
		
#if !RP_NODE /* Calculate the mean thermal conductivity if host or serial*/
mean_lambda_unsaturated = summ_lambda_tot/vol_compt_tot;
Message("mean thermal conductivity of the unsaturated phase (gas) is %f (W/m.K) \n", mean_lambda_unsaturated);
#endif
/*Give the result of the global mean to the compute nodes*/
host_to_node_real_1(mean_lambda_unsaturated);

Again, please tell me if I am wrong.

3. Use volume based averaging instead of arithmatic. So, prefer C_K_L * C_VOLUME and then divide by total volume. Not an error but improvement.

4. Your values for cut-off fractions are very strict. Use 0.9 and 0.1 instead of 0.999999 and 0.000001. Numerically, even 1 is represented as 0.999999. And even with 0.9 and 0.1, liquid and gas will be present only near the interface. So, you will get correct results.

I followed your recommendation. Thank you.

5. You are determining face indices in cell loop and then fetching same cell information from each face. This is making it go into a very big loop. To fetch gradients, simply use C_T_G. I could not understand all the equations that you are using but temperature gradient can be fetched using C_T_G.

I am aware of the command C_T_G however this command gives the the gradient in the cell, not the gradient at the cell face which I need. Unfortunately, the command F_T_G does not exist.

I realized my first post was probably trying to address too many issues. Let's focus on the infinite computation. As I said before, I put most of my calculation as comment to try solving the issue. Here are the details:

This code is causing the infinite computation after the message "Hybrid Initialization is done":

#include "udf.h"
DEFINE_MASS_TRANSFER(mass_transfer, cell, mix_thread, from_index, from_species_index, to_index, to_species_index)

{   
real m_to_v=0.0;
real summ_lambda = 0.0;
real summ_lambda_tot=0.0;
real mean_lambda_unsaturated = 0.0;
real summ_q=0.0;
real summ_q_cis=0.0;
real vol_compt=0.0;
real vol_compt_tot=0.0;
real L_heat = 1158.1;       
real Vol=0.0;
real Vol_cis=0.0;
#if !RP_HOST /*Following variables are only used in nodes or in serial*/
real q_c_i_in_int_2=0.0;
real q_c_i_in_int_2_cis=0.0;
real NV_VEC(q_c_i_in_int_1);
real NV_VEC(q_c_i_in_int_1_cis);
real NV_VEC(gradT);
real NV_VEC(gradT_cis);
real NV_VEC(area_vec);
real NV_VEC(area_vec_cis);
real threshold_two_phase=0.1;
real threshold_one_phase=0.9;
real NV_VEC(A_cs);
real NV_VEC(A_cis);
real NV_VEC(pos);
real NV_VEC(pos_cis);
real c0_centroid[ND_ND];
real c0_cis_centroid[ND_ND];
real c1_centroid[ND_ND];
real c1_cis_centroid[ND_ND];
face_t f;
face_t f_cis;
cell_t c0, c0_cis;
cell_t c1, c1_cis;
Thread *gas,*liq, *ft, *ft_cis;
int n;
int m;

liq = THREAD_SUB_THREAD(mix_thread, from_index);
gas = THREAD_SUB_THREAD(mix_thread, to_index);
/* Compute the mean thermal conductivity of the gas phase (unsaturated
phase during the boiling process) 
begin_c_loop_int(cell,gas) 
{
if (C_VOF(cell,gas) > threshold_one_phase) /*because I need gas only*/
{
summ_lambda += C_K_L(cell,gas)*C_VOLUME(cell,gas);  
vol_compt+=C_VOLUME(cell,gas);
}
}
end_c_loop_int(cell, gas)
#endif
/*Without the 5 following lines , the code works normally, computation stops
after the initialization and then I am able to launch the calculation until
convergence. Obviously it does not add anything to the calculation as it
always returns 0.0  */
#if RP_NODE 
PRF_GSYNC();    
summ_lambda_tot = PRF_GRSUM1(summ_lambda);
vol_compt_tot = PRF_GRSUM1(vol_compt);
#endif
/*end of the five faulty lines, problem is the same with or without
PRF_GSYNC, it is also the same with only PRF_GSYNC*/
#if !RP_NODE
return(m_to_v);
#endif

I also tried to pass the variables to the node as recommended by Vinerm but result is the same.
I hope this give a better understanding of the issue.
I added a new version of the complete code if needed, please keep in mind that for the test I mentioned above, lines 68 to 173 as well as lines 175 and 176 were not active.
mass_transfer_2.c

Again, Thank you for your help.

[vinerm]

Modify the code as follows

#if !RP_HOST /* Keep its location wherever it is */

#if RP_NODE
PRF commands
#endif

#endif /* End !RP_HOST here, after including PRF commands */

[fonk]

Good Morning,

Thank you for the fast answer.
I tried your solution, without success unfortunately. This is the code I tried:

#include "udf.h"
DEFINE_MASS_TRANSFER(mass_transfer, cell, mix_thread, from_index, from_species_index, to_index, to_species_index)

{   
real m_to_v=0.0;
real summ_lambda = 0.0;
real summ_lambda_tot=0.0;
real mean_lambda_unsaturated = 0.0;
real summ_q=0.0;
real summ_q_cis=0.0;
real vol_compt=0.0;
real vol_compt_tot=0.0;
real L_heat = 1158.1;			
real Vol=0.0;
real Vol_cis=0.0;
#if !RP_HOST 
real q_c_i_in_int_2=0.0;
real q_c_i_in_int_2_cis=0.0;
real NV_VEC(q_c_i_in_int_1);
real NV_VEC(q_c_i_in_int_1_cis);
real NV_VEC(gradT);
real NV_VEC(gradT_cis);
real NV_VEC(area_vec);
real NV_VEC(area_vec_cis);
real threshold_two_phase=0.1;
real threshold_one_phase=0.9;
real NV_VEC(A_cs);
real NV_VEC(A_cis);
real NV_VEC(pos);
real NV_VEC(pos_cis);
real c0_centroid[ND_ND];
real c0_cis_centroid[ND_ND];
real c1_centroid[ND_ND];
real c1_cis_centroid[ND_ND];
face_t f;
face_t f_cis;
cell_t c0, c0_cis;
cell_t c1, c1_cis;
Thread *gas,*liq, *ft, *ft_cis;
int n;
int m;
liq = THREAD_SUB_THREAD(mix_thread, from_index); 	 
gas = THREAD_SUB_THREAD(mix_thread, to_index);		
begin_c_loop_int(cell,gas) 
{
if (C_VOF(cell,gas) > threshold_one_phase) 
{
summ_lambda += C_K_L(cell,gas)*C_VOLUME(cell,gas);  
vol_compt+=C_VOLUME(cell,gas);
}
}		
end_c_loop_int(cell, gas)
#if RP_NODE 
summ_lambda_tot = PRF_GRSUM1(summ_lambda);
vol_compt_tot = PRF_GRSUM1(vol_compt);
#endif /* RP_NODE*/
#endif /* !RP_HOST*/
#if !RP_NODE  
Message("calculated mass transfer flux is %f (kg/(m^3.s)) \n",m_to_v);
return(m_to_v);
#endif
}

I also used PRF_GSYNC(); just before PRF_GRSUM1.
As a result, I get infinite computation after initialization with PRF_GSYNC(); and Fluent shutdown after initialization with PRF_GRSUM1 only.

I managed to get the error message before Fluent shutdown:

C:\PROGRA~1\ANSYSI~1\v190\fluent\fluent19.0.0\win6 4\3ddp_node\fl_mpil900.exe: Rank 0:8: MPI_Allreduce: Message truncated
MPI Application rank 8 exited before MPI_Finalize() with status 14
The fl process could not be started

I could not find any information with this specific status.
I hope this will help identify the issue.

Thank you for the help.
Kevin.

[vinerm]

AllReduce is essentially the function that is at the background of PRF_ functionality. Does this work in Serial? Check it in Serial because it appears the issue is with the compiler?

[fonk]

Good morning.

Thank you for your answer.

I launched the initialization with the full code with Fluent in serial (parallel 0 processes).

It worked although not as I thought it would. It seems that DEFINE_MASS_TRANSFER is called by Fluent within a cell loop. My code therefore includes unnecessary loops inside an already existing loop.

Indeed messages I programmed to get only once per iteration kept popped up after initialization(it did stop after a while but there is enough to fill up all the console window as high as I can scroll).

I guess it might be causing the issue with the PRF commands.

My calculations requires the mean thermal conductivity of one phase so I will have to re-think the code. I will probably set the thermal conductivity as a global variable and use a DEFINE_ADJUST to calculate the value before DEFINE_MASS_TRANSFER occurs.

Thank you again.
Kevin.

[vinerm]

Yes, there are many functions, such _PROPERTIES, _MASS_TRANSFER, _SOURCE, etc. that are executed per cell. So, you don't need to have a cell loop.

[AlexanderZ]

since you don't read/write to file or use rp_variables
for code above you don't need HOST at all, remove it

you need only this

#if RP_NODE 
summ_lambda_tot = PRF_GRSUM1(summ_lambda);
vol_compt_tot = PRF_GRSUM1(vol_compt);
#endif /* RP_NODE*/

to print these values use Message0(""); - information from node 0

[fonk]

Good morning,

Thank you for the help you both provided.

As I said previously I wrote a DEFINE_ADJUST macro to compute the mean thermal conductivity I needed as a global variable. This works without issue.

#include "udf.h"
#include "math.h"
	/*mean_lambda_unsaturated declared as global variable*/
real mean_lambda_unsaturated;
real threshold_one_phase  = 0.9;
DEFINE_ADJUST(volume_average_thermal_conductivity,mix_domain)
{
real vol_compt=0.0;
real summ_lambda=0.0;
#if !RP_HOST
real vof_liq;
real vof_gas;
cell_t cell;
Thread *t = Lookup_Thread(mix_domain,8); /*zone_ID given in boundary condition for interior-fluid_ammonia*/
Domain *gas_domain=Get_Domain(2); /*Assume gas domain ID is 2, check in phases menu in Fluent. Careful: 3 is liquid NH3 */
Domain *liq_domain=Get_Domain(3);
int gas_phase_domain_index=PHASE_DOMAIN_INDEX(gas_domain);
int liq_phase_domain_index=PHASE_DOMAIN_INDEX(liq_domain);
Thread *gas_thread = THREAD_SUB_THREAD(t,gas_phase_domain_index);
Thread *liq_thread = THREAD_SUB_THREAD(t,liq_phase_domain_index);

begin_c_loop_int(cell,t)
{
if (C_VOF(cell,gas_thread) > threshold_one_phase)
{
summ_lambda += C_K_L(cell,t)*C_VOLUME(cell,t);
vol_compt += C_VOLUME(cell,t);
}
}
end_c_loop_int(cell,t)
#if RP_NODE
summ_lambda=PRF_GRSUM1(summ_lambda);
vol_compt=PRF_GRSUM1(vol_compt);
#endif
#endif
node_to_host_real_2(summ_lambda,vol_compt);
#if !RP_NODE
if (vol_compt != 0)
{
mean_lambda_unsaturated = summ_lambda/vol_compt;	
}
else
{
mean_lambda_unsaturated = 28.343e-3;
}
#endif
host_to_node_real_1(mean_lambda_unsaturated);
}

The macro DEFINE_MASS_TRANSFER that follows works when using serial mode but not parallel mode. I get the following error then Fluent shutdown.

Node 9: Process 9104: Received signal SIGSEGV.

999999: mpt_accept: error: accept failed: No such file or directory (I got the same message 11 times followed by)

MPI_Application rank 0 exited before MPI_Finalize() with status 2
the fl process could not be started.

I could identify the line that cause this issue as the second c_face_loop macro.

Here is the code that I used to try to identify the issue:

DEFINE_MASS_TRANSFER(evaporation_mass_transfer, cell, mixture_thread, from_phase_index, from_species_index, to_phase_index, to_species_index)
{
/*DEFINE_MASS_TRANSFER is called in a cell loop*/
#if !RP_HOST
int n;
face_t f, f_cis;
Thread *gas_thread, *liq_thread, *ft, *ft_cis;
cell_t c0, c1, c_calc, c0_cis, c1_cis;
real c0_centroid[ND_ND];
real c0_cis_centroid[ND_ND];
real c1_centroid[ND_ND];
real c1_cis_centroid[ND_ND];
real NV_VEC(pos);
real NV_VEC(pos_cis);
real NV_VEC(A_cs);
real NV_VEC(A_cis);
real NV_VEC(area_vec);
real NV_VEC(area_vec_cis);
real NV_VEC(gradT);
real NV_VEC(gradT_cis);
real vol_cis;
real vol_compt_cis;
real q_c_i_in;
real NV_VEC(q_c_i_in_int_1);
real q_c_i_in_int_2;
real q_c_i_in_cis;
real NV_VEC(q_c_i_in_int_1_cis);
real q_c_i_in_int_2_cis;
real q_c_i_in_int_3_cis;
real L_heat = 1158.1e3;
real l_to_g;
int m;
liq_thread=THREAD_SUB_THREAD(mixture_thread,from_phase_index);
gas_thread=THREAD_SUB_THREAD(mixture_thread,to_phase_index);
/*reset l_to_g, q_c_i_in (imported heat in the cell at the interface),
q_c_i_in_cis(sum of the imported heat in all the neighbouring cells) */ 
l_to_g=0.0;
q_c_i_in=0.0;
q_c_i_in_cis=0.0;
vol_compt_cis=0.0;
/*Check if interface is present*/
if (C_VOF(cell,liq_thread)<threshold_one_phase && C_VOF(cell,gas_thread)<threshold_one_phase) 
{
c_face_loop(cell,mixture_thread,n)
{/*Calculation of the imported heat in the cells at the interface gas-liq*/
f = C_FACE(cell, mixture_thread, n);	/* get face n of cell */
ft = C_FACE_THREAD(cell, mixture_thread, n);	/* get face thread of face f */
NV_VS(q_c_i_in_int_1 , = , q_c_i_in_int_1 , * ,0);
q_c_i_in_int_2=0.0;
if PRINCIPAL_FACE_P(f,ft)
{
c0 = F_C0(f, ft);   /* get c0, cell at one side of the face f */
c1 = F_C1(f, ft);   /* get c1, cell at the other side of face f */
if (c1 != c0)
{
F_AREA(A_cs, f, ft);	/* get the vector area of the face (magnitude =
surface area),direction and sign correspond to the normal vector that point
out of the domain */
NV_VS(area_vec, = , A_cs , / , -NV_MAG(A_cs));	/* compute the normal vector of the face using the vector area */
C_CENTROID(c0_centroid, c0, mixture_thread);	/* get the centroid of cell c0 */
C_CENTROID(c1_centroid, c1, mixture_thread);	/* get the centroid of cell c1 */
NV_VV(pos, =, c1_centroid, -, c0_centroid);    /* get the position vector for gradient calculation */
/* calculation of the temperature gradient at the cell face */
NV_VS(gradT, =, pos, *, ((C_T(c1, mixture_thread) - C_T(c0, mixture_thread))/(NV_MAG2(pos))));	
NV_VS(q_c_i_in_int_1, =, gradT, *, (-1)*mean_lambda_unsaturated*NV_MAG(A_cs));  
q_c_i_in_int_2 = NV_DOT(q_c_i_in_int_1,area_vec); 

if (q_c_i_in_int_2 < 0)
{
q_c_i_in_int_2 = 0;
}

q_c_i_in += q_c_i_in_int_2; /*+= to account for the heat coming from all faces*/

/*Calculation of the imported heat in neighbour cells of the saturated phase
(liquid) 
set c_calc as the neighbour cell*/
if (cell == c0)
{
c_calc=c1;
}
else
{
c_calc=c0;
}
q_c_i_in_int_3_cis=0.0;/*reset q_c_i_in_int_3_cis: q_in the neighbouring
cell corresponding to face n*/
if (C_VOF(c_calc,liq_thread) >= threshold_one_phase) /*Check if c_calc is liquid phase*/
{/*then make calculations*/
vol_compt_cis += C_VOLUME(c_calc,mixture_thread);
c_face_loop(c_calc,mixture_thread,m)
{
Message("c_face_loop\n");
}/*end face loop on c_calc*/
}/*enf if C_VOF(c_calc)*/
}/*end if c1!=c0*/
}/*end if principal face(main cell)*/
}/*end face loop on main cell*/
l_to_g=q_c_i_in/(L_heat*C_VOLUME(cell,mixture_thread));
}/*end if main cell is interface*/
return(l_to_g);
#endif
}/*end of function go to neext cell*/

If I remove the second face loop the code works.

Thank you,
Kevin.

PS: I join the full code if needed:adjust.c

[vinerm]

If the problem is with the second loop, then you should start debugging one line at a time. You are also using -NV_MAG. Better would be to write it as (-1)*NV_MAG since NV_MAG is a macro. -NV_MAG may or may not work. But I doubt if that is the reason for segmentation violation. So, best would be to debug one line at a time.

Why do you have a nested c_face_loop? It is not that it is not allowed, but what could be the purpose of this. Since DEFINE_MASS_TRANSFER goes over each cell, and it appears that the code might be going to adjacent cells due to dual loop, the data might be overwritten when DEFINE macro goes over the adjacent cell and then modifies the data in the cell currently being addressed by the macro. I think you can avoid double loops. Furthermore, c1 is never equal to c0. c1, if exists, is always different from c0. However, there are situations when c1 does not exist. So, you should check if c1 != NULL before using it. This could lead to segmentation violation.

[fonk]

Dear Vinerm,

Thank you very much. The problem was caused, as you said, by Fluent attempting to read c1 when it is not allocated. I could solve it by adding if (!BOUNDARY_FACE_THREAD_P(face_thread)) before the calculation. (I could not use c1 != NULL because non allocation is different from void).

Regarding the dual loop:
I compute the mass transfer of the cell n as a function of both the heat imported in the cell n Qn,in and the heat imported in the cell next to it if it is a cell with only liquid Qn,adj,in.

Here is a picture:
fluent_grid.PNG
where the arrows represents the imported heat Q. I need to access the temperature at the cell adjacent to the adjacent cell for those calculation which explain the dual loop. (Q is computed as a function of the temperature gradient at a face)

However the dual loop did not appear to cause any issue. I had no error during the calculation. Now I have to figure out what is wrong with the model as I got energy residuals between 1e20 and 1e15 but I will post my problem on the Fluent Multiphase forum if I cannot figure out myself.

Again, thank you for helping me debugging the code.
Kevin.

[vinerm]

The problem with dual loop is not its running; it will run. However, it may lead to non-conservation of energy and might be leading to divergence of energy equation as you mentioned. It's not that you should not use dual loop, however, you may require some extra checks and bounds to ensure that energy conservation is not violated. I haven't gone through the details of the code (it appeared too involved to me without the availability of the mathematical equations) but you should check for bugs in the dual loop; the bugs won't be syntactical but logical.