[mnaufalazwar]

Hi guys! Im solving incompressible flow in a 2D axisymmetrical pipe with heat transfer in steady state, but i have problem with my work.

First I modify simpleFoam solver and I name it simpleThermFoam, this is how I modified the solver :

createFields.H :

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,
dimensionedVector("U",dimensionSet(0,1,-1,0,0,0,0),vector::zero)
);

#include "createPhi.H"


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

singlePhaseTransportModel laminarTransport(U, phi);

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

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

Info<< "Reading transportProperties\n" << endl;
IOdictionary transportProperties
(
IOobject
(
"transportProperties",
runTime.constant(),
mesh,
IOobject::MUST_READ_IF_MODIFIED,
IOobject::NO_WRITE
)
);

Info<< "Reading density rho\n" << endl;
dimensionedScalar rho
(
transportProperties.lookup("rho")
);

Info<< "Reading heat capacity Cp\n" << endl;
dimensionedScalar Cp
(
transportProperties.lookup("Cp")
);

Info<< "Reading konductivity k\n" << endl;
dimensionedScalar k
(
transportProperties.lookup("k")
);

simpleThermFoam.C :

#include "fvCFD.H"
#include "singlePhaseTransportModel.H"
#include "turbulentTransportModel.H"
#include "simpleControl.H"
#include "fvIOoptionList.H"

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

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

simpleControl simple(mesh);

#include "createFields.H"
#include "createMRF.H"
#include "createFvOptions.H"
#include "initContinuityErrs.H"

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

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"
}

laminarTransport.correct();
turbulence->correct();

// --- Temperature field
solve
(
fvm::ddt(rho * Cp, T)
+ rho * Cp * fvm::div(phi, T)
- fvm::laplacian(k, T)
);

runTime.write();

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

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

return 0;
}

and then I modify totalTemperature boundary condition and I name it convectiveHeatFlux, this is how I modify the boundary condition :

convectiveHeatFluxFvPatchScalarField.H :

#ifndef convectiveHeatFluxFvPatchScalarField_H
#define convectiveHeatFluxFvPatchScalarField_H

//#include "fixedValueFvPatchFields.H"
#include "fvPatchFields.H"
#include "fixedGradientFvPatchFields.H"

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

namespace Foam
{

/*---------------------------------------------------------------------------*\
Class convectiveHeatFluxFvPatchScalarField Declaration
\*---------------------------------------------------------------------------*/

class convectiveHeatFluxFvPatchScalarField
:
public fixedGradientFvPatchScalarField
{
// Private data

//- Name of the temperature field
word TName_;
public:

//- Runtime type information
TypeName("convectiveHeatFlux");


// Constructors

//- Construct from patch and internal field
convectiveHeatFluxFvPatchScalarField
(
const fvPatch&,
const DimensionedField<scalar, volMesh>&
);

//- Construct from patch, internal field and dictionary
convectiveHeatFluxFvPatchScalarField
(
const fvPatch&,
const DimensionedField<scalar, volMesh>&,
const dictionary&
);

//- Construct by mapping given convectiveHeatFluxFvPatchScalarField
// onto a new patch
convectiveHeatFluxFvPatchScalarField
(
const convectiveHeatFluxFvPatchScalarField&,
const fvPatch&,
const DimensionedField<scalar, volMesh>&,
const fvPatchFieldMapper&
);

//- Construct as copy
convectiveHeatFluxFvPatchScalarField
(
const convectiveHeatFluxFvPatchScalarField&
);

//- Construct and return a clone
virtual tmp<fvPatchScalarField> clone() const
{
return tmp<fvPatchScalarField>
(
new convectiveHeatFluxFvPatchScalarField(*this)
);
}

//- Construct as copy setting internal field reference
convectiveHeatFluxFvPatchScalarField
(
const convectiveHeatFluxFvPatchScalarField&,
const DimensionedField<scalar, volMesh>&
);

//- Construct and return a clone setting internal field reference
virtual tmp<fvPatchScalarField> clone
(
const DimensionedField<scalar, volMesh>& iF
) const
{
return tmp<fvPatchScalarField>
(
new convectiveHeatFluxFvPatchScalarField(*this, iF)
);
}

// Evaluation functions

//- Update the coefficients associated with the patch field
virtual void updateCoeffs();


//- Write
virtual void write(Ostream&) const;
};


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

} // End namespace Foam

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

#endif

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

convectiveHeatFluxFvPatchScalarField.C :

#include "convectiveHeatFluxFvPatchScalarField.H"
#include "addToRunTimeSelectionTable.H"
#include "fvPatchFieldMapper.H"
#include "volFields.H"
#include "surfaceFields.H"

// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //

