Merge branch 'master' of /home/noisy2/OpenFOAM/OpenFOAM-dev/

This commit is contained in:
andy
2008-05-09 13:04:42 +01:00
6 changed files with 253 additions and 312 deletions

View File

@ -75,6 +75,9 @@ int main(int argc, char *argv[])
# include "correctPhi.H" # include "correctPhi.H"
} }
// Keep the absolute fluxes for use in ddtPhiCorr
surfaceScalarField phiAbs("phiAbs", phi);
// Make the fluxes relative to the mesh motion // Make the fluxes relative to the mesh motion
fvc::makeRelative(phi, U); fvc::makeRelative(phi, U);
@ -95,7 +98,11 @@ int main(int argc, char *argv[])
U = rAU*UEqn.H(); U = rAU*UEqn.H();
phi = (fvc::interpolate(U) & mesh.Sf()); phi = (fvc::interpolate(U) & mesh.Sf());
//+ fvc::ddtPhiCorr(rAU, U, phi);
if (ddtPhiCorr)
{
phi += fvc::ddtPhiCorr(rAU, U, phiAbs);
}
adjustPhi(phi, U, p); adjustPhi(phi, U, p);

View File

@ -74,7 +74,7 @@ int main(int argc, char *argv[])
Info<< "Time = " << runTime.timeName() << nl << endl; Info<< "Time = " << runTime.timeName() << nl << endl;
// Make the fluxes absolute // Make the fluxes absolute
if (mesh.changing() && correctPhi) if (mesh.changing())
{ {
fvc::makeAbsolute(phi, U); fvc::makeAbsolute(phi, U);
} }
@ -99,8 +99,11 @@ int main(int argc, char *argv[])
#include "correctPhi.H" #include "correctPhi.H"
} }
// Keep the absolute fluxes for use in ddtPhiCorr
surfaceScalarField phiAbs("phiAbs", phi);
// Make the fluxes relative to the mesh motion // Make the fluxes relative to the mesh motion
if (mesh.changing() && correctPhi) if (mesh.changing())
{ {
fvc::makeRelative(phi, U); fvc::makeRelative(phi, U);
} }

View File

