[hiwactb]

Hi all!

I'm doing simulations at supercritical pressures. As the specific heat changes dramatically with temperature, I have to use UDF by defining the macro named define_specific_heat. However, i cant get accurate results as the temperature hardly changed along the flow direction. I don't know why and wonder if the udf for specific heat is incorrect. Does anyone can help me check my code? Thank you very much.

The code im using is as follows:


#include"udf.h"
DEFINE_SPECIFIC_HEAT(superfluid_cp, T, Tref, h, yi)
{
real cp;
if (369 <= T < 430)
cp = 3*pow(10, -5)*(T-369)*(T-369)*(T-369)+2*pow(10, -4)*(T-369)*(T-369)+3.8841*(T-369)+2235.7;
elseif (430 <= T < 471)
cp = 0.00001*(T-430)*(T-430)*(T-430)+0.0003*(T-430)*(T-430)+4.1015*(T-430)+2479.4;
elseif (471 <= T < 507)
cp = 1*pow(10, -6)*(T-471)*(T-471)*(T-471)*(T-471)-0.0002*(T-471)*(T-471)*(T-471)+0.0121*(T-471)*(T-471)+3.9779*(T-471)+2648.6;
elseif (507 <= T < 559)
cp = 1*pow(10, -6)*pow((T-507),4)-0.00006*(T-507)*(T-507)*(T-507)-0.0005*(T-507)*(T-507)+4.3283*(T-507)+2800.8;
elseif (559 <= T < 598)
cp = 7*pow(10, -6)*(T-559)*(T-559)*(T-559)*(T-559)-0.00008*(T-559)*(T-559)*(T-559)+0.012*(T-559)*(T-559)+4.6167*(T-559)+3026.7;
elseif (598 <= T < 632)
cp = 2*pow(10, -7)*(T-598)*(T-598)*(T-598)*(T-598)*(T-598)+2*pow(10, -6)*(T-598)*(T-598)*(T-598)*(T-598)+0.0005*(T-598)*(T-598)*(T-598)+0.0354*(T-598)*(T-598)+6.305*(T-598)+3233.3;
elseif (632 <= T < 643)
cp = 0.0003*(T-632)*(T-632)*(T-632)*(T-632)+0.003*(T-632)*(T-632)*(T-632)+0.227*(T-632)*(T-632)+12.47*(T-632)+3522.3;
elseif (643 <= T < 654)
cp = 9*pow(10, -5)*(T-643)*(T-643)*(T-643)*(T-643)*(T-643)-0.0006*(T-643)*(T-643)*(T-643)*(T-643)+0.029*(T-643)*(T-643)*(T-643)+0.523*(T-643)*(T-643)+20.273*(T-643)+3695.5;
elseif (654 <= T < 665)
cp = 6*pow(10, -5)* pow(T-654, 6)-0.00005*pow(T-654, 5)+0.0065*pow(T-654, 4)+0.142*pow(T-654, 3)+2.2238*pow(T-654, 2)+46.029*(T-654)+4026.3;
elseif (665 <= T < 672)
cp = -1.2999*pow(T-665, 4)+9.5373*pow(T-665, 3)+5.6801*pow(T-665, 2)+256.28*(T-665)+5180.7;
elseif (672 <= T < 680)
cp = 0.0053*pow(T-672, 6)-0.1823*pow(T-672, 5)+2.0009*pow(T-672, 4)-2.6145*pow(T-672, 3)-77.333*pow(T-672, 2)+81.763*(T-672)+7453.3;
elseif (680 <= T < 692)
cp = -0.0007*pow(T-680, 5)+0.0443*pow(T-680, 4)-1.2394*pow(T-680, 3)+21.253*pow(T-680, 2)-256.32*(T-680)+5429.2;
elseif (692 <= T < 710)
cp = -0.00007*pow(T-692, 5)+0.0039*pow(T-692, 4)-0.1134*pow(T-692, 3)+2.5342*pow(T-692, 2)-49.579*(T-692)+4012.7;
elseif (710 <= T < 731)
cp = -0.000002*pow(T-710, 5)+0.0002*pow(T-710, 4)-0.0093*pow(T-710, 3)+0.3742*pow(T-710, 2)-12.225*(T-710)+3560.3;
elseif (731 <= T < 760)
cp = -0.00000008*pow(T-731, 5)+0.00002*pow(T-731, 4)-0.0015*pow(T-731, 3)+0.0993*pow(T-731, 2)-3.8505*(T-731)+3410.4;
elseif (760 <= T < 793)
cp = 0.000003*pow(T-760, 4)-0.0003*pow(T-760, 3)+0.0301*pow(T-760, 2)-0.5887*(T-760)+3355.4;
else
cp = 0.0000005*pow(T-793, 4)-0.0001*pow(T-793, 3)+0.0108*pow(T-793, 2)+0.6326*(T-793)+3359.5;

*h = cp * (T-Tref);
return cp;
}

