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

how to loop over cells in a customised way ( UDM usage )

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

Like Tree1Likes
  • 1 Post By blackmask

Reply
 
LinkBack Thread Tools Display Modes
Old   July 10, 2013, 05:21
Default how to loop over cells in a customised way ( UDM usage )
  #1
Senior Member
 
Ananthakrishnan.A.S
Join Date: Feb 2012
Location: Mumbai (Bombay), India
Posts: 157
Rep Power: 8
Ananthakrishnan is on a distinguished road
Send a message via Skype™ to Ananthakrishnan
Hello everybody,

Thanks a lot for your help in advance. First of all, the title of this thread could be misleading, because this is the best i could think of for describing my problem!!!

I have shown some parts of my program below. Its function is to split the flow domain based on temperature and store it in a UDM through "reaction ID". But now , I need to loop over the cells in a particular reaction ID and i dont know how to do that. I would really appreciate if anyone can give any sort of help regarding this.. Thanks a lot in advance..

please let me know if i was not clear in my explanations.

#include "udf.h"
enum
{
CHEM_REACTION_ID,
N_CHEM_REACTION_UDM
};
#define C_CHEM_REACTION_ID(C,T) C_UDMI(C,T,CHEM_REACTION_ID)
#define PDF_SPECIES_NOT_FOUND -1

DEFINE_ON_DEMAND(set_chemical_reaction_regions)
{
Domain *d;
Thread *ct;
cell_t c;
real cen[ND_ND],x,y,z;
int reaction_ID;
int species_index;
real mass_fraction;

d = Get_Domain(ROOT_DOMAIN_ID);
species_index = find_pdf_species_index("o2");

thread_loop_c(ct,d)
if(FLUID_THREAD_P(ct))
begin_c_loop(c,ct)
{
/* Temperature segregation */
if (C_T(c,ct) < 1000.0)
reaction_ID = 12;
else if (C_T(c,ct) < 1500.0)
reaction_ID = 9;
else if (C_T(c,ct) < 2000.0)
reaction_ID = 6;
else if (C_T(c,ct) < 2250.0)
reaction_ID = 3;
else
reaction_ID = 0;
C_CHEM_REACTION_ID(c,ct) = reaction_ID;
}
end_c_loop(c,ct)
Message("Done.\n");
}
Ananthakrishnan is offline   Reply With Quote

Old   July 10, 2013, 05:42
Default
  #2
Senior Member
 
Join Date: Aug 2011
Posts: 315
Rep Power: 12
blackmask will become famous soon enough
Either you loop over all the cells and find those with particular reaction id(s), or you record the cell id c and thread id ct for particular reaction id(s) with arrays in DEFINE_ON_DEMAND and then loop over the arrays.
blackmask is offline   Reply With Quote

Old   July 10, 2013, 06:35
Default
  #3
Senior Member
 
Ananthakrishnan.A.S
Join Date: Feb 2012
Location: Mumbai (Bombay), India
Posts: 157
Rep Power: 8
Ananthakrishnan is on a distinguished road
Send a message via Skype™ to Ananthakrishnan
Hi,
Thanks a lot for your tip.. Is it somethinglike this

option 1
thread_loop_c(ct,d)
begin_c_loop(c,ct)
if C_CHEM_REACTION_ID(c,ct) = 1;
....
else if C_CHEM_REACTION_ID(c,ct) = 2;
....

or do i have to use reaction_ID in place of C_CHEM_REACTION_ID(c,ct). I am really confused how this works..

option 2

Am i not already recording the cell id and thread id for each reaction ID

I am really sorry but I will be really thankful if you can provide me a head start by writing a couple of lines of codes..

Thanks a lot,
Ananthakrishnan is offline   Reply With Quote

Old   July 10, 2013, 06:52
Default
  #4
Senior Member
 
Join Date: Aug 2011
Posts: 315
Rep Power: 12
blackmask will become famous soon enough
Then let discuss option 1 only. It goes like this
Code:
Domain *d;
Thread *ct;
cell_t c;

d = Get_Domain(ROOT_DOMAIN_ID);

thread_loop_c(ct,d) {
    if(FLUID_THREAD_P(ct))
      begin_c_loop(c,ct)
        {
          
  switch ( C_CHEM_REACTION_ID(c,ct) ) {
case 3:
//do something 
break;
case 6:
//do something
break;
...
case default:

break;
}
        }
      end_c_loop(c,ct)
}
Ananthakrishnan likes this.
blackmask is offline   Reply With Quote

Old   July 10, 2013, 07:00
Default
  #5
Senior Member
 
Ananthakrishnan.A.S
Join Date: Feb 2012
Location: Mumbai (Bombay), India
Posts: 157
Rep Power: 8
Ananthakrishnan is on a distinguished road
Send a message via Skype™ to Ananthakrishnan
Thanks a lot,