Foam::convectiveHeatFluxFvPatchScalarField::convec tiveHeatFluxFvPatchScalarField
(
const fvPatch& p,
const DimensionedField<scalar, volMesh>& iF
)
:
fixedGradientFvPatchScalarField(p, iF),
TName_("T")
{}


Foam::convectiveHeatFluxFvPatchScalarField::convec tiveHeatFluxFvPatchScalarField
(
const convectiveHeatFluxFvPatchScalarField& ptf,
const fvPatch& p,
const DimensionedField<scalar, volMesh>& iF,
const fvPatchFieldMapper& mapper
)
:
fixedGradientFvPatchScalarField(ptf, p, iF, mapper),
TName_(ptf.TName_)
{}


Foam::convectiveHeatFluxFvPatchScalarField::convec tiveHeatFluxFvPatchScalarField
(
const fvPatch& p,
const DimensionedField<scalar, volMesh>& iF,
const dictionary& dict
)
:
fixedGradientFvPatchScalarField(p, iF),
TName_(dict.lookupOrDefault<word>("T", "T"))
{
fvPatchField<scalar>:perator=(patchInternalField ());
gradient() = 0.0;
}


Foam::convectiveHeatFluxFvPatchScalarField::convec tiveHeatFluxFvPatchScalarField
(
const convectiveHeatFluxFvPatchScalarField& tppsf
)
:
fixedGradientFvPatchScalarField(tppsf),
TName_(tppsf.TName_)
{}


Foam::convectiveHeatFluxFvPatchScalarField::convec tiveHeatFluxFvPatchScalarField
(
const convectiveHeatFluxFvPatchScalarField& tppsf,
const DimensionedField<scalar, volMesh>& iF
)
:
fixedGradientFvPatchScalarField(tppsf, iF),
TName_(tppsf.TName_)
{}


// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //

void Foam::convectiveHeatFluxFvPatchScalarField::update Coeffs()
{
if (updated())
{
return;
}

const fvPatchField<scalar>& T =
patch().lookupPatchField<volScalarField, scalar>(TName_);

const dictionary& transportPropertiesIsolant = db().lookupObject<IOdictionary> ("transportProperties");

dimensionedScalar Tsur (transportPropertiesIsolant.lookup("Tsur"));
dimensionedScalar h (transportPropertiesIsolant.lookup("h"));
dimensionedScalar k (transportPropertiesIsolant.lookup("k"));

gradient() = h.value() / k.value() * (Tsur.value() - T);

fixedGradientFvPatchScalarField::updateCoeffs();
}


void Foam::convectiveHeatFluxFvPatchScalarField::write( Ostream& os) const
{
fvPatchScalarField::write(os);
writeEntryIfDifferent<word>(os, "T", "T", TName_);
writeEntry("value", os);
}


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

namespace Foam
{
makePatchTypeField
(
fvPatchScalarField,
convectiveHeatFluxFvPatchScalarField
);
}

and this is how I make my case :

0/p :

