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/)
-   -   How to use a UDF to set the volume fraction in the cells next to a wall? (https://www.cfd-online.com/Forums/fluent-udf/152350-how-use-udf-set-volume-fraction-cells-next-wall.html)

DF15 April 29, 2015 06:58

How to use a UDF to set the volume fraction in the cells next to a wall?
 
Hi,

I have an issue regarding a UDF to be used in a multiphase problem where water and air flows around a rigid body.

In my results from a simulation I have an unphysical layer of air in the first cell layer on a wall which I know should be completely wet, i.e. the volume fraction of air should be zero. I would therefore like to create a UDF which loops through all cells adjacent to this wall, and if the volume fraction of air is below a certain limit I want to set the volume fraction to 0.

I have constructed a UDF with a DEFINE_ADJUST macro in which I loop through all faces on the wall (this should work because I've taken it from another UDF used to compute the area of the same wall). I think the problem is that I have two phases and I want to evaluate the volume fraction just for one phase (air). I have not fully understood how to use subthreads and subdomains in multiphase simulations.

When I run a time step, I get the following error:
"Received signal SIGSEGV"

Can anyone help me with this problem? The current UDF looks as follows. It crashes when entering the second #if !RP_HOST loop.

Best regards,
David

Code:

#include "udf.h"

DEFINE_ADJUST(remove_air,domain)
{
    /* Cut-off volume fraction */
    real VOF_limit=0.3;            /* Volume fraction limit (of air) */
   
    /* Get the thread pointer*/
    #if !RP_HOST
        Thread *t;            /* Mixture level thread */
        Thread *t_air;            /* Phase level thread */
        face_t f;            /* Face thread */
        cell_t c;            /* Cell thread */
        domain = Get_Domain(1);    /* Get the domain using ANSYS Fluent utility */
    #endif
     
    /* Set thread id of wall (from Fluent) */
    int surface_thread_id=12;        /* Wall thread id */
    host_to_node_int_1(surface_thread_id);    /* Pass thread id to nodes */
    int phase_domain_index = 0;        /* Air phase index (primary phase) */
    host_to_node_int_1(phase_domain_index);    /* Pass domain index to nodes */

    /* Loop through nodes */
    #if !RP_HOST
        /* Find threads of mixture and air */
        t = Lookup_Thread(domain,surface_thread_id);
        t_air = THREAD_SUB_THREAD(t, phase_domain_index);   

        /* Loop through all faces on wall */
        begin_f_loop(f,t)   
        {   
            /* Get thread for adjacent cell */
            c = F_C0(f,t);
   
            /* Reset volume fraction if it is below the limit */
            if (C_VOF(c,t_air) <= VOF_limit)
            {
                C_VOF(c,t_air)=0;
            }
        }
        end_f_loop(f,t)

    #endif
}


`e` April 29, 2015 07:49

Declare your cell thread variable, "cell_t c", with the other declarations (not halfway through your code, this placement might technically be legal with ANSI C but there's nothing else obviously causing the malloc error).

There's a number of lines in your code, try to narrow down the cause of the error. Either reduce the lines of code or add Message() to print information to the screen about where the process reaches before crashing.

DF15 April 29, 2015 08:07

Thank you for your advice, I am not very experienced in programming so my UDF is not structured very well.

However, I now declared the cell-thread earlier, and I also used host_to_node to pass the phase_domain_index to the nodes as I did with the surface_thread_id. Now, it seems like the crash occurs when entering the node-loop; it never reaches the line "t = Lookup_Thread(domain,surface_thread_id);" but it does reach just before this loop starts. I have used the same loop earlier without problems, so I don't see why it crashes when entering the node-loop now... The error I get is just

"Received signal SIGSEGV"

The code in my first post is updated to the current version.

`e` April 29, 2015 08:39

Are you sure the code falls over at initialising "t" with Lookup_Thread? Did you check by removing the other lines of code? Did this section of code work with the air thread pointer?

The cell thread, "c", might be different between the two phases (I'm not familiar with the Euler-Euler multiphase model in Fluent). If Fluent is searching for a cell which is not within a thread then I'd expect a SIGSEGV error. For example your code:

Code:

if (C_VOF(c,t_air) <= VOF_limit)
{
    C_VOF(c,t_air)=0;
}

may need to be along the lines of:

Code:

if (C_VOF(c_air,t_air) <= VOF_limit)
{
    C_VOF(c_air,t_air)=0;
}


DF15 April 29, 2015 09:00

I put Message-lines in the code, and it prints the message just before "#if !RP_HOST" but not the message just before "t = Lookup_Thread(domain,surface_thread_id);" The Lookup_thread thing worked in an earlier case where I looped over the surface faces and calculated the total area.

Now when I added a Message line to the first #if - part with the thread declarations, I see that it does not go in there. However, it does not crash there but continues to the next #if. What could this be?

UPDATE: If it crashes before the part of the code that you referred to about the cell thread, could it still be because something is wrong there?

DF15 April 29, 2015 10:59

UPDATE: Now I have removed the parallelization-stuff since I only want to replace the volume fraction value. Now, when I want to loop through all the faces on the wall, I do not manage to get the loop working. It gets to message 1, but not to message 2, so it wont enter the loop.

Have I constructed the face-loop correctly, or should there be another thread t that I run the loop though? The part of the UDF containing these things looks as follows.

/David

Code:

        t = Lookup_Thread(domain,surface_thread_id);
        t_air = THREAD_SUB_THREAD(t, phase_domain_index);
        Message("\n 1");
       
        /* Loop through all faces on wall */
        begin_f_loop(f,t)
        {   
            Message("\n 2");
            /* Get thread for adjacent cell */
            c = F_C0(f,t);
            ...


`e` April 29, 2015 23:29

Quote:

Originally Posted by DF15 (Post 544306)
Now when I added a Message line to the first #if - part with the thread declarations, I see that it does not go in there. However, it does not crash there but continues to the next #if. What could this be?

UPDATE: If it crashes before the part of the code that you referred to about the cell thread, could it still be because something is wrong there?

If it's not declaring those threads at the start then that's the problem you need to solve. I don't understand why that block of code isn't being called, perhaps try isolating this block of code. Another approach of simplifying your code would be to look at single phase first because using multiple phases (different threads) might be complicating things and I don't have experience with this area (but nothing seems obviously out of place).

Quote:

Originally Posted by DF15 (Post 544325)
Have I constructed the face-loop correctly, or should there be another thread t that I run the loop though? The part of the UDF containing these things looks as follows.

The face loop looks OK but perhaps use a principal face check to ensure you're only adjusting each face once (read the parallelisation section in the UDF manual).

DF15 April 30, 2015 07:56

Hello again! Thank you `e` for your help, I am really thankful! I have now worked further on my UDF and tried to look at a previous, similar UDF which worked fine. I have come up with a working test code, which is run at each iteration without any errors. To test the face loop, I am now summing the viscous force in x-direction and compare to the report from Fluent and they agree. Also, I figured I should set the volume fraction of water to 1 if I set that of air to 0, so I have made a phase thread for water as well.

Before I try to overwrite the volume fraction in the cells where it is below a certain limit, I wanted to see if I could use the C_VOF to get any values from it. So I created a variable volume_fraction_air which I print at every face if it's within a certain range to see if the volume fractions add up to 1. When I look in Fluent's terminal, I see that the water volume fractions are not 1-volume_fraction air, and they are sometimes negative. Any ideas what this could be? Here are some different lines that are printed:

Volume fraction air: 0.386194
Volume fraction water: -0.000000
Volume fraction air: 0.384285
Volume fraction water: 0.000000
Volume fraction air: 0.384660
Volume fraction water: 30067229385242264364586082335177459086538981525062 95365192962563293835631977532052698749741028725707 767808.000000

The code looks as follows:

Code:

#include "udf.h"

DEFINE_ADJUST(remove_air,domain)
{
    /* Cut-off volume fraction */
    real VOF_limit=0.3;            /* Volume fraction limit (of air) */

    /* Create the thread pointers for the wall */
    Thread *t;            /* Mixture level thread */
    Thread *t_air;            /* Phase level thread */
    Thread *t_water;        /* Phase level thread */
    face_t f;            /* Face thread */
    cell_t c;            /* Cell thread */
    domain = Get_Domain(1);    /* Get the domain */

    /* Set thread id's of wall (from Fluent) and phase domain indices */
    int surface_thread_id=12;        /* Wall thread id */
    host_to_node_int_1(surface_thread_id);    /* Pass thread id to nodes */
    int phase_domain_index_air = 0;        /* Air phase index */
    int phase_domain_index_water = 1;    /* Water phase index */
    host_to_node_int_2(phase_domain_index_air, phase_domain_index_water);    /* Pass domain index to nodes */


    /* Calculate viscous resistance */
    real R_V=0.0;
    real volume_fraction_air;
    real volume_fraction_water;

    /* Node calculations */
    #if !RP_HOST
        /* Find mixture level thread and phase level threads */
        t = Lookup_Thread(domain,surface_thread_id);
        t_air = THREAD_SUB_THREAD(t, phase_domain_index_air);
        t_water = THREAD_SUB_THREAD(t, phase_domain_index_water);

        /* Loop through faces on wall */
      begin_f_loop(f,t)
        {
        if (PRINCIPAL_FACE_P(f,t))
        {
            /* Get thread for adjacent cell */
            c = F_C0(f,t);

            volume_fraction_air=C_VOF(c,t_air);
            volume_fraction_water=C_VOF(c,t_water);
            /* Message("Volume fraction air %f\n", volume_fraction_air); */

            if (volume_fraction_air >= 0.38 && volume_fraction_air <= 0.39)
            {
                Message("Volume fraction air:    %f\n", volume_fraction_air);
                Message("Volume fraction water: %f\n", volume_fraction_water);
            }
           
            R_V -= F_STORAGE_R_N3V(f,t,SV_WALL_SHEAR)[0];
        }
        }
        end_f_loop(f,t)

        /* Sum values from nodes */
        R_V = PRF_GRSUM1(R_V);
    #endif
    node_to_host_real_3(R_V, volume_fraction_air, volume_fraction_water);

    /* Print variables to validate UDF */
    #if !RP_NODE
        Message("Viscous force: %f [N]\n", -R_V);
        Message("Volume fraction air: %f\n", volume_fraction_air);
        Message("Volume fraction water: %f\n", volume_fraction_water);
    #endif
}


`e` April 30, 2015 21:29

Let's walk through the code. DEFINE_ADJUST is executed at each iteration before the solver is called and you've declared a real variable "volume_fraction_water" on all nodes and on the host process. I assume your messages you have posted is from the host process (or serial in serial mode) because of the absence of the space. Why is there no viscous force reported? Are the volume fractions reported at each iteration or for every face for every iteration?

Each node has its own partition and therefore set of faces to look after. You're overwriting the "volume_fraction_water" variable every time through the face loop. If for whatever reason a node has no faces from this boundary associated with it (its partition is elsewhere in the domain) then you're sending whatever was previously in memory at the declared space in memory for "volume_fraction_water" (rubbish, which would explain your huge unexpected number).

DF15 May 1, 2015 02:40

Quote:

Originally Posted by `e` (Post 544570)
Let's walk through the code. DEFINE_ADJUST is executed at each iteration before the solver is called and you've declared a real variable "volume_fraction_water" on all nodes and on the host process. I assume your messages you have posted is from the host process (or serial in serial mode) because of the absence of the space.

Ok, so when I declare a real variable outside any #if RP_NODE or #if RP_HOST, then it is declared on each node and on the host. The messages I posted were written manually because I couldn't find a was to copy-paste from the Fluent terminal. What do you mean with the absence of space?

Quote:

Originally Posted by `e` (Post 544570)
Why is there no viscous force reported? Are the volume fractions reported at each iteration or for every face for every iteration?

These were just some lines I took from the middle of the long row of messages. At the end of all the messages the viscous force is printed once (and it is correct). The volume fractions are reported for every adjacent cell for every iteration, which I suppose is what I want since I would finally like to run through all the cells at the wall in every iteration and replace the volume fraction if some condition is fulfilled.

Quote:

Originally Posted by `e` (Post 544570)
Each node has its own partition and therefore set of faces to look after. You're overwriting the "volume_fraction_water" variable every time through the face loop. If for whatever reason a node has no faces from this boundary associated with it (its partition is elsewhere in the domain) then you're sending whatever was previously in memory at the declared space in memory for "volume_fraction_water" (rubbish, which would explain your huge unexpected number).

Alright, I think I understand that in general there might be incorrect volume fractions reported if there are no wall-faces in a specific node. But, since it prints the volume fraction of air within my specified range and I try to print the volume fraction of water in the same cell, shouldn't there be realistic numbers in that case?

Since the viscous force is reported correctly, I assume that the face loop runs through the desired cell faces (i.e. only those on the wall). But this does not seem to be the case judging from the volume fraction report. If a node has no faces on the wall associated to it, shouldn't the face loop just skip this node since it has no cells in the mixture level thread t? Could it be that I have used the face-loop or the C_VOF function incorrectly, i.e. with the wrong type of threads as input? I have tried to look at the first example found here: http://jullio.pe.kr/fluent6.1/help/html/udf/node115.htm , and they seem to use the same threads for C_VOF...

`e` May 1, 2015 03:16

Yes, I believe DEFINE_ADJUST is called on all of the nodes including the host and therefore these processors will execute the lines within your code block. Select the text and then right click > copy (Ctrl + C is reserved for interrupting Fluent; similar to a terminal).

You have two sets of messages of the following form:

Code:

Message("Volume fraction air:    %f\n", volume_fraction_air);
One with a bunch of spaces before the number (as seen above) and the second set is without these spaces. I could only assume your output was from the second set of messages because there was no bunch of spaces: "absence of the space". It's unclear where you're printing the "volume_fraction_water" variable. You've now said that it's successfully printing for every adjacent cell for every iteration and those messages are a short excerpt (OK). However, the garbage number must be from the second set of messages because there's no case where the first set of messages is called with an uninitialised "volume_fraction_water" variable.

You've initialised your viscous force variable, "R_V", with zero at the place of declaring said variable. Therefore any nodes or host process which executes this code has a default value of zero assigned to "R_V". Whereas your other two volume fractions has no such initialisation and therefore returns garbage for nodes and the host which have no relevant faces assigned to them.

DF15 May 1, 2015 07:41

Ok, I see what you mean. Initializing the volume fraction variables with a value is maybe better, I guess I can set them to 1 so that the if-statement is not entered in case they are not overwritten by a new value in the present cell. Do you see what I mean?

I have tried to simplify the code a bit, here is the current one. It runs without problems, and I can see that the whole code is run through on each of the nodes I have. For some nodes, no volume fractions are printed which should indicate that the volume fraction is not within the specified range OR the partition associated to the node does not contain a part of the wall. Though, I still cannot see why the volume fraction of water can be something else than 1-volume_fraction_air, can you?

The "garbage numbers" are now gone, but there are still some (small) negative volume fractions of water for the corresponding cell (I hope that they correspond at least). Maybe it isn't as easy as defining the sub threads for the two phases as I do, could it be that the phase level threads point to different cells or something? A typical terminal output is (I can't right click to copy on this Linux system...).

Volume fraction air: 0.299322
Volume fraction water: -0.000000
Volume fraction air: 0.287581
Volume fraction water: 0.000000
Volume fraction air: 0.287909
Volume fraction water: 0.000000
Volume fraction air: 0.287016
Volume fraction water: 0.000041

Code:

#include "udf.h"

DEFINE_ADJUST(remove_air,domain)
{
Message("0");

    /* Cut-off volume fraction */
    real VOF_limit=0.3;        /* Volume fraction limit (of air) */
   
    /* Create the thread pointers for the wall */
    Thread *t;            /* Mixture level thread */
    Thread *t_air;            /* Phase level thread */
    Thread *t_water;        /* Phase level thread */
    face_t f;            /* Face thread */
    cell_t c;            /* Cell thread */
    domain = Get_Domain(1);    /* Get the domain */
   
    /* Calculate viscous resistance */
    real R_V=0.0;
    real volume_fraction_air=1.0;
    real volume_fraction_water=1.0;

    /* Node calculations */
    #if !RP_HOST
        /* Find mixture level thread and phase level threads */
        t = Lookup_Thread(domain,12);        /* Surface thread ID from Fluent */
        t_air = THREAD_SUB_THREAD(t, 0);    /* Primary phase */
        t_water = THREAD_SUB_THREAD(t, 1);    /* Secondary phase */
       
        /* Loop through all faces on surface */
        begin_f_loop(f,t)
        {
        if (PRINCIPAL_FACE_P(f,t))
        {
            /* Get thread for adjacent cell */
            c = F_C0(f,t);
           
            /* Get volume fractions of air and water */
            volume_fraction_air=C_VOF(c,t_air);
            volume_fraction_water=C_VOF(c,t_water);

            /* Update volume fractions */
            if (volume_fraction_air >= 0.287 && volume_fraction_air <= 0.3)
            {
                Message("Volume fraction air:      %f\n", volume_fraction_air);
                Message("Volume fraction water:  %f\n", volume_fraction_water);
            }

            /* Calculate viscous force (to test UDF) */
            R_V -= F_STORAGE_R_N3V(f,t,SV_WALL_SHEAR)[0];
        }
        }
        end_f_loop(f,t)
       
        /* Sum values from each node */
        R_V = PRF_GRSUM1(R_V);
    #endif
    /* Send viscous force from nodes to host */
    node_to_host_real_1(R_V);
   

    /* Print results (to test UDF) */
    #if !RP_NODE
        Message("Viscous force: %f [N]\n", -R_V);
    #endif
}


`e` May 2, 2015 02:34

I haven't used VOF with Fluent before and have little knowledge of this area. The way you've called the volume fractions with C_VOF appear to agree with how the developers have shown it in the UDF manual.

It seems strange that the two volume fractions don't add to unity; is there another phase? Perhaps try setting a constant volume fraction for your entire domain and boundary for each phase and then run this script to check if your code is reading the correct volume fraction values.

How small are the negative numbers, could you check with using the "%e" format in message? If they're on the order of numeric precision then I wouldn't worry about it. If you're concerned about these negative values, you could place a lower limit on the volume fraction which is very small but above the numeric precision.

Are you suggesting the cell threads, "c", could be unique between the phases? I wondered about this before. Perhaps check the Fluent documentation if you're not receiving the correct volume fractions after initialising them as a constant (suggestion above).

DF15 May 4, 2015 08:14

Quote:

Originally Posted by `e` (Post 544698)
It seems strange that the two volume fractions don't add to unity; is there another phase?

No, there are only two phases - air and water.

Quote:

Originally Posted by `e` (Post 544698)
How small are the negative numbers, could you check with using the "%e" format in message? If they're on the order of numeric precision then I wouldn't worry about it. If you're concerned about these negative values, you could place a lower limit on the volume fraction which is very small but above the numeric precision.

Some of them are in the order of 0.0001, which I suppose is above the numeric precision? I run Fluent in double precision if that's related to the numeric precision.

Quote:

Originally Posted by `e` (Post 544698)
Are you suggesting the cell threads, "c", could be unique between the phases? I wondered about this before. Perhaps check the Fluent documentation if you're not receiving the correct volume fractions after initialising them as a constant (suggestion above).

I can only find that the cell and face threads are related to cells and faces, and that the different phases are accounted for through the phase level threads (t_air and t_water in this case).

I will try to initialize the volume fraction in a simple manner, but according to the checks I have made so far, the volume fractions are not reported correctly with my use of C_VOF...

Katja September 6, 2016 10:07

Hello,
Sorry to ask my question here, may be you can help me with this problem. I am simulating a multiphase flow where one of the solid phases should stay in the domain while the other fluid phase can leave at the outlet. How I can restrict this phase from leaving at the outlet? Do I need a UDF for that or is there a setting in Fluent where I can simply do that?

Thank you!

Rubblevg March 14, 2017 10:12

Quote:

Originally Posted by Katja (Post 616798)
Hello,
Sorry to ask my question here, may be you can help me with this problem. I am simulating a multiphase flow where one of the solid phases should stay in the domain while the other fluid phase can leave at the outlet. How I can restrict this phase from leaving at the outlet? Do I need a UDF for that or is there a setting in Fluent where I can simply do that?

Thank you!

I have exactly the same problem as you! Have you been able to fix it yet? Does anyone else know how to solve this?

sewgolam September 7, 2017 09:00

volume fraction
 
Dear friends,

i also need the folume fraction of my secondary phase in each cell. I tryed to make a uDF but couldnot get it work. can you please help

Mel B December 14, 2017 02:03

hello,
i want to cut-off air volume fraction when it is greater than 0.01 for able to an realistic results in eulerian model.
i need to write UDF.

sand zone ID=5 and 1th secondary phase (air) ID=3

my UDF;

DEFINE_ADJUST(mixture_domain, phase_domain_index)
{
int phase_domain_index = PHASE_DOMAIN_INDEX(subdomain);
cell_t cell;
Thread *cell_thread;
thread **pt
Domain *subdomain;
Domain *mixture_domain;
real VOF_limit=0.01;

/* loop over all subdomains (phases) in the superdomain (mixture) */
sub_domain_loop(subdomain, mixture_domain,phase_domain_index)
{
/* loop if air phase */
if (DOMAIN_ID(subdomain) == 3)

/* loop over all cell threads in the air phase domain */
thread_loop_c(cell_thread,subdomain)
{
/* loop if air phase */
if (THREAD_ID(cell_thread) == 5)
{
/* loop over all cells in air phase cell threads */
begin_c_loop_all (cell,cell_thread)
{
C_VOF(c,pt[3]);
if ( C_VOF(c,pt[3])=>VOF_limit

/* set volume fraction to 0 for air phase */
C_VOF(c,pt[3]) = 0.;
}
end_c_loop_all (cell,cell_thread)
}

}
}
}



I don't know the correctness of the code I wrote
also it gives parse error and declared variable during interpretting
someone that can help me?

obscureed December 14, 2017 05:50

Hi Mel B,

It might have been better to start a new thread, with a link to this thread if you thought it was particularly relevant.

You should try to fix minor compilation bugs before you post code. It will help you (and everyone) if you can get an editor to auto-indent the code, to make sure brackets etc line up. (Maybe you've done this, but it doesn't show up if you don't quote code as code in the forum post.) This should point you to some errors:
-- "thread **pt" is missing a semicolon.
-- The line "if ( C_VOF(c,pt[3])=>VOF_limit" is a mess -- you must mean ">=" and the brackets need fixing. (As an aside, the line before it does nothing.)
-- The line "if (DOMAIN_ID(subdomain) == 3)" needs some "{ ... }" braces.

In the bigger picture, poking values into solution variables is generally a very bad idea. You are basically fighting against the solver: the solver puts values into the solution variables, then you put in different values, then repeat. You never know which values are present at any instant. Also, the solver sometimes compares its new values with previous values, and may diverge if there are discrepancies that it does not expect.

Really, the only safe way to manipulate solution variables is via source terms. This should remind you that you're doing something "extra" to the physics: you can't just tell the solver to engineer a value somewhere without giving the solver a reason to do so. (That reason is the source term. You don't have to explain how the source term happened; as far as the solver is concerned, the extra mass or whatever is just teleported in/out of that location.) So, for example, if you code up a source term that teleports secondary phase out of any cell that achieves C_VOF > 0.01, then of course you should not expect mass conservation. You also need to think: is material actually removed by the source term, or is it transmuted to primary phase? (Probably the latter. If you have a compressible fluid, or a density-based solver, this is more complicated.) And then you need to think: is this what you intended?

Good luck!
Ed

Mel B December 15, 2017 01:46

firstly i would like thank you
obscureed, you right.
my model is consist three phase for open channel flow.(sediment transport in open channel)
(usually researchers use only two phases in the similar flow problems. I think this is not appropriate for open channel.)
therefore I must also use the air phase in addition to the water and sand phase.
primary phase is water, secondary phases are air and sand.
but there is several problem in my solution results.
should not be air phase in realistic flow field of experimental sand bed. maybe it is should be ignored value for air volume fraction. but
there is air volume fraction in sand zone as a usual result of the solution with three phases in numerical solution.
but volume fraction of air phase in sand zone is reaching large values(%10) in my solution.
i want to cut-off air volume fraction. maybe this idea is entirely unreasonable. i should check my problem.

maybe the value of air phase is a condition resulting from inappropriate initial conditions.
packing limit values is %63 for monodisperse spherical particle. this is causes the water volume fraction to be large.
but sand zone in the experiment has different spherical particles, so packing limit value is greater than 63%.
can we reach results about that model is not appropriate for sediment transport?


All times are GMT -4. The time now is 02:07.