CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > Software User Forums > OpenFOAM

a problem with wallGradT

Register Blogs Community New Posts Updated Threads Search

Reply
 
LinkBack Thread Tools Search this Thread Display Modes
Old   December 9, 2013, 14:00
Default wallGradT and wrong answer
  #1
Member
 
Arjang Behnoud
Join Date: Oct 2012
Posts: 63
Rep Power: 14
AmirBaqa1987 is on a distinguished road
Hi everyone

I want to get normal gradient of temperature at wall boundary in laminar incompressible flow.I've modified wallGradU utility as someone in other posts said.

I want to know why wallGradT remains fixed in some cells?and how I can make it correct?

my geometry is pipe with 1 meter length and 0.003 meter radius and the blockMeshDict is like:

Code:
convertToMeters 1;
vertices
(
    (0 0 0)      //0
    (0 0 0.001)
    (0 0.001 0.001)      //2
    (0 0.001 0)
    (0 0 0.003)   //4
    (0 0.0021213 0.0021213)
    (0 0.003 0 )      //6

    (1 0 0)      //7
    (1 0 0.001)
    (1 0.001 0.001)      //9
    (1 0.001 0)
    (1 0 0.003)   //11
    (1 0.0021213 0.0021213)
    (1 0.003 0 )      //13
    

);

blocks
(
    // inlet block
    hex (7 8 9 10 0 1 2 3) (9 9 300) simpleGrading (1 1 1)
    hex (8 11 12 9 1 4 5 2) (51 9 300) simpleGrading (0.03075 1 1)
    hex (10 9 12 13 3 2 5 6) (9 51 300) simpleGrading (1 0.03075 1)
);

edges
(
    arc 4 5 (0 0.001148 0.00277)
    arc 5 6 (0 0.00277 0.001148)
    arc 11 12 (1 0.001148 0.00277)
    arc 12 13 (1 0.00277 0.001148)
);
boundary
(
    inlet
    {
        type patch;
        faces
        (
            (0 1 2 3)
            (1 4 5 2)
            (2 5 6 3)
        );
    }

    outlet
    {
        type patch;
        faces
        (
            (7 10 9 8)
            (8 9 12 11)
            (9 10 13 12)
        );
    }

    side1
    {
        type symmetryPlane;
        neighbourPatch side2;
        faces
        (
            (0 7 8 1)
            (1 8 11 4) 
        );

       transform rotational;
        rotationAxis (1 0 0);
        rotationCentre (0 0 0);
    }

    side2
    {
        type symmetryPlane;
        neighbourPatch side1;
        faces
        (
            (0 3 10 7)
            (3 6 13 10)
        );

        transform rotational;
        rotationAxis (1 0 0);
        rotationCentre (0 0 0);
    }

    walls
    {
        type wall;
        faces
        (
            (4 11 12 5)
           (5 12 13 6)
        );
    }

   
);
and the boundary condition for T is like:
Code:
dimensions      [0 0 0 1 0 0 0];

internalField   uniform 293;

boundaryField
{
    inlet
    {
        type            fixedValue;
        value           uniform 293;
    }
    outlet
    {
        type            zeroGradient;
    }
    
    walls
    {
        type            fixedValue;
        value           uniform 373;
    }

  side2
    {
        type            symmetryPlane;
        
    }
    side1
    {
        type            symmetryPlane;
        
    }  
}
and the fvScheme is like:
Code:
ddtSchemes
{
    default         steadyState;
}

gradSchemes
{
    default         Gauss linear;
    grad(p)         Gauss linear;
    grad(U)         Gauss linear;
}

divSchemes
{
    default         none;
    div(phi,U)      Gauss limitedLinearV 1;
    div(phi,k)      Gauss limitedLinear 1;
    div(phi,epsilon) Gauss limitedLinear 1;
    div(phi,R)      Gauss limitedLinear 1;
    div(R)          Gauss linear;
    div(phi,nuTilda) Gauss limitedLinear 1;
    div((nuEff*dev(T(grad(U))))) Gauss linear;
    div(phi,T)    Gauss linear;
   
}