I will try it out and if necessary will come back to you..!!

Thanks once again..
Ananthakrishnan is offline   Reply With Quote

Old   July 10, 2013, 11:13
Default
  #6
Senior Member
 
Ananthakrishnan.A.S
Join Date: Feb 2012
Location: Mumbai (Bombay), India
Posts: 157
Rep Power: 8
Ananthakrishnan is on a distinguished road
Send a message via Skype™ to Ananthakrishnan
Hi,
I am getting the error "switch expression not integral".. any idea??

thanks in advance
Ananthakrishnan is offline   Reply With Quote

Old   July 10, 2013, 19:42
Default
  #7
Senior Member
 
Join Date: Aug 2011
Posts: 315
Rep Power: 12
blackmask will become famous soon enough
switch ( C_CHEM_REACTION_ID(c,ct) )
=>
switch ( (int) C_CHEM_REACTION_ID(c,ct) )
blackmask is offline   Reply With Quote

Old   July 11, 2013, 04:38
Default
  #8
Senior Member
 
Ananthakrishnan.A.S
Join Date: Feb 2012
Location: Mumbai (Bombay), India
Posts: 157
Rep Power: 8
Ananthakrishnan is on a distinguished road
Send a message via Skype™ to Ananthakrishnan
Hi,
Thanks,

The error disappeared but i dont think it is looping over the cells, what i can think of is C_CHEM_REACTION_ID(c,ct) does not contain the reaction IDs (3, 6, etc..) . I just tried to display the temperature in that thread and it returned zero.

