[awacs]

Dear foamers,

I'd like to create a piecewise viscosity model described as following:

(1) when T<Tg, viscosity is equal to a constant:

nu = nu0;

(2) when T>=Tg, viscosity varies with temperature, pressure and strain rate:

nu= (B_*exp(B1_*pd)*exp(Tb_/T))/(RHO_*(scalar(1)+pow((B_*exp(B1_*pd)*exp(Tb_/T))*strainRate()/t1_, scalar(1)-n_)));

where B, B1, Tb, t1, n are material constants.

How to define this piecewise function in calcNu() in user-defined viscosity model?

Best regards,
Jitao

[awacs]

The modefied NewCrossArrhenius.C:

#include "NewCrossArrhenius.H"
//#include "twoPhaseMixture.H"
#include "addToRunTimeSelectionTable.H"
#include "surfaceFields.H"

// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //

namespace Foam
{
namespace viscosityModels
{
defineTypeNameAndDebug(NewCrossArrhenius, 0);

addToRunTimeSelectionTable
(
viscosityModel,
NewCrossArrhenius,
dictionary
);
}
}


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

Foam::tmp<Foam::volScalarField>
Foam::viscosityModels::NewCrossArrhenius::calcNu() const
{

const volScalarField& T=U_.mesh().lookupObject<volScalarField>("T");
const volScalarField& pd=U_.mesh().lookupObject<volScalarField>("pd");


forAll(T, cell)
{
if ( T [cell]<Tg_)

return INF_;

else

return (B_*exp(B1_*pd)*exp(Tb_/T))/(RHO_*(scalar(1)+pow((B_*exp(B1_*pd)*exp(Tb_/T))*strainRate()/t1_, scalar(1)-n_)));
}


}


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

Foam::viscosityModels::NewCrossArrhenius::NewCross Arrhenius
(
const word& name,
const dictionary& viscosityProperties,
const volVectorField& U,
const surfaceScalarField& phi
)
:
viscosityModel(name, viscosityProperties, U, phi),
NewCrossArrheniusCoeffs_(viscosityProperties.subDi ct(typeName + "Coeffs")),
B_(NewCrossArrheniusCoeffs_.lookup("B")),
B1_(NewCrossArrheniusCoeffs_.lookup("B1")),
Tb_(NewCrossArrheniusCoeffs_.lookup("Tb")),
t1_(NewCrossArrheniusCoeffs_.lookup("t1")),
n_(NewCrossArrheniusCoeffs_.lookup("n")),
RHO_(NewCrossArrheniusCoeffs_.lookup("RHO")),
INF_(NewCrossArrheniusCoeffs_.lookup("INF")),
Tg_(NewCrossArrheniusCoeffs_.lookup("Tg")),


nu_
(
IOobject
(
name,
U_.time().timeName(),
U_.db(),
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
calcNu()
)
{}


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

bool Foam::viscosityModels::NewCrossArrhenius::read
(
const dictionary& viscosityProperties
)
{
viscosityModel::read(viscosityProperties);

NewCrossArrheniusCoeffs_ = viscosityProperties.subDict(typeName + "Coeffs");

NewCrossArrheniusCoeffs_.lookup("B") >> B_;
NewCrossArrheniusCoeffs_.lookup("B1") >> B1_;
NewCrossArrheniusCoeffs_.lookup("Tb") >> Tb_;
NewCrossArrheniusCoeffs_.lookup("t1") >> t1_;
NewCrossArrheniusCoeffs_.lookup("n") >> n_;
NewCrossArrheniusCoeffs_.lookup("RHO") >> RHO_;
NewCrossArrheniusCoeffs_.lookup("INF") >> INF_;
NewCrossArrheniusCoeffs_.lookup("Tg") >> Tg_;

return true;
}


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

The compilation of this visocosity model ended up with following errors :

viscosityModels/NewCrossArrhenius/NewCrossArrhenius.C:81: error: no match for ‘operator<’ in ‘((const Foam::volScalarField*)T)->Foam::GeometricField<double,Foam::fvPatchField, Foam::volMesh>::<anonymous>.Foam:imensionedField <double, Foam::volMesh>::<anonymous>.Foam::Field<double>::< anonymous>.Foam::List<double>::<anonymous>.Foam::U List<T>:perator[] [with T = double](cell) < ((const Foam::viscosityModels::NewCrossArrhenius*)this)->Foam::viscosityModels::NewCrossArrhenius::Tg_’
viscosityModels/NewCrossArrhenius/NewCrossArrhenius.C:83: error: conversion from ‘const Foam::dimensionedScalar’ to non-scalar type ‘Foam::tmp<Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh> >’ requested

[dgadensg]

Hi..

I have just recently had similar problems in implementing a frozen layer. I did not succeed as it showed som strange results, however i did get a working compilation.

When comparing a cell value (T[cell]) with a constant (Tg_) i found that it was neccessary to do it with the suffix .value() on the constant, as shown below:

...

if ( T [cell]<Tg_.value())
...

Hope that works!

[awacs]

Originally Posted by dgadensg

Hi..

I have just recently had similar problems in implementing a frozen layer. I did not succeed as it showed som strange results, however i did get a working compilation.

When comparing a cell value (T[cell]) with a constant (Tg_) i found that it was neccessary to do it with the suffix .value() on the constant, as shown below:

...

if ( T [cell]<Tg_.value())
...

Hope that works!

Hi Dan,

Thank you very much. I have compiled this Piecewise viscosity model. When running a case using this model, the calculation always stopped due divergence. And the claculated temperature field is strange.

Have you successed in implementing the frozen layer? Please give me some suggections. Thanks in advance.

Kind regards,
Jitao