laplacianSchemes
{
    default         none;
    laplacian(nuEff,U) Gauss linear corrected;
    laplacian((1|A(U)),p) Gauss linear corrected;
    laplacian(DkEff,k) Gauss linear corrected;
    laplacian(DepsilonEff,epsilon) Gauss linear corrected;
    laplacian(DREff,R) Gauss linear corrected;
    laplacian(DnuTildaEff,nuTilda) Gauss linear corrected;
    laplacian(alpha,T)                 Gauss linear corrected;
    laplacian(1,p)                 Gauss linear corrected;
    
}

interpolationSchemes
{
    default         linear;
    interpolate(U)  linear;
}

snGradSchemes
{
    default         Uncorrected;
}

fluxRequired
{
    default         no;
    p               ;
}


// ************************************************************************* //
by considering the above conditions the amount of wallGradT should be constantly decreased through pipe because the amount of the temperature must be increased constantly . but as you can see in the attached files the amount of wallGradT in some cells remains fixed and do not decrease .

can any body help me?
thanks.
arjang
Attached Files
File Type: xls Untitled 1.xls (19.0 KB, 22 views)

Last edited by AmirBaqa1987; December 10, 2013 at 07:10.
AmirBaqa1987 is offline   Reply With Quote

Old   December 10, 2013, 02:45
Default
  #2
Member
 
Arjang Behnoud
Join Date: Oct 2012
Posts: 63
Rep Power: 14
AmirBaqa1987 is on a distinguished road
Nobody can help?
AmirBaqa1987 is offline   Reply With Quote

Old   December 10, 2013, 17:33
Default
  #3
Retired Super Moderator
 
Bruno Santos
Join Date: Mar 2009
Location: Lisbon, Portugal
Posts: 10,982
Blog Entries: 45
Rep Power: 128
wyldckat is a name known to allwyldckat is a name known to allwyldckat is a name known to allwyldckat is a name known to allwyldckat is a name known to allwyldckat is a name known to all
Greetings Arjang,

Well, knowing which posts you are referring to would help, in order to assess if something might have gone wrong.

In addition, knowing the "U" and "p" boundary conditions would also help, as well as knowing which solver you have used and if has converged.

And a picture of what your are seeing in ParaView would also help.

Best regards,
Bruno
__________________
wyldckat is offline   Reply With Quote

Old   December 11, 2013, 03:43
Default
  #4
Member
 
Arjang Behnoud
Join Date: Oct 2012
Posts: 63
Rep Power: 14
AmirBaqa1987 is on a distinguished road
Quote:
Originally Posted by wyldckat View Post
Greetings Arjang,

Well, knowing which posts you are referring to would help, in order to assess if something might have gone wrong.

In addition, knowing the "U" and "p" boundary conditions would also help, as well as knowing which solver you have used and if has converged.

And a picture of what your are seeing in ParaView would also help.

Best regards,
Bruno
dear Bruno
thanks for your attention.
I'm want help about wallGradT.

Nu = wallGradT*D/(Twall - Tbulk)

for calculating Nusselt number I've calculated bulk Temperature and also normal gradient of temperature at wall boundary.
as you can see in the attached file,the amount of wallGradT in some cells remains fixed but the amount of Tbulk has increased continuously and this fact caused the Nusselt diagram to have oscillation.(please look at attached files) but It must not have this behaviour.
I've used newSimpleFoam solver which I've added energy equation to. and the solution has converged( but the manner of convergence of T equation is not good,I mean the initial residuals of T did not constantly decrease).

TEqn:

Code:
fvScalarMatrix TEqn
 (
    fvm::div(phi,T)
   -fvm::laplacian(alpha,T)
 );

TEqn.solve();
conditions for U ;
Code:
dimensions      [0 1 -1 0 0 0 0];

internalField   uniform (0 0 0);

boundaryField
{
    inlet
    {
        type            fixedValue;
        value           uniform (0.10088 0 0);
    }
    outlet
    {
        type            zeroGradient;
    }
   side1
    {
        type           symmetryPlane;
       
    }
    side2
    {
        type           symmetryPlane;
       
    }
    walls
    {
        type            fixedValue;
        value           uniform (0 0 0);
    }
 
}
and for P:
Code:
dimensions      [0 2 -2 0 0 0 0];