Any help would be great. thanks again.

[mprinkey]

Does the UDF even compile? elseif is not C. The syntax is "else if"

http://www.cprogramming.com/tutorial/c/lesson2.html

I cut&pasted your code and turned it into a standalone program, changed the "elseif"s to "else if"s, and it seems to work:

#include <stdio.h>
#include <math.h>
#define real double
real specificHeat(real T)
{
real cp;
if (369 <= T < 430)
cp = 3*pow(10, -5)*(T-369)*(T-369)*(T-369)+2*pow(10, -4)*(T-369)*(T-369)+3.8841*(T-369)+2235.7;
else if (430 <= T < 471)
cp = 0.00001*(T-430)*(T-430)*(T-430)+0.0003*(T-430)*(T-430)+4.1015*(T-430)+2479.4;
else if (471 <= T < 507)
cp = 1*pow(10, -6)*(T-471)*(T-471)*(T-471)*(T-471)-0.0002*(T-471)*(T-471)*(T-471)+0.0121*(T-471)*(T-471)+3.9779*(T-471)+2648.6;
else if (507 <= T < 559)
cp = 1*pow(10, -6)*pow((T-507),4)-0.00006*(T-507)*(T-507)*(T-507)-0.0005*(T-507)*(T-507)+4.3283*(T-507)+2800.8;
else if (559 <= T < 598)
cp = 7*pow(10, -6)*(T-559)*(T-559)*(T-559)*(T-559)-0.00008*(T-559)*(T-559)*(T-559)+0.012*(T-559)*(T-559)+4.6167*(T-559)+3026.7;
else if (598 <= T < 632)
cp = 2*pow(10, -7)*(T-598)*(T-598)*(T-598)*(T-598)*(T-598)+2*pow(10, -6)*(T-598)*(T-598)*(T-598)*(T-598)+0.0005*(T-598)*(T-598)*(T-598)+0.0354*(T-598)*(T-598)+6.305*(T-598)+3233.3;
else if (632 <= T < 643)
cp = 0.0003*(T-632)*(T-632)*(T-632)*(T-632)+0.003*(T-632)*(T-632)*(T-632)+0.227*(T-632)*(T-632)+12.47*(T-632)+3522.3;
else if (643 <= T < 654)
cp = 9*pow(10, -5)*(T-643)*(T-643)*(T-643)*(T-643)*(T-643)-0.0006*(T-643)*(T-643)*(T-643)*(T-643)+0.029*(T-643)*(T-643)*(T-643)+0.523*(T-643)*(T-643)+20.273*(T-643)+3695.5;
else if (654 <= T < 665)
cp = 6*pow(10, -5)* pow(T-654, 6)-0.00005*pow(T-654, 5)+0.0065*pow(T-654, 4)+0.142*pow(T-654, 3)+2.2238*pow(T-654, 2)+46.029*(T-654)+4026.3;
else if (665 <= T < 672)
cp = -1.2999*pow(T-665, 4)+9.5373*pow(T-665, 3)+5.6801*pow(T-665, 2)+256.28*(T-665)+5180.7;
else if (672 <= T < 680)
cp = 0.0053*pow(T-672, 6)-0.1823*pow(T-672, 5)+2.0009*pow(T-672, 4)-2.6145*pow(T-672, 3)-77.333*pow(T-672, 2)+81.763*(T-672)+7453.3;
else if (680 <= T < 692)
cp = -0.0007*pow(T-680, 5)+0.0443*pow(T-680, 4)-1.2394*pow(T-680, 3)+21.253*pow(T-680, 2)-256.32*(T-680)+5429.2;
else if (692 <= T < 710)
cp = -0.00007*pow(T-692, 5)+0.0039*pow(T-692, 4)-0.1134*pow(T-692, 3)+2.5342*pow(T-692, 2)-49.579*(T-692)+4012.7;
else if (710 <= T < 731)
cp = -0.000002*pow(T-710, 5)+0.0002*pow(T-710, 4)-0.0093*pow(T-710, 3)+0.3742*pow(T-710, 2)-12.225*(T-710)+3560.3;
else if (731 <= T < 760)
cp = -0.00000008*pow(T-731, 5)+0.00002*pow(T-731, 4)-0.0015*pow(T-731, 3)+0.0993*pow(T-731, 2)-3.8505*(T-731)+3410.4;
else if (760 <= T < 793)
cp = 0.000003*pow(T-760, 4)-0.0003*pow(T-760, 3)+0.0301*pow(T-760, 2)-0.5887*(T-760)+3355.4;
else
cp = 0.0000005*pow(T-793, 4)-0.0001*pow(T-793, 3)+0.0108*pow(T-793, 2)+0.6326*(T-793)+3359.5;

//*h = cp * (T-Tref);
return cp;
}

