[sooraj546]

I have an annular disc whose radius goes from 1.3 m to 2.7 m with 48 mm thickness and angle subtended 316 degree. I need to write a UDF for heat generation as a function of angle. I will use the DEFINE_SOURCE function but I am not able to understand how to in corporate the angle variation in disc without the radius.

Rotation angle of the hearth Heat generation rate (kW/m3 )
50 0
100 0.28
125 0.285
150 0.25
175 0.21
200 0.12
225 0.07
250 0.03
275 0.01
300 0.01
from this data i have formulated the equation by curve fitting

Q = (0.0001*x^5) - (0.0044*x^4) + (0.0565*x^3) - (0.349*x^2) + (0.9717*x) - 0.6707 this is the variation of the heat source with angular change Q is the heat generated x is the angle in radians.




for(x=0;x<317;x++)
{
Q = (0.0001*x^5) - (0.0044*x^4) + (0.0565*x^3) - (0.349*x^2) + (0.9717*x) - 0.6707
}

how to covert from Cartesian to cylindrical co-ordinate i have no idea
Kindly help, I have searched a lot but not able to figure out the solution.

[obscureed]

Hi Sooraj546,

A few comments:
(1) The polynomial you have supplied is not a good fit when you supply the polynomial coefficients at low precision. (Try it in a spreadsheet!) If you are going to use a polynomial fit, you will need to supply many decimal places for each coefficient. You should also check that the polynomial fit gives good results when the angle is lower than 50deg or higher than 300deg.
(2) Note that the source term in Fluent needs to be in [W/m3]. (I steer clear of the density-based solver when using energy sources.)
(3) You want to know "how to covert from Cartesian to cylindrical co-ordinate" -- have you tried a Google search for those words? This should help a lot. To program this, one important trick is to use atan2(y,x) instead of atan(y/x) -- atan2 deals with all the combinations of positive, negative and zero for x and y. See https://www.tutorialspoint.com/c_sta...tion_atan2.htm for example. (If you try to test the calculations in Excel, be aware that Excel's atan2 has arguments in alphabetical order, atan2(x,y), unlike almost any other implementation.) Note that atan2 returns an angle in the range [-M_PI,M_PI], so you should probably add (2*M_PI) to negative values.
(4) You will be evaluating the source term cell by cell in a UDF of type DEFINE_SOURCE, so there is no need for the loop "for(x=0;x<317;x++)" that you mention (which uses degrees, by the way!). You can obtain the vector of coordinates x using the command "C_CENTROID(x,c,t);".

Good luck! Ed

[sooraj546]

thanks a lot sir, i was searching for it in UDF if any macro etc are there .

[sooraj546]

So if i change the co-ordinate system from Cartesian to cylindrical system in meshing part as shown in the attached diagram

then apply this for loop for angular variation of heat source will it work

for(x=0;x<317;x++), x is angular displacement in radians
{
Q = f(x)
}


will it work??????

[obscureed]

Adding that coordinate system to the meshing set-up will not achieve anything in Fluent. (Also, try to get an update if you can -- the software at version 18 is genuinely better than version 15. Version 19 is coming soon.) I would advise you to go back to the plan of looking up Cartesian coordinates using C_CENTROID, and then calculating an angle from something like "angle = atan2(x[1],x[0]); if(angle<0.) { angle += 2.*M_PI; }"

Regarding the loop, see my previous point (4). If you have a loop command in your DEFINE_SOURCE, you are probably doing it wrong (unless, for example, you decide to evaluate your source term by linear interpolation on the data, which would not be a bad idea).

[sooraj546]

THANKS A LOT SIR.
I will do as per ur advise and c whats happening

[sooraj546]

In the attachment there is a pellet layer shown, i want to attach this UDF on that layer its a 3-D annular disc like zone. so this code will be ok
geometry-RHF1.png

#include "udf.h"

#define pi 3.14159265

DEFINE_SOURCE(CO_source, cell, thread, dS, eqn)
{

real x[ND_ND];
real angle, source;

C_CENTROID(x, cell, thread);
angle = atan2(x[1],x[0]);/* in rad*/

if(angle<0.)
{
Nangle = angle+(2.*PI);
source = -(6E-05*(pow(angle,6)))+(0.002*(pow(angle,5)))-(0.0279*(pow(angle,4)))+(0.2032*(pow(angle,3)))-(0.8114*(pow(angle,2)))+(1.6473*angle)-1.012;
}
else
{
source = -(6E-05*(pow(angle,6)))+(0.002*(pow(angle,5)))-(0.0279*(pow(angle,4)))+(0.2032*(pow(angle,3)))-(0.8114*(pow(angle,2)))+(1.6473*angle)-1.012;
}
dS[eqn] = 0.0;
return source;
}

[obscureed]