internalField   uniform 0;

boundaryField
{
    inlet
    {
        type            zeroGradient;
    }
    outlet
    {
        type            fixedValue;
        value           uniform 0;
    }
    
    walls
    {
        type            zeroGradient;
    }
  side2
    {
        type            symmetryPlane;
        
    }
    side1
    {
        type            symmetryPlane;
        
    }
}
Thanks again.
sincerely yours,
Arjang

wallGradT&Tbulk.xls

Nu.jpg
Attached Images
File Type: jpg Tbulk.jpg (28.9 KB, 30 views)
AmirBaqa1987 is offline   Reply With Quote

Old   December 15, 2013, 11:44
Default
  #5
Retired Super Moderator
 
Bruno Santos
Join Date: Mar 2009
Location: Lisbon, Portugal
Posts: 10,982
Blog Entries: 45
Rep Power: 128
wyldckat is a name known to allwyldckat is a name known to allwyldckat is a name known to allwyldckat is a name known to allwyldckat is a name known to allwyldckat is a name known to all
Hi Arjang,

Is there any chance you can share the complete case, instead of only parts of it?
I ask this because at least the following information is still missing:
  1. What are the fluid properties being used?
  2. What is the Reynolds number at work here? More specifically, are you 100% certain the flow is laminar? I say this because 0.1 m/s vs 0.003 m of radius in a pipe, does not feel very laminar to me, so it really depends on the fluid in question.
  3. How exactly are you sampling these data points? I ask this because there are only 101 points and there are 300 cells along the pipe.
    • Attached is an image that depicts the sampling positions at 0.57 and 0.58 m (2 consecutive points where the values are identical) and they are placed in very complicated locations.
    • The locations are complicated, because they are located right between 2 cells each, which means that the calculated value is in fact interpolated from the centers of the cells in question, which can introduce some numerical error. Wait, actually, they are located between 4 cells, which can introduce some additional error.
  4. How or where exactly does "wallGradT" come from?
Please provide as much information as possible, as explained here: http://www.cfd-online.com/Forums/ope...-get-help.html - otherwise this is just a guessing game and people here on the forum don't have time to play games


Because right now, my guess is that this is a sampling problem (insufficient data points or wrongly located points), or the case is incorrectly prepared.


Best regards,
Bruno
Attached Images
File Type: jpg Sampling positions at 0.57 and 0.58.jpg (20.5 KB, 27 views)
__________________
wyldckat is offline   Reply With Quote

Old   December 16, 2013, 05:21
Default
  #6
Member
 
Arjang Behnoud
Join Date: Oct 2012
Posts: 63
Rep Power: 14
AmirBaqa1987 is on a distinguished road
Hi Dear Bruno
yes why not.

and I've developed solver like below:
createFields
Code:
 Info<< "Reading transportProperties\n" << endl;

    IOdictionary transportProperties
    (
        IOobject
        (
            "transportProperties",
            runTime.constant(),
            mesh,
            IOobject::MUST_READ_IF_MODIFIED,
            IOobject::NO_WRITE
        )
    );

   Info<< "Reading field p\n" << endl;
    volScalarField p
    (
        IOobject
        (
            "p",
            runTime.timeName(),
            mesh,
            IOobject::MUST_READ,
            IOobject::AUTO_WRITE
        ),
        mesh
    );

    Info<< "Reading field U\n" << endl;
    volVectorField U
    (
        IOobject
        (
            "U",
            runTime.timeName(),
            mesh,
            IOobject::MUST_READ,
            IOobject::AUTO_WRITE
        ),
        mesh
    );

    Info<< "Reading field T\n" << endl;
    volScalarField T
    (
        IOobject
        (
            "T",
            runTime.timeName(),
            mesh,
            IOobject::MUST_READ,
            IOobject::AUTO_WRITE
        ),
        mesh
    ); 

     dimensionedScalar nu
    (
        transportProperties.lookup("nu")
    );

   
    dimensionedScalar rho
    (
        transportProperties.lookup("rho")
    );

    
    dimensionedScalar c
    (
        transportProperties.lookup("c")
    );

    dimensionedScalar alpha
    (
        transportProperties.lookup("alpha")
    );

    
    dimensionedScalar k
    (
        transportProperties.lookup("k")
    );

    #include "createPhi.H"


    label pRefCell = 0;
    scalar pRefValue = 0.0;
    setRefCell(p, mesh.solutionDict().subDict("SIMPLE"), pRefCell, pRefValue);

    singlePhaseTransportModel laminarTransport(U, phi);

    autoPtr<incompressible::RASModel> turbulence
    (
        incompressible::RASModel::New(U, phi, laminarTransport)
    );

    IObasicSourceList sources(mesh);