int main (void)
{
    real T;
    for (T=370.0;T<790.0;T+=10)
        printf("T: %f, cp: %f\n",T,specificHeat(T));
    return 0;
}

Results:

mike@anubis:/r/home/mike/test> ./specificHeat 
T: 370.000000, cp: 2239.584330
T: 380.000000, cp: 2278.489230
T: 390.000000, cp: 2317.632130
T: 400.000000, cp: 2357.193030
T: 410.000000, cp: 2397.351930
T: 420.000000, cp: 2438.288830
T: 430.000000, cp: 2480.183730
T: 440.000000, cp: 2523.216630
T: 450.000000, cp: 2567.567530
T: 460.000000, cp: 2613.416430
T: 470.000000, cp: 2660.943330
T: 480.000000, cp: 2710.328230
T: 490.000000, cp: 2761.751130
T: 500.000000, cp: 2815.392030
T: 510.000000, cp: 2871.430930
T: 520.000000, cp: 2930.047830
T: 530.000000, cp: 2991.422730
T: 540.000000, cp: 3055.735630
T: 550.000000, cp: 3123.166530
T: 560.000000, cp: 3193.895430
T: 570.000000, cp: 3268.102330
T: 580.000000, cp: 3345.967230
T: 590.000000, cp: 3427.670130
T: 600.000000, cp: 3513.391030
T: 610.000000, cp: 3603.309930
T: 620.000000, cp: 3697.606830
T: 630.000000, cp: 3796.461730
T: 640.000000, cp: 3900.054630
T: 650.000000, cp: 4008.565530
T: 660.000000, cp: 4122.174430
T: 670.000000, cp: 4241.061330
T: 680.000000, cp: 4365.406230
T: 690.000000, cp: 4495.389130
T: 700.000000, cp: 4631.190030
T: 710.000000, cp: 4772.988930
T: 720.000000, cp: 4920.965830
T: 730.000000, cp: 5075.300730
T: 740.000000, cp: 5236.173630
T: 750.000000, cp: 5403.764530
T: 760.000000, cp: 5578.253430
T: 770.000000, cp: 5759.820330
T: 780.000000, cp: 5948.645230

[hiwactb]

Thank you very much. It is else if in my code, elseif is a typewriting mistake. I can build and load the udf for cp but the temperature variation along the tube is rather small. Maybe I should check other conditions like the boundary condition or the solution method. Thanks again.

[mprinkey]

Temperature-dependent specific heat is not going to INDUCE temperature variations. Your UDF will do the opposite, in fact. As you add energy to a cell, its temperature goes up. But the increasing value cp with temperature leads to a SMALLER temperature increase being associated with that energy increase versus a case with constant specific heat.