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

Adding temperature in icoFoam. Not working properly!

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

Reply
 
LinkBack Thread Tools Search this Thread Display Modes
Old   December 17, 2016, 17:56
Default Adding temperature in icoFoam. Not working properly!
  #1
New Member
 
Join Date: May 2016
Posts: 25
Rep Power: 10
mgab is on a distinguished road
Hi,
i modified the icoFoam solver to solve the temperature following this wiki:
https://openfoamwiki.net/index.php/H...ure_to_icoFoam
I run the cavity example and looks OK.

Now i'm trying to simulate a mixing elbow in which the two flows have different temperatures and the tube is adiabatic.
(I know icoFoam is only laminar but i'm learning this and then i'll move to simpleFoam to account the turbulence)

The simulation runs without errors but the results are incorrect (see attached image)

This is the solver:
Code:
/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
    \\  /    A nd           | Copyright (C) 2011-2016 OpenFOAM Foundation
     \\/     M anipulation  |
-------------------------------------------------------------------------------
License
    This file is part of OpenFOAM.

    OpenFOAM is free software: you can redistribute it and/or modify it
    under the terms of the GNU General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.

    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
    for more details.

    You should have received a copy of the GNU General Public License
    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.

Application
    icoFoamT

Description
    Transient solver for incompressible, laminar flow of Newtonian fluids with temperature field.

\*---------------------------------------------------------------------------*/

#include "fvCFD.H"
#include "pisoControl.H"

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

int main(int argc, char *argv[])
{
    #include "setRootCase.H"
    #include "createTime.H"
    #include "createMesh.H"

    pisoControl piso(mesh);

    #include "createFields.H"
    #include "initContinuityErrs.H"

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

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

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

        #include "CourantNo.H"

        // Momentum predictor

        fvVectorMatrix UEqn
        (
            fvm::ddt(U)
          + fvm::div(phi, U)
          - fvm::laplacian(nu, U)
        );

        if (piso.momentumPredictor())
        {
            solve(UEqn == -fvc::grad(p));
        }

        // --- PISO loop
        while (piso.correct())
        {
            volScalarField rAU(1.0/UEqn.A());
            volVectorField HbyA(constrainHbyA(rAU*UEqn.H(), U, p));
            surfaceScalarField phiHbyA
            (
                "phiHbyA",
                fvc::flux(HbyA)
              + fvc::interpolate(rAU)*fvc::ddtCorr(U, phi)
            );

            adjustPhi(phiHbyA, U, p);

            // Update the pressure BCs to ensure flux consistency
            constrainPressure(p, U, phiHbyA, rAU);

            // Non-orthogonal pressure corrector loop
            while (piso.correctNonOrthogonal())
            {
                // Pressure corrector

                fvScalarMatrix pEqn
                (
                    fvm::laplacian(rAU, p) == fvc::div(phiHbyA)
                );

                pEqn.setReference(pRefCell, pRefValue);

                pEqn.solve(mesh.solver(p.select(piso.finalInnerIter())));

                if (piso.finalNonOrthogonalIter())
                {
                    phi = phiHbyA - pEqn.flux();
                }
            }

            #include "continuityErrs.H"

            U = HbyA - rAU*fvc::grad(p);
            U.correctBoundaryConditions();
        }

        fvScalarMatrix TEqn
        (
            fvm::ddt(T)
            + fvm::div(phi, T)
            - fvm::laplacian(DT, T)
        );

        TEqn.solve();

        runTime.write();

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

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

    return 0;
}


// ************************************************************************* //
Code:
Info<< "Reading transportProperties\n" << endl;

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

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

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

Info<< "Reading field p\n" << endl;
volScalarField p
(
    IOobject
    (
        "p",
        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
);

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


#include "createPhi.H"


label pRefCell = 0;
scalar pRefValue = 0.0;
setRefCell(p, mesh.solutionDict().subDict("PISO"), pRefCell, pRefValue);
mesh.setFluxRequired(p.name());
And this is the temperature BC:
Code:
dimensions      [0 0 0 1 0 0 0];

internalField   uniform 0;

boundaryField
{
    wall
    {
        type            zeroGradient;
    }

    inletSmall
    {
        type            fixedValue;
        value            uniform 40;
    }

    inletBig
    {
        type            fixedValue;
        value            uniform 20;
    }

    outlet
    {
        type            zeroGradient;
    }
}
Does anyone know how to help me?
Thanks!
Attached Images
File Type: png T.png (12.9 KB, 11 views)
mgab is offline   Reply With Quote

Old   December 19, 2016, 07:12
Default
  #2
Senior Member
 
Agustín Villa
Join Date: Apr 2013
Location: Alcorcón
Posts: 314
Rep Power: 15
agustinvo is on a distinguished road
Hi,

at what time did you take that image? If you are working in unsteady, and you know your mass flows from the two inlets, you can have an idea about when you should compare. It seems it is still a bit soon to compare.
agustinvo is offline   Reply With Quote

Old   December 19, 2016, 07:28
Default
  #3
New Member
 
Join Date: May 2016
Posts: 25
Rep Power: 10
mgab is on a distinguished road
The image is taken at time=1 second.
The fact is that there is no variation from t=0
mgab is offline   Reply With Quote

Old   December 19, 2016, 07:30
Default
  #4
Senior Member
 
Agustín Villa
Join Date: Apr 2013
Location: Alcorcón
Posts: 314
Rep Power: 15
agustinvo is on a distinguished road
And the peaks of temperature you see close to the inlets? Are they related with the mesh? What are the dimensions and the velocity you have?
agustinvo is offline   Reply With Quote

Old   December 19, 2016, 07:38
Default
  #5
New Member
 
Join Date: May 2016
Posts: 25
Rep Power: 10
mgab is on a distinguished road
Those are the boundary conditions. If I change the visualization to T with the cube the domain if completely uniform.
Inlet velocities are 0.4 m/s for the big inlet and 1.2 m/s for the small one.
Attached it's the velocity field
Attached Images
File Type: jpg v.jpg (24.7 KB, 10 views)
mgab is offline   Reply With Quote

Reply

Tags
icofoam, icofoam problem, solver, temperature

Thread Tools Search this Thread
Search this Thread:

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


Similar Threads
Thread Thread Starter Forum Replies Last Post
How I can introduce my power heat (W) in chtMultiRegionFoam? aminem OpenFOAM Pre-Processing 32 August 29, 2019 03:23
adding temperature to simpleFoam waters OpenFOAM Programming & Development 80 May 30, 2017 22:30
chtMultiRegionSimpleFoam samiam1000 OpenFOAM Running, Solving & CFD 39 March 31, 2016 09:43
Adding porous zones in icoFoam solver josephn OpenFOAM Running, Solving & CFD 4 March 7, 2015 01:28
Adding temperature field to InterFoam yapalparvi OpenFOAM Running, Solving & CFD 8 October 14, 2009 21:18


All times are GMT -4. The time now is 17:57.