TEqn.h
Code:
fvScalarMatrix TEqn
 (
   (fvm::div(phi,T)
  -fvm::Sp(fvc::div(phi),T))
  -fvm::laplacian(alpha,T)
 );

TEqn.solve();
heatSimpleFoam.c
Code:
int main(int argc, char *argv[])
{
    #include "setRootCase.H"
    #include "createTime.H"
    #include "createMesh.H"
    #include "createFields.H"
    #include "initContinuityErrs.H"

    simpleControl simple(mesh);

    // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

    Info<< "\nStarting time loop\n" << endl;

    while (simple.loop())
    {
        Info<< "Time = " << runTime.timeName() << nl << endl;

        // --- Pressure-velocity SIMPLE corrector
        {
            #include "UEqn.H"
            #include "pEqn.H"
            #include "TEqn.H"
        }

        turbulence->correct();

        runTime.write();

        Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
            << "  ClockTime = " << runTime.elapsedClockTime() << " s"
            << nl << endl;
    }

    Info<< "End\n" << endl;

    return 0;
}
wallGradT code
Code:
int main(int argc, char *argv[])
{
    timeSelector::addOptions();
    #include "setRootCase.H"
#   include "createTime.H"
    instantList timeDirs = timeSelector::select0(runTime, args);
#   include "createMesh.H"

    forAll(timeDirs, timeI)
    {
        runTime.setTime(timeDirs[timeI], timeI);
        Info<< "Time = " << runTime.timeName() << endl;

        IOobject Theader
        (
            "T",
            runTime.timeName(),
            mesh,
            IOobject::MUST_READ
        );

        // Check T exists
        if (Theader.headerOk())
        {
            mesh.readUpdate();

            Info<< "    Reading T" << endl;
            volScalarField T(Theader, mesh);

            Info<< "    Calculating wallGradT" << endl;

            volScalarField wallGradT
            (
                IOobject
                (
                    "wallGradT",
                    runTime.timeName(),
                    mesh,
                    IOobject::NO_READ,
                    IOobject::AUTO_WRITE
                ),
                mesh,
                dimensionedScalar
                (
                    "wallGradT",
                    T.dimensions()/dimLength,
                    0
                )
            );

            forAll(wallGradT.boundaryField(), patchi)
            {
                wallGradT.boundaryField()[patchi] =
                    mag(T.boundaryField()[patchi].snGrad());
            }

            wallGradT.write();
        }
        else
        {
            Info<< "    No T" << endl;
        }
    }

    Info<< "End" << endl;

    return 0;
}
- Re number here is 600.
- I've plot these data over line in paraFoam with following coordination:
point1: (0 0 0.003)
point2: (1 0 0.003)

and the resolution is 100 (this cause 101 points being plotted).

Thanks
Best Regards.
Attached Images
File Type: jpg geometry.jpg (79.3 KB, 19 views)
Attached Files
File Type: gz complete case.tar.gz (5.5 KB, 12 views)

Last edited by AmirBaqa1987; December 16, 2013 at 06:33.
AmirBaqa1987 is offline   Reply With Quote

Old   December 16, 2013, 07:39
Default
  #7
Super Moderator
 
Tobi's Avatar
 
Tobias Holzmann
Join Date: Oct 2010
Location: Bad Wörishofen
Posts: 2,711
Blog Entries: 6
Rep Power: 52
Tobi has a spectacular aura aboutTobi has a spectacular aura aboutTobi has a spectacular aura about
Send a message via ICQ to Tobi Send a message via Skype™ to Tobi
Hello Arjang,