@ -5,12 +5,12 @@
volVectorField HU = UEqn.H(); volVectorField HU = UEqn.H();
U = rAU*HU; U = rAU*HU;
surfaceScalarField phiU surfaceScalarField phiU("phiU", (fvc::interpolate(U) & mesh.Sf()));
(
"phiU", if (ddtPhiCorr)
(fvc::interpolate(U) & mesh.Sf()) {
//+ fvc::ddtPhiCorr(rAU, rho, U, phiAbs) phiU += fvc::ddtPhiCorr(rAU, rho, U, phiAbs);
); }
phi = phiU + phi = phiU +
( (

View File

@ -25,3 +25,9 @@
{ {
nOuterCorr = readInt(piso.lookup("nOuterCorrectors")); nOuterCorr = readInt(piso.lookup("nOuterCorrectors"));
} }
bool ddtPhiCorr = false;
if (piso.found("ddtPhiCorr"))
{
ddtPhiCorr = Switch(piso.lookup("ddtPhiCorr"));
}

View File

@ -1,4 +1,4 @@
/*---------------------------------------------------------------------------*\ /*---------------------------------------------------------------------------* \
========= | ========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | \\ / O peration |
@ -375,50 +375,30 @@ EulerDdtScheme<Type>::fvcDdtPhiCorr
( (
const volScalarField& rA, const volScalarField& rA,
const GeometricField<Type, fvPatchField, volMesh>& U, const GeometricField<Type, fvPatchField, volMesh>& U,
const fluxFieldType& phi const fluxFieldType& phiAbs
) )
{ {
dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT(); dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
IOobject ddtIOobject IOobject ddtIOobject
( (
"ddtPhiCorr(" + rA.name() + ',' + U.name() + ',' + phi.name() + ')', "ddtPhiCorr(" + rA.name() + ',' + U.name() + ',' + phiAbs.name() + ')',
mesh().time().timeName(), mesh().time().timeName(),
mesh() mesh()
); );
if (mesh().moving()) tmp<fluxFieldType> phiCorr =
{ phiAbs.oldTime() - (fvc::interpolate(U.oldTime()) & mesh().Sf());
return tmp<fluxFieldType>
(
new fluxFieldType
(
ddtIOobject,
mesh(),
dimensioned<typename flux<Type>::type>
(
"0",
rDeltaT.dimensions()*rA.dimensions()*phi.dimensions(),
pTraits<typename flux<Type>::type>::zero
)
)
);
}
else
{
tmp<fluxFieldType> phiCorr =
phi.oldTime() - (fvc::interpolate(U.oldTime()) & mesh().Sf());
return tmp<fluxFieldType> return tmp<fluxFieldType>
(
new fluxFieldType
( (
new fluxFieldType ddtIOobject,
( fvcDdtPhiCoeff(U.oldTime(), phiAbs.oldTime(), phiCorr())
ddtIOobject, *fvc::interpolate(rDeltaT*rA)*phiCorr
fvcDdtPhiCoeff(U.oldTime(), phi.oldTime(), phiCorr()) )
*fvc::interpolate(rDeltaT*rA)*phiCorr );
)
);
}
} }
@ -429,7 +409,7 @@ EulerDdtScheme<Type>::fvcDdtPhiCorr
const volScalarField& rA, const volScalarField& rA,
const volScalarField& rho, const volScalarField& rho,
const GeometricField<Type, fvPatchField, volMesh>& U, const GeometricField<Type, fvPatchField, volMesh>& U,
const fluxFieldType& phi const fluxFieldType& phiAbs
) )
{ {
dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT(); dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
@ -437,111 +417,94 @@ EulerDdtScheme<Type>::fvcDdtPhiCorr
IOobject ddtIOobject IOobject ddtIOobject
( (
"ddtPhiCorr(" "ddtPhiCorr("
+ rA.name() + ',' + rho.name() + ',' + U.name() + ',' + phi.name() + ')', + rA.name() + ','
+ rho.name() + ','
+ U.name() + ','
+ phiAbs.name() + ')',
mesh().time().timeName(), mesh().time().timeName(),
mesh() mesh()
); );
if (mesh().moving()) if
(
U.dimensions() == dimVelocity
&& phiAbs.dimensions() == dimVelocity*dimArea
)
{ {
return tmp<fluxFieldType> return tmp<fluxFieldType>
( (
new fluxFieldType new fluxFieldType
( (
ddtIOobject, ddtIOobject,
mesh(), rDeltaT
dimensioned<typename flux<Type>::type> *fvcDdtPhiCoeff(U.oldTime(), phiAbs.oldTime())
*(
fvc::interpolate(rA*rho.oldTime())*phiAbs.oldTime()
- (fvc::interpolate(rA*rho.oldTime()*U.oldTime())
& mesh().Sf())
)
)
);
}
else if
(
U.dimensions() == dimVelocity
&& phiAbs.dimensions() == rho.dimensions()*dimVelocity*dimArea
)
{
return tmp<fluxFieldType>
(
new fluxFieldType
(
ddtIOobject,
rDeltaT
*fvcDdtPhiCoeff
( (
"0", U.oldTime(),
rDeltaT.dimensions()*rA.dimensions() phiAbs.oldTime()/fvc::interpolate(rho.oldTime())
*rho.dimensions()*phi.dimensions(), )
pTraits<typename flux<Type>::type>::zero *(
fvc::interpolate(rA*rho.oldTime())
*phiAbs.oldTime()/fvc::interpolate(rho.oldTime())
- (
fvc::interpolate
(
rA*rho.oldTime()*U.oldTime()
) & mesh().Sf()
)
)
)
);
}
else if
(
U.dimensions() == rho.dimensions()*dimVelocity
&& phiAbs.dimensions() == rho.dimensions()*dimVelocity*dimArea
)
{
return tmp<fluxFieldType>
(
new fluxFieldType
(
ddtIOobject,
rDeltaT
*fvcDdtPhiCoeff(rho.oldTime(), U.oldTime(), phiAbs.oldTime())
*(
fvc::interpolate(rA)*phiAbs.oldTime()
- (fvc::interpolate(rA*U.oldTime()) & mesh().Sf())
) )
) )
); );
} }
else else
{ {
if FatalErrorIn
( (
U.dimensions() == dimVelocity "EulerDdtScheme<Type>::fvcDdtPhiCorr"
&& phi.dimensions() == dimVelocity*dimArea ) << "dimensions of phiAbs are not correct"
) << abort(FatalError);
{
return tmp<fluxFieldType>
(
new fluxFieldType
(
ddtIOobject,
rDeltaT*fvcDdtPhiCoeff(U.oldTime(), phi.oldTime())
*(
fvc::interpolate(rA*rho.oldTime())*phi.oldTime()
- (fvc::interpolate(rA*rho.oldTime()*U.oldTime())
& mesh().Sf())
)
)
);
}
else if
(
U.dimensions() == dimVelocity
&& phi.dimensions() == rho.dimensions()*dimVelocity*dimArea
)
{
return tmp<fluxFieldType>
(
new fluxFieldType
(
ddtIOobject,
rDeltaT
*fvcDdtPhiCoeff
(
U.oldTime(),
phi.oldTime()/fvc::interpolate(rho.oldTime())
)
*(
fvc::interpolate(rA*rho.oldTime())
*phi.oldTime()/fvc::interpolate(rho.oldTime())
- (
fvc::interpolate
(
rA*rho.oldTime()*U.oldTime()
) & mesh().Sf()
)
)
)
);
}
else if
(
U.dimensions() == rho.dimensions()*dimVelocity
&& phi.dimensions() == rho.dimensions()*dimVelocity*dimArea
)
{
return tmp<fluxFieldType>
(
new fluxFieldType
(
ddtIOobject,
rDeltaT
*fvcDdtPhiCoeff(rho.oldTime(), U.oldTime(), phi.oldTime())
*(
fvc::interpolate(rA)*phi.oldTime()
- (fvc::interpolate(rA*U.oldTime()) & mesh().Sf())
)
)
);
}
else
{
FatalErrorIn
(
"EulerDdtScheme<Type>::fvcDdtPhiCorr"
) << "dimensions of phi are not correct"
<< abort(FatalError);
return fluxFieldType::null(); return fluxFieldType::null();
}
} }
} }

View File

@ -531,43 +531,160 @@ backwardDdtScheme<Type>::fvcDdtPhiCorr
mesh() mesh()
); );
if (mesh().moving()) scalar deltaT = deltaT_();
scalar deltaT0 = deltaT0_(U);
scalar coefft = 1 + deltaT/(deltaT + deltaT0);
scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
scalar coefft0 = coefft + coefft00;
return tmp<fluxFieldType>
(
new fluxFieldType
(
ddtIOobject,
rDeltaT*fvcDdtPhiCoeff(U.oldTime(), phi.oldTime())
*(
fvc::interpolate(rA)
*(
coefft0*phi.oldTime()
- coefft00*phi.oldTime().oldTime()
)
- (
fvc::interpolate
(
rA*
(
coefft0*U.oldTime()
- coefft00*U.oldTime().oldTime()
)
) & mesh().Sf()
)
)
)
);
}
template<class Type>
tmp<typename backwardDdtScheme<Type>::fluxFieldType>
backwardDdtScheme<Type>::fvcDdtPhiCorr
(
const volScalarField& rA,
const volScalarField& rho,
const GeometricField<Type, fvPatchField, volMesh>& U,
const fluxFieldType& phiAbs
)
{
dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
IOobject ddtIOobject
(
"ddtPhiCorr("
+ rA.name() + ','
+ rho.name() + ','
+ U.name() + ','
+ phiAbs.name() + ')',
mesh().time().timeName(),
mesh()
);
scalar deltaT = deltaT_();
scalar deltaT0 = deltaT0_(U);
scalar coefft = 1 + deltaT/(deltaT + deltaT0);
scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
scalar coefft0 = coefft + coefft00;
if
(
U.dimensions() == dimVelocity
&& phiAbs.dimensions() == dimVelocity*dimArea
)
{ {
return tmp<fluxFieldType> return tmp<fluxFieldType>
( (
new fluxFieldType new fluxFieldType
( (
ddtIOobject, ddtIOobject,
mesh(), rDeltaT*fvcDdtPhiCoeff(U.oldTime(), phiAbs.oldTime())
dimensioned<typename flux<Type>::type> *(
( coefft0*fvc::interpolate(rA*rho.oldTime())
"0", *phiAbs.oldTime()
rDeltaT.dimensions()*rA.dimensions()*phi.dimensions(), - coefft00*fvc::interpolate(rA*rho.oldTime().oldTime())
pTraits<typename flux<Type>::type>::zero *phiAbs.oldTime().oldTime()
- (
fvc::interpolate
(
rA*
(
coefft0*rho.oldTime()*U.oldTime()
- coefft00*rho.oldTime().oldTime()
*U.oldTime().oldTime()
)
) & mesh().Sf()
)
) )
) )
); );
} }
else else if
(
U.dimensions() == dimVelocity
&& phiAbs.dimensions() == rho.dimensions()*dimVelocity*dimArea
)
{ {
scalar deltaT = deltaT_();
scalar deltaT0 = deltaT0_(U);
scalar coefft = 1 + deltaT/(deltaT + deltaT0);
scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
scalar coefft0 = coefft + coefft00;
return tmp<fluxFieldType> return tmp<fluxFieldType>
( (
new fluxFieldType new fluxFieldType
( (
ddtIOobject, ddtIOobject,
rDeltaT*fvcDdtPhiCoeff(U.oldTime(), phi.oldTime()) rDeltaT
*fvcDdtPhiCoeff
(
U.oldTime(),
phiAbs.oldTime()/fvc::interpolate(rho.oldTime())
)
*(
fvc::interpolate(rA*rho.oldTime())
*(
coefft0*phiAbs.oldTime()
/fvc::interpolate(rho.oldTime())
- coefft00*phiAbs.oldTime().oldTime()
/fvc::interpolate(rho.oldTime().oldTime())
)
- (
fvc::interpolate
(
rA*rho.oldTime()*
(
coefft0*U.oldTime()
- coefft00*U.oldTime().oldTime()
)
) & mesh().Sf()
)
)
)
);
}
else if
(
U.dimensions() == rho.dimensions()*dimVelocity
&& phiAbs.dimensions() == rho.dimensions()*dimVelocity*dimArea
)
{
return tmp<fluxFieldType>
(
new fluxFieldType
(
ddtIOobject,
rDeltaT
*fvcDdtPhiCoeff(rho.oldTime(), U.oldTime(), phiAbs.oldTime())
*( *(
fvc::interpolate(rA) fvc::interpolate(rA)
*( *(
coefft0*phi.oldTime() coefft0*phiAbs.oldTime()
- coefft00*phi.oldTime().oldTime() - coefft00*phiAbs.oldTime().oldTime()
) )
- ( - (
fvc::interpolate fvc::interpolate
@ -583,170 +700,15 @@ backwardDdtScheme<Type>::fvcDdtPhiCorr
) )
); );
} }
}
template<class Type>
tmp<typename backwardDdtScheme<Type>::fluxFieldType>
backwardDdtScheme<Type>::fvcDdtPhiCorr
(
const volScalarField& rA,
const volScalarField& rho,
const GeometricField<Type, fvPatchField, volMesh>& U,
const fluxFieldType& phi
)
{
dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
IOobject ddtIOobject
(
"ddtPhiCorr("
+ rA.name() + ',' + rho.name() + ',' + U.name() + ',' + phi.name() + ')',
mesh().time().timeName(),
mesh()
);
if (mesh().moving())
{
return tmp<fluxFieldType>
(
new fluxFieldType
(
ddtIOobject,
mesh(),
dimensioned<typename flux<Type>::type>
(
"0",
rDeltaT.dimensions()*rA.dimensions()
*rho.dimensions()*phi.dimensions(),
pTraits<typename flux<Type>::type>::zero
)
)
);
}
else else
{ {
scalar deltaT = deltaT_(); FatalErrorIn
scalar deltaT0 = deltaT0_(U);
scalar coefft = 1 + deltaT/(deltaT + deltaT0);
scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
scalar coefft0 = coefft + coefft00;
if
( (
U.dimensions() == dimVelocity "backwardDdtScheme<Type>::fvcDdtPhiCorr"
&& phi.dimensions() == dimVelocity*dimArea ) << "dimensions of phiAbs are not correct"
) << abort(FatalError);
{
return tmp<fluxFieldType>
(
new fluxFieldType
(
ddtIOobject,
rDeltaT*fvcDdtPhiCoeff(U.oldTime(), phi.oldTime())
*(
coefft0*fvc::interpolate(rA*rho.oldTime())
*phi.oldTime()
- coefft00*fvc::interpolate(rA*rho.oldTime().oldTime())
*phi.oldTime().oldTime()
- (
fvc::interpolate
(
rA*
(
coefft0*rho.oldTime()*U.oldTime()
- coefft00*rho.oldTime().oldTime()
*U.oldTime().oldTime()
)
) & mesh().Sf()
)
)
)
);
}
else if
(
U.dimensions() == dimVelocity
&& phi.dimensions() == rho.dimensions()*dimVelocity*dimArea
)
{
return tmp<fluxFieldType>
(
new fluxFieldType
(
ddtIOobject,
rDeltaT
*fvcDdtPhiCoeff
(
U.oldTime(),
phi.oldTime()/fvc::interpolate(rho.oldTime())
)
*(
fvc::interpolate(rA*rho.oldTime())
*(
coefft0*phi.oldTime()
/fvc::interpolate(rho.oldTime())
- coefft00*phi.oldTime().oldTime()
/fvc::interpolate(rho.oldTime().oldTime())
)
- (
fvc::interpolate
(
rA*rho.oldTime()*
(
coefft0*U.oldTime()
- coefft00*U.oldTime().oldTime()
)
) & mesh().Sf()
)
)
)
);
}
else if
(
U.dimensions() == rho.dimensions()*dimVelocity
&& phi.dimensions() == rho.dimensions()*dimVelocity*dimArea
)
{
return tmp<fluxFieldType>
(
new fluxFieldType
(
ddtIOobject,
rDeltaT
*fvcDdtPhiCoeff(rho.oldTime(), U.oldTime(), phi.oldTime())
*(
fvc::interpolate(rA)
*(
coefft0*phi.oldTime()
- coefft00*phi.oldTime().oldTime()
)
- (
fvc::interpolate
(
rA*
(
coefft0*U.oldTime()
- coefft00*U.oldTime().oldTime()
)
) & mesh().Sf()
)
)
)
);
}
else
{
FatalErrorIn
(
"backwardDdtScheme<Type>::fvcDdtPhiCorr"
) << "dimensions of phi are not correct"
<< abort(FatalError);
return fluxFieldType::null(); return fluxFieldType::null();
}
} }
} }