thread_loop_c(ct,d) {
if(FLUID_THREAD_P(ct))
begin_c_loop(c,ct)
{

switch ( (int) C_CHEM_REACTION_ID(c,ct) ) {
case 3:
Message("\n temp is %g",C_T(c,ct));
break;
}
}
end_c_loop(c,ct)
Ananthakrishnan is offline   Reply With Quote

Old   July 11, 2013, 07:40
Default
  #9
Senior Member
 
Ananthakrishnan.A.S
Join Date: Feb 2012
Location: Mumbai (Bombay), India
Posts: 157
Rep Power: 8
Ananthakrishnan is on a distinguished road
Send a message via Skype™ to Ananthakrishnan
Hi again,
Instead of switch case i tried
if ( C_CHEM_REACTION_ID(c,ct) < 3 or 4 or whatever)

and it worked..
I guess it is not necessary for C_CHEM_REACTION_ID(c,ct) to have integers..No clue why!! but i am getting some logical answers but still i dont know if i can depend on it or not!!!!
Ananthakrishnan is offline   Reply With Quote

Old   July 11, 2013, 08:06
Default
  #10
Senior Member
 
Join Date: Aug 2011
Posts: 315
Rep Power: 12
blackmask will become famous soon enough
The type of C_CHEM_REACTION_ID(c,ct) should be real so that it should be casted to integer in the switch expression. The 'lround' function is more accurate but I think 'int' is enough for small integers. You can print the reaction ids and temperature for cells to see whether it works as expected, i.e.
id = 12, T < 1000
id = 9, T < 1500
id = 6, T < 2000
id = 3, T < 2250
id = 0, otherwise
blackmask is offline   Reply With Quote

Old   July 11, 2013, 11:41
Default
  #11
Senior Member
 
Ananthakrishnan.A.S
Join Date: Feb 2012
Location: Mumbai (Bombay), India
Posts: 157
Rep Power: 8
Ananthakrishnan is on a distinguished road
Send a message via Skype™ to Ananthakrishnan
Hi,
thanks a lot..

yup..it works fine..checked.. I am sticking with the If loop..
I am really sorry to bother you even more but it will be really nice if you can help me out with two other doubts

1. do you know how fluent calculates the velocity magnitude..I tried a couple of formulas but i dont get anything similar to Fluent..

2.This one is regarding my second part of the program

I also need to find out the mass flow rate between each of the reaction IDs.. I am able to find out the over all mass flow with the program below but i am unable to identify the mass flow separately between two pairs..

thread_loop_f(ft,domain)
if(THREAD_TYPE(ft) == THREAD_F_INTERIOR)
begin_f_loop(f,ft)
{

id0 = C_CHEM_REACTION_ID(F_C0(f,ft) THREAD_T0(ft));

id1 = C_CHEM_REACTION_ID(F_C1(f,ft) THREAD_T1(ft));
if(id0 != id1)
{
mass_FLUX += F_FLUX(ft,ct);
//is there a way to keep track of the id0,id1 so that i can know what is the mass flux between two given reaction IDs..
}

end_f_loop(f,ft)

Thank a lot once again...
Ananthakrishnan is offline   Reply With Quote

Old   July 11, 2013, 21:51
Default
  #12
Senior Member
 
Join Date: Aug 2011
Posts: 315
Rep Power: 12
blackmask will become famous soon enough
1. Should not the magnitude be the square root of its components for each phase? Is there any chance that the phase does not match in your comparison?

2. I think that the UDM should also be defined for the face threads so that you can use
Code:
 C_CHEM_REACTION_ID(f, ft) 
to indicate whether the reaction ids of adjacent cells differs or not. You can do this value immediately after the assignment of reaction ids for cells.
blackmask is offline   Reply With Quote

Old   July 15, 2013, 04:02
Default
  #13
Senior Member
 
Ananthakrishnan.A.S
Join Date: Feb 2012
Location: Mumbai (Bombay), India
Posts: 157
Rep Power: 8
Ananthakrishnan is on a distinguished road
Send a message via Skype™ to Ananthakrishnan
Hi,

Thanks a lot for the reply..
I tried to work out the logic like you mentioned but couldnt get going. Sorry about my ignorance but once again i will really happy if you can give me an idea similar to switch program!!!

and regarding the switch case option..
i am not using the switch anymore. I am using somthing like this below

#define MAX_REACTION_ID 50
DEFINE_ON_DEMAND(report_reaction_region_averages)
{
cell_t c;
Thread *ct;
Domain *d;
int rid;
real tot_temp[MAX_REACTION_ID];
real tot_vol[MAX_REACTION_ID];
real tot_mass[MAX_REACTION_ID];
real mass, vol;

d = Get_Domain(ROOT_DOMAIN_ID);
/* Initialise totals */
for (rid=0;rid<MAX_REACTION_ID;rid++)
{
tot_temp[rid] = 0.0;
tot_volume[rid] = 0.0;
tot_mass[rid] = 0.0;
}
thread_loop_c(ct,d)
if(FLUID_THREAD_P(ct))
begin_c_loop(c,ct)
{
rid = C_CHEM_REACTION_ID(c,ct)
if(rid < MAX_REACTION_ID)
{
vol = C_VOLUME(c,ct);
mass = vol*C_R(c,ct);
tot_vol += vol;
tot_mass += mass;
tot_temp[rid] += mass * C_T(c,ct);
}
}
end_c_loop(c,ct)
/* Report averages */
for (rid=0;rid<MAX_REACTION_ID;rid++)
{
if(tot_mass[rid] > 0.0)
{
Message("Average Temperature for Reaction %d is %f\n",i,tot_temp[rid]/tot_mass[rid]);
}
else
Message("Reaction %d does not occur anywhere\n");
}
}

Thanks,
Ananthakrishnan is offline   Reply With Quote

Old   July 15, 2013, 19:49
Default
  #14
Senior Member
 
Join Date: Aug 2011
Posts: 315
Rep Power: 12
blackmask will become famous soon enough
Suppose that id0, id1 < 100, what I am trying to do is to store the combination of id0,id1 in someway and retrieve them later,
Code:
#include <math.h>
#include "udf.h"

thread_loop_f(ft,domain)
{
        if(THREAD_TYPE(ft) == THREAD_F_INTERIOR)
                begin_f_loop(f,ft)
                {

                        id0 = C_CHEM_REACTION_ID(F_C0(f,ft), THREAD_T0(ft));

                        id1 = C_CHEM_REACTION_ID(F_C1(f,ft), THREAD_T1(ft));

                        C_CHEM_REACTION_ID(f, ft) = id0*100+id1;

                }   
        end_f_loop(f,ft)
                ;   
}

thread_loop_f(ft,domain)
{
        if(THREAD_TYPE(ft) == THREAD_F_INTERIOR)
                begin_f_loop(f,ft)
                {

                        id01 = C_CHEM_REACTION_ID(f, ft);

                        id0 = lround(id01/100);

                        id1 = lround(id01 - 100*id0);


                        if (id0 != id1)
                        {   
                                mass_FLUX += F_FLUX(ft,ct);
                        }
                        //
                }   
        end_f_loop(f,ft)
                ;   
}
I assume the storage for UDM in thread "ft" has been allocated so that I did not check it.
blackmask is offline   Reply With Quote

Reply

Thread Tools
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 On
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
No layers in a small gap bobburnquist OpenFOAM Native Meshers: snappyHexMesh and Others 6 August 26, 2015 09:38
SnappyHexMesh for internal Flow vishwa OpenFOAM Native Meshers: snappyHexMesh and Others 23 August 6, 2014 03:50
Import netgen mesh to OpenFOAM hsieh Open Source Meshers: Gmsh, Netgen, CGNS, ... 32 September 13, 2011 05:50
NACA0012 geometry/design software needed Franny Main CFD Forum 13 July 7, 2007 15:57
physical boundary error!! kris CD-adapco 2 August 3, 2005 00:32


All times are GMT -4. The time now is 18:34.