just a few question and hints.

  1. You are using the SIMPLE algorithm. You mixed up some fvSolution files, aren´t you?
  2. In my opinion you should use PBiCG instead of BICCG
  3. Use PBiCG instead of smoothSolver
  4. You use a underrelaxation factor for T but I hope that you know that this does not be used in your calculation due to the fact, that your are not using a relaxation for your TEqn like:
Code:
fvScalarMatrix TEqn
 (
   (fvm::div(phi,T)
  -fvm::Sp(fvc::div(phi),T))   // Red brackets can be removed
  -fvm::laplacian(alpha,T)
 );
 
TEqn.solve();
 
TEqn.relax();          // therefor your Relaxation factor is used
5. As you mentioned that your case has no smooth convergence for T i think this is a manner of the not used underrelaxation in your equation.
6. I think your case stops if you reach your very low residual control Limits:
Code:
{
p 1e-2;        // 1e-5
U 1e-3;        // 1e-6
T 1e-3;        // 1e-6
}

I think if you Change your solvers and add the Relaxation to your TEqn you will get accurate and better results.
Try it and let us know.

Regards Tobi
Tobi is offline   Reply With Quote

Old   December 18, 2013, 04:02
Default
  #8
Member
 
Arjang Behnoud
Join Date: Oct 2012
Posts: 63
Rep Power: 14
AmirBaqa1987 is on a distinguished road
Hi Tobias
thanks for your Tips.

I did your suggestions I mean:
1.using PBiCG instead of BICCG for T.
2.adding TEqn.relax(); inTEqn.H

now the convergence of T is smooth but the oscillations in Nu diagram exists yet.

Danke Schon .
Arjang
AmirBaqa1987 is offline   Reply With Quote

Old   December 19, 2013, 07:52
Default
  #9
Super Moderator
 
Tobi's Avatar
 
Tobias Holzmann
Join Date: Oct 2010
Location: Bad Wörishofen
Posts: 2,711
Blog Entries: 6
Rep Power: 52
Tobi has a spectacular aura aboutTobi has a spectacular aura aboutTobi has a spectacular aura about
Send a message via ICQ to Tobi Send a message via Skype™ to Tobi
Quote:
Originally Posted by wyldckat View Post
Hi Arjang,

Is there any chance you can share the complete case, instead of only parts of it?
I ask this because at least the following information is still missing:

  1. How exactly are you sampling these data points? I ask this because there are only 101 points and there are 300 cells along the pipe.
    • Attached is an image that depicts the sampling positions at 0.57 and 0.58 m (2 consecutive points where the values are identical) and they are placed in very complicated locations.
    • The locations are complicated, because they are located right between 2 cells each, which means that the calculated value is in fact interpolated from the centers of the cells in question, which can introduce some numerical error. Wait, actually, they are located between 4 cells, which can introduce some additional error.
Because right now, my guess is that this is a sampling problem (insufficient data points or wrongly located points), or the case is incorrectly prepared.


Best regards,
Bruno

Hi,

did you check the hints Bruno mentioned?
At the moment it seems that it is a sampling Problem like he said.

Additionally. Check the regions with These "Nu" oscillations. I think Bruno is telling you the correct Problem because if you interpolate the massflow of an inlet and outlet you will not get the same values (Interpolation Problem with paraview).


Regards
Tobi
Tobi is offline   Reply With Quote

Reply

Tags
eshdict


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


Similar Threads
Thread Thread Starter Forum Replies Last Post
UDF compiling problem Wouter Fluent UDF and Scheme Programming 6 June 6, 2012 05:43
Gambit - meshing over airfoil wrapping (?) problem JFDC FLUENT 1 July 11, 2011 06:59
natural convection problem for a CHT problem Se-Hee CFX 2 June 10, 2007 07:29
Adiabatic and Rotating wall (Convection problem) ParodDav CFX 5 April 29, 2007 20:13
Is this problem well posed? Thomas P. Abraham Main CFD Forum 5 September 8, 1999 15:52


All times are GMT -4. The time now is 00:16.