This looks like a good start, but you have some debugging to do and then some geometry. Some hints:
(a) "Nangle"
(b) When you come to choose between pi, PI and M_PI, choose the one that is built in.
(c) You need to work out for yourself which values (either 0, 1 or 2) to put in "angle = atan2(x[???],x[???]);" -- this depends on the geometry. This should also prompt you to consider where angle equals 0, and which direction angle sweeps round. You will probably need to apply some kind of offset to make it start where you want.
(d) Rearrange your "if" and brackets so that you only type the polynomial once.
(e) If you insist on using a polynomial, a more efficient way to evaluate it is ((((((-6E-05*angle+0.002)*angle-0.0279)*angle+0.2032)*angle-0.8114)*angle+1.6473)*angle-1.012). (https://en.wikipedia.org/wiki/Horner%27s_method) For style points, define the coefficients in a const array and loop through them -- it will make it easier to put in new coefficients without typing errors.
(f) But, seriously, look again at my earlier point (1) -- try it in a spreadsheet! I can't even tell whether you've acted on point (2).

[sooraj546]

Sir, thanks a lot for all your suggestions

replies

point-1 -->you are absolutely correct sir. i have tried it in spread sheet its coming somewhat close to the values but not accurate, i guess its fitting error. Before 50 and after 300 i have changed the code like shown below.

point-2-->When you come to choose between pi, PI and M_PI, choose the one that is built in-------i cant find the representation of pi in ansys i searched a lot as u said i believe its M_PI as shown in the code.

point-3-->Note that the source term in Fluent needs to be in [W/m3]. (I steer clear of the density-based solver when using energy sources.)---------yes sir all the values of y are already in W/m3 and x in degrees -----------so obviously y = f(x) will be in w/m3 i guess and i am using a pressure based solver

point-4-->i am changing this statement ------angle = atan2(x[2],x[0]); i am attaching a the 3D geometry with this ---------the 3D geometry with origin as centre and angular variation along z and x axis thats y atan2(x[2],x[0])? am i correct sir

I have a doubt sir-if i rotate the geometry and places its one edge of the yellow colored zone (shown in attached diagram) on X-axis where the heat source is to be incorporated --------will that start measuring the angle from x-axis in clockwise or anti-clock wise direction ? or its measuring along the geometry as from cell to cell along cell thread etc-----------my geometry has about 316 span so that first will be zero as its on x-axis etc

According to my exp------ measurements are done in anticlockwise direction i.e the progression of heat source from 0 to 316 along the yellow zone is in anticlockwise 0-50 is the preheat zone etc.

#include "udf.h"
#include "mem.h"

DEFINE_SOURCE(Heat_source, cell, thread, dS, eqn)
{
real M_PI
real x[ND_ND];
real angle;
real source;
real Nangle;
source = 0.0;
C_CENTROID(x, cell, thread);
angle = atan2(x[2],x[0]);

if(angle<=0.)
{
angle = (angle+(2.*M_PI));

if (angle<=0.872664626)
{
source=0;
}
else if (angle>=0.872664626) && angle<=5.235987756)
{
source = ((((((354.9*angle-12576)*angle+178009)*angle-pow(10,6))*angle+
(5*pow(10,6)))*angle-(pow(10,7)))*angle)+(6*pow(10,6));
}
else (angle>5.235987756)
;{
source=0;
}
}

else if (angle>0 && angle<=0.872664626)
{
source=0;
}
else if (angle>=0.872664626 && angle<=5.235987756)
{
source = ((((((354.9*angle-12576)*angle+178009)*angle-pow(10,6))*angle+
(5*pow(10,6)))*angle-(pow(10,7)))*angle)+(6*pow(10,6));
}
else (angle>5.235987756)
;{
source=0;
}

dS[eqn] = 0.0;

return source;
}


THANKS A LOT IN ADVANCE

[obscureed]

You have got this far, so I think you should be able to work out for yourself where the angles are. A good principle is not to trust your own calculations, but to actually see the correct inputs in place. (If necessary, store the value of the source term and some intermediate calculation results in User-Defined Memories, and plot contours to show that they are as expected.)

Some last comments (in increasing order of severity), and then you are on your own:
-- You have fixed my point (a) in one respect, but made it worse in another.
-- See my point (d) again.
-- To replace 6e6 with "6*pow(10,6)" is pointless (and slow, and difficult to read). C understands 6e6.
-- "else (angle > 300)" is wrong. "; {" is also wrong -- I suspect that you added semicolons until it compiled without errors. This is not a good way to make C do what you intend.

The really important one is point (1). You might be shocked at how bad that polynomial is ***when some of the coefficients are defined with one significant digit of precision***. When you showed the spreadsheet generating a nice curve, the spreadsheet was calculating the coefficients and using them with double precision. What I wanted you to do was to evaluate the polynomial using the the coefficients as displayed (as you are proposing to do in the UDF). This is a massive error. Try it on a calculator! I could show you a graph, but the polynomial is so enormously positive that the small negative values would not be visible below the axis. For example, at angle=300, the answer is intended to be approximately zero, but the polynomial as stated equals +2.3e17. Polynomials are not robust fitting functions.

[sooraj546]

yes sir u r absolutely correct, thanks a lot for that. I will adapt the polynomial accordingly

[sooraj546]

Dear sir,

i have changed the code accordingly except the polynomial part which i will change too but my concern here is now this code is getting compiled but floating point error is coming i dont no why - without UDF the problem is running fine and i am getting good results as well. floating point error is mainly associated with a boundary condition errors. In my case the boundary condition before and after attaching are the same so it doesnt make sence . Kindly help me out sir

#include "udf.h"
#include "mem.h"

DEFINE_SOURCE(CO_source, cell, thread, dS, eqn)
{

real x[ND_ND];
real angle,c,source;

C_CENTROID(x, cell, thread);

angle = atan2(x[2],x[0]);

if(angle<=0.)
{
angle = (angle+(2.*M_PI));
}

if (angle>=0.872664626 && angle<=5.235987756)
{
c = (angle*180)/M_PI;
source = (((((-(6e-14)*c+(6e-11))*c-(3e-8))*c+(6e-6))*c-(0.0008))*c+(0.0593))*c-(1.5632);

}
else
{
source=0;
}

dS[eqn] = 0.0;
return source;
}

[obscureed]

Perhaps I have not explained this forcefully enough: the error in the polynomial is huge: literally, several orders of magnitude. Until you fix this, it is not worth looking for any other issues.

[sooraj546]

Oh Understood sir thanks a lot i will rectify it and get back to you