CFD Online Discussion Forums

CFD Online Discussion Forums (https://www.cfd-online.com/Forums/)
-   Fluent UDF and Scheme Programming (https://www.cfd-online.com/Forums/fluent-udf/)
-   -   enquire explanation of the code (https://www.cfd-online.com/Forums/fluent-udf/161793-enquire-explanation-code.html)

mimi0201 October 28, 2015 20:59

enquire explanation of the code
 
Here is an example from ANSYS help, tracking particle's reflections at walls.

Code:

/* reflect boundary condition for inert particles */
 #include "udf.h"
 DEFINE_DPM_BC(bc_reflect,p,t,f,f_normal,dim)
 {
    real alpha; /* angle of particle path with face normal */
    real vn=0.;
    real nor_coeff = 1.;
    real tan_coeff = 0.3;
    real normal[3];
    int i, idim = dim;
    real NV_VEC(x);

for (i=0; i<idim; i++)
      normal[i] = f_normal[i];

    if(p->type==DPM_TYPE_INERT)
      {
        alpha = M_PI/2. - acos(MAX(-1.,MIN(1.,NV_DOT(normal,P_VEL(p))/
              MAX(NV_MAG(P_VEL(p)),DPM_SMALL))));
        if ((NNULLP(t)) && (THREAD_TYPE(t) == THREAD_F_WALL))
              F_CENTROID(x,f,t);

        /* calculate the normal component, rescale its magnitude by
          the coefficient of restitution and subtract the change */

        /* Compute normal velocity. */
        for(i=0; i<idim; i++)
          vn += P_VEL(p)[i]*normal[i];

        /* Subtract off normal velocity. */
          for(i=0; i<idim; i++)
            P_VEL(p)[i] -= vn*normal[i];

        /* Apply tangential coefficient of restitution. */
          for(i=0; i<idim; i++)
            P_VEL(p)[i] *= tan_coeff;

        /* Add reflected normal velocity. */
          for(i=0; i<idim; i++)
            P_VEL(p)[i] -= nor_coeff*vn*normal[i];

        /* Store new velocity in P_VEL0 of particle */
        for(i=0; i<idim; i++)
          P_VEL0(p)[i] = P_VEL(p)[i];

        return PATH_ACTIVE;
      }
  return PATH_ABORT;
 }

In the code, what exactly does the statement do as shown below?

Code:

for (i=0; i<idim; i++)
      normal[i] = f_normal[i];

and

if ((NNULLP(t)) && (THREAD_TYPE(t) == THREAD_F_WALL))
              F_CENTROID(x,f,t);


upeksa October 29, 2015 04:10

for (i=0; i<idim; i++)
normal[i] = f_normal[i];


int idim=dim, where dim is argument in DEFINE_DPM_BC, which stands for the dimension of the flow problem. The value is 2 in 2d, for 2d-axisymmetric and 2d-axisymmetric-swirling flow, while it is 3 in 3d flows.

You are looping through the components for the vector "normal" and setting an equivalence between "normal" and "f_normal", which this last one is another argument of DEFINE_DPM_BC that contains the unit vector which is normal to the face (the boundary condition where DEFINE_DPM_BC is hooked to).

if ((NNULLP(t)) && (THREAD_TYPE(t) == THREAD_F_WALL))
F_CENTROID(x,f,t);


(Copied from the udf guide): You can use the NULLP and NNULLP functions to check whether storage has been allocated for userdefined scalars. NULLP returns TRUE if storage is not allocated, and NNULLP returns TRUE if storage is allocated.
If you are not using UDS, I guess you can omit it.

(THREAD_TYPE(t) == THREAD_F_WALL) is used to check if you are in a wall

F_CENTROID(x,f,t) sets the array "x" as the coordinates of the face centroid of the wall.

Cheers

mimi0201 October 29, 2015 18:43

Quote:

Originally Posted by upeksa (Post 570806)
for (i=0; i<idim; i++)
normal[i] = f_normal[i];


int idim=dim, where dim is argument in DEFINE_DPM_BC, which stands for the dimension of the flow problem. The value is 2 in 2d, for 2d-axisymmetric and 2d-axisymmetric-swirling flow, while it is 3 in 3d flows.

You are looping through the components for the vector "normal" and setting an equivalence between "normal" and "f_normal", which this last one is another argument of DEFINE_DPM_BC that contains the unit vector which is normal to the face (the boundary condition where DEFINE_DPM_BC is hooked to).

if ((NNULLP(t)) && (THREAD_TYPE(t) == THREAD_F_WALL))
F_CENTROID(x,f,t);


(Copied from the udf guide): You can use the NULLP and NNULLP functions to check whether storage has been allocated for userdefined scalars. NULLP returns TRUE if storage is not allocated, and NNULLP returns TRUE if storage is allocated.
If you are not using UDS, I guess you can omit it.

(THREAD_TYPE(t) == THREAD_F_WALL) is used to check if you are in a wall

F_CENTROID(x,f,t) sets the array "x" as the coordinates of the face centroid of the wall.

Cheers

Thank you so much for the reply, and i have several questions.
1. So u mean dim is argument in DEFINE_DPM_BC, but why don't we just use dim then, why we need to define dim as idim?
2. Similarly, what is the purpose to define normal[i] = f_normal[i]?
3. What does the statement "F_CENTROID(x,f,t)" do? it calls the controid if statement "((NNULLP(t)) && (THREAD_TYPE(t) == THREAD_F_WALL))" is true, but why do we need the centroid?

I'm new to UDF so I have a lot of silly questions. Really appreciate!

pakk October 30, 2015 07:54

Quote:

Originally Posted by mimi0201 (Post 570923)
Thank you so much for the reply, and i have several questions.
1. So u mean dim is argument in DEFINE_DPM_BC, but why don't we just use dim then, why we need to define dim as idim?

I see no reason in this code. It would have been easier to use dim.
Quote:

2. Similarly, what is the purpose to define normal[i] = f_normal[i]?
I see no purpose in this code. It would have been easier and more clear to use f_normal everywhere.
Quote:

3. What does the statement "F_CENTROID(x,f,t)" do? it calls the controid if statement "((NNULLP(t)) && (THREAD_TYPE(t) == THREAD_F_WALL))" is true, but why do we need the centroid?
Here it gives 'x' the coordinates of the center of the cell where the particle is in. But 'x' is never used again, so this is not really useful to do.
Quote:

I'm new to UDF so I have a lot of silly questions. Really appreciate!
These are definitely not silly questions, I think they are good questions for somebody new to UDF. This example has some unnecessary lines of code, and you have seen that they have no purpose. I think you understand UDFs better than some other people!

mimi0201 October 31, 2015 21:49

Quote:

Originally Posted by pakk (Post 571010)
I see no reason in this code. It would have been easier to use dim.

I see no purpose in this code. It would have been easier and more clear to use f_normal everywhere.

Here it gives 'x' the coordinates of the center of the cell where the particle is in. But 'x' is never used again, so this is not really useful to do.

These are definitely not silly questions, I think they are good questions for somebody new to UDF. This example has some unnecessary lines of code, and you have seen that they have no purpose. I think you understand UDFs better than some other people!

So as you explained, codes shown below are all redundant if I just want to track the particle hitting on the wall:
Code:

dim=idim;
normal[i] = f_normal[i];
((NNULLP(t)) && (THREAD_TYPE(t) == THREAD_F_WALL))";
F_CENTROID(x,f,t);

But which code shows that, when particle hitting on the wall then calculate the normal velocity etc.?
Cheers!

pakk November 2, 2015 05:32

This code is a boundary condition. So the code is only run when a particle has reached a wall. And the steps that are taken then, are in the code. Look at the comments!
Code:

        /* Compute normal velocity. */
        /* Subtract off normal velocity. */
        /* Apply tangential coefficient of restitution. */
        /* Add reflected normal velocity. */
        /* Store new velocity in P_VEL0 of particle */


mimi0201 November 2, 2015 17:45

Quote:

Originally Posted by pakk (Post 571415)
This code is a boundary condition. So the code is only run when a particle has reached a wall. And the steps that are taken then, are in the code. Look at the comments!
Code:

        /* Compute normal velocity. */
        /* Subtract off normal velocity. */
        /* Apply tangential coefficient of restitution. */
        /* Add reflected normal velocity. */
        /* Store new velocity in P_VEL0 of particle */


Oh ok! That was a bit tricky,haha. I think I get it now, i was stuck in the beginning of the code for a long time.


All times are GMT -4. The time now is 22:38.