FoamFile
{
version 2.0;
format ascii;
class volScalarField;
object p;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions [0 2 -2 0 0 0 0];

internalField uniform 0;

boundaryField
{
axis
{
type symmetry;
}

inlet
{
type zeroGradient;
}

pipe
{
type zeroGradient;
}

outlet
{
type fixedValue;
value uniform 0;
}

front
{
type wedge;
}

back
{
type wedge;
}
}

0/T :

FoamFile
{
version 2.0;
format ascii;
class volScalarField;
object T;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions [0 0 0 1 0 0 0];

internalField uniform 20;

boundaryField
{
inlet
{
type fixedValue;
value uniform 20;
}

outlet
{
type zeroGradient;
}

pipe
{
type convectiveHeatFlux;
}

axis
{
type symmetry;
}

front
{
type wedge;
}

back
{
type wedge;
}
}

0/U :

FoamFile
{
version 2.0;
format ascii;
class volVectorField;
object U;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions [0 1 -1 0 0 0 0];

internalField uniform (0 0 0);

boundaryField
{
axis
{
type symmetry;
}

inlet
{
type fixedValue;
value uniform (0.1 0 0);

}

pipe
{
type fixedValue;
value uniform (0 0 0);
}

outlet
{
type zeroGradient;
}

front
{
type wedge;
}

back
{
type wedge;
}
}

blockMeshDict :

convertToMeters 0.01;

vertices
(
(0 0 0)
(0 0.9990482216 -0.0436193874)
(0 0.9990482216 0.0436193874)
(10 0 0)
(10 0.9990482216 -0.0436193874)
(10 0.9990482216 0.0436193874)
);

blocks
(
hex (0 1 2 0 3 4 5 3) (30 1 300) simpleGrading (0.1 1 1)
);

edges
(
);

boundary
(
axis
{
type symmetry;
faces
(
(0 3 3 0)
);
}

inlet
{
type patch;
faces
(
(0 0 2 1)
);
}
pipe
{
type wall;
faces
(
(2 5 4 1)
);
}
outlet
{
type patch;
faces
(
(3 4 5 3)
);
}

front
{
type wedge;
faces
(
(0 3 5 2)
);
}

back
{
type wedge;
faces
(
(0 1 4 3)
);
}
);

mergePatchPairs
(
);

fvSchemes :

FoamFile
{
version 2.0;
format ascii;
class dictionary;
location "system";
object fvSchemes;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

ddtSchemes
{
default steadyState;
}

gradSchemes
{
default Gauss linear;
}

divSchemes
{
default none;
div(phi,U) bounded Gauss linearUpwind grad(U);
div(phi,k) bounded Gauss limitedLinear 1;
div(phi,epsilon) bounded Gauss limitedLinear 1;
div(phi,omega) bounded Gauss limitedLinear 1;
div(phi,v2) bounded Gauss limitedLinear 1;
div((nuEff*dev2(T(grad(U))))) Gauss linear;
div(nonlinearStress) Gauss linear;
}

laplacianSchemes
{
default Gauss linear corrected;
}

interpolationSchemes
{
default linear;
}

snGradSchemes
{
default corrected;
}

wallDist
{
method meshWave;
}

fvSolution :

FoamFile
{
version 2.0;
format ascii;
class dictionary;
location "system";
object fvSolution;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

solvers
{
p
{
solver PCG;
preconditioner DIC;
tolerance 1e-06;
relTol 0;
}

U
{
solver PBiCG;
preconditioner DILU;
tolerance 1e-05;
relTol 0;
}

T
{
solver PBiCG;
preconditioner DILU;
tolerance 1e-06;
relTol 0;
}
}

SIMPLE
{
nNonOrthogonalCorrectors 0;
//pRefCell 0;
//pRefValue 0;
consistent yes;

residualControl
{
p 1e-6;
U 1e-2;
T 1e-6;
}
}

relaxationFactors
{
fields
{
p 0.3;
T 0.7;
}
equations
{
U 0.7; // 0.9 is more stable but 0.95 more convergent
}
}

controlDict :

FoamFile
{
version 2.0;
format ascii;
class dictionary;
location "system";
object controlDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

application simpleThermFoam;

startFrom latestTime;

startTime 0;

stopAt endTime;

endTime 1000;

deltaT 1;

writeControl timeStep;

writeInterval 50;

purgeWrite 0;

writeFormat ascii;

writePrecision 6;

writeCompression off;

timeFormat general;

timePrecision 6;

runTimeModifiable true;

libs ( "lconvectiveHeatFlux.so" );

but when I run the program the terminal said this :

/*---------------------------------------------------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: 3.0.1 |
| \\ / A nd | Web: www.OpenFOAM.org |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
Build : 3.0.1-d8a290b55d28
Exec : simpleThermFoam
Date : Aug 20 2018
Time : 15:19:12
Host : "DESKTOP-GDC6B34"
PID : 1268
Case : /home/naufal/OpenFOAM/naufal-3.0.1/run/Ex8
nProcs : 1
sigFpe : Enabling floating point exception trapping (FOAM_SIGFPE).
fileModificationChecking : Monitoring run-time modified files using timeStampMaster
allowSystemOperations : Allowing user-supplied system call operations

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Create time

Create mesh for time = 0


SIMPLE: convergence criteria
field p tolerance 1e-06
field U tolerance 0.01
field T tolerance 1e-06

Reading field p

Reading field U

--> FOAM Warning :
From function GeometricField<Type, PatchField, GeoMesh>::readIfPresent()
in file /home/naufal/OpenFOAM/OpenFOAM-3.0.1/src/OpenFOAM/lnInclude/GeometricField.C at line 106
read option IOobject::MUST_READ or MUST_READ_IF_MODIFIED suggests that a read constructor for field U would be more appropriate.
Reading/calculating face flux field phi

Selecting incompressible transport model Newtonian
Selecting turbulence model type laminar
Reading field T

Reading transportProperties

Reading density rho

Reading heat capacity Cp

Reading konductivity k

No MRF models present

No finite volume options present


Starting time loop

Time = 1



--> FOAM FATAL ERROR:

gradientInternalCoeffs cannot be called for a calculatedFvPatchField
on patch inlet of field U in file "/home/naufal/OpenFOAM/naufal-3.0.1/run/Ex8/0/U"
You are probably trying to solve for a field with a default boundary condition.

From function calculatedFvPatchField<Type>::gradientInternalCoef fs() const
in file fields/fvPatchFields/basic/calculated/calculatedFvPatchField.C at line 199.

FOAM exiting

Does anyone know how to solve this error? thanks a lot.