diff --git a/applications/solvers/compressible/rhoSimpleFoam/pEqn.H b/applications/solvers/compressible/rhoSimpleFoam/pEqn.H index cefee48969..6c64056be2 100644 --- a/applications/solvers/compressible/rhoSimpleFoam/pEqn.H +++ b/applications/solvers/compressible/rhoSimpleFoam/pEqn.H @@ -54,7 +54,8 @@ else { fvScalarMatrix pEqn ( - fvm::laplacian(rho*rAU, p) == fvc::div(phiHbyA) + fvc::div(phiHbyA) + - fvm::laplacian(rho*rAU, p) ); pEqn.setReference(pRefCell, pRefValue); @@ -63,7 +64,7 @@ else if (simple.finalNonOrthogonalIter()) { - phi = phiHbyA - pEqn.flux(); + phi = phiHbyA + pEqn.flux(); } } } @@ -88,5 +89,10 @@ if (closedVolume) rho = thermo.rho(); rho = max(rho, rhoMin); rho = min(rho, rhoMax); -rho.relax(); + +if (!simple.transonic()) +{ + rho.relax(); +} + Info<< "rho max/min : " << max(rho).value() << " " << min(rho).value() << endl; diff --git a/applications/solvers/compressible/rhoSimpleFoam/rhoSimplecFoam/pEqn.H b/applications/solvers/compressible/rhoSimpleFoam/rhoSimplecFoam/pEqn.H index 379e720c04..44babc3a6b 100644 --- a/applications/solvers/compressible/rhoSimpleFoam/rhoSimplecFoam/pEqn.H +++ b/applications/solvers/compressible/rhoSimpleFoam/rhoSimplecFoam/pEqn.H @@ -3,8 +3,6 @@ rho = max(rho, rhoMin); rho = min(rho, rhoMax); rho.relax(); -volScalarField p0(p); - volScalarField AU(UEqn().A()); volScalarField AtU(AU - UEqn().H1()); @@ -22,7 +20,7 @@ if (simple.transonic()) surfaceScalarField phid ( "phid", - fvc::interpolate(psi*HbyA) & mesh.Sf() + fvc::interpolate(psi)*(fvc::interpolate(HbyA) & mesh.Sf()) ); surfaceScalarField phic @@ -31,6 +29,8 @@ if (simple.transonic()) fvc::interpolate(rho/AtU - rho/AU)*fvc::snGrad(p)*mesh.magSf() ); + HbyA -= fvc::grad(p)*(1.0/AU - 1.0/AtU); + fvScalarMatrix pEqn ( fvm::div(phid, p) @@ -58,11 +58,15 @@ else surfaceScalarField phiHbyA ( "phiHbyA", - fvc::interpolate(rho*HbyA) & mesh.Sf() + fvc::interpolate(rho)*fvc::interpolate(HbyA) & mesh.Sf() ); - closedVolume = adjustPhi(phi, U, p); - phi += fvc::interpolate(rho/AtU - rho/AU)*fvc::snGrad(p)*mesh.magSf(); + closedVolume = adjustPhi(phiHbyA, U, p); + + phiHbyA += + fvc::interpolate(rho/AtU - rho/AU)*fvc::snGrad(p)*mesh.magSf(); + + HbyA -= fvc::grad(p)*(1.0/AU - 1.0/AtU); fvScalarMatrix pEqn ( @@ -88,9 +92,7 @@ else // Explicitly relax pressure for momentum corrector p.relax(); -U = HbyA - (fvc::grad(p0)*(1.0/AU - 1.0/AtU) + fvc::grad(p)/AtU); -//U = HbyA - fvc::grad(p)/AU; - +U = HbyA - fvc::grad(p)/AtU; U.correctBoundaryConditions(); // For closed-volume cases adjust the pressure and density levels diff --git a/src/OpenFOAM/matrices/LduMatrix/LduMatrix/SolverPerformance.C b/src/OpenFOAM/matrices/LduMatrix/LduMatrix/SolverPerformance.C index 62ce85abb7..5f307718eb 100644 --- a/src/OpenFOAM/matrices/LduMatrix/LduMatrix/SolverPerformance.C +++ b/src/OpenFOAM/matrices/LduMatrix/LduMatrix/SolverPerformance.C @@ -99,8 +99,16 @@ void Foam::SolverPerformance::print { for(direction cmpt=0; cmpt::nComponents; cmpt++) { - os << solverName_ << ": Solving for " - << word(fieldName_ + pTraits::componentNames[cmpt]); + if (pTraits::nComponents == 1) + { + os << solverName_ << ": Solving for " << fieldName_; + + } + else + { + os << solverName_ << ": Solving for " + << word(fieldName_ + pTraits::componentNames[cmpt]); + } if (singular_[cmpt]) { diff --git a/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixATmul.C b/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixATmul.C index db2b6d8bbb..61ba90d830 100644 --- a/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixATmul.C +++ b/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixATmul.C @@ -295,4 +295,36 @@ Foam::tmp Foam::lduMatrix::residual } +Foam::tmp Foam::lduMatrix::H1() const +{ + tmp tH1 + ( + new scalarField(lduAddr().size(), 0.0) + ); + + if (lowerPtr_ || upperPtr_) + { + scalarField& H1_ = tH1(); + + scalar* __restrict__ H1Ptr = H1_.begin(); + + const label* __restrict__ uPtr = lduAddr().upperAddr().begin(); + const label* __restrict__ lPtr = lduAddr().lowerAddr().begin(); + + const scalar* __restrict__ lowerPtr = lower().begin(); + const scalar* __restrict__ upperPtr = upper().begin(); + + register const label nFaces = upper().size(); + + for (register label face=0; face Foam::lduMatrix::H1() const -{ - tmp tH1 - ( - new scalarField(lduAddr().size(), 0.0) - ); - - if (lowerPtr_ || upperPtr_) - { - scalarField& H1_ = tH1(); - - scalar* __restrict__ H1Ptr = H1_.begin(); - - const label* __restrict__ uPtr = lduAddr().upperAddr().begin(); - const label* __restrict__ lPtr = lduAddr().lowerAddr().begin(); - - const scalar* __restrict__ lowerPtr = lower().begin(); - const scalar* __restrict__ upperPtr = upper().begin(); - - register const label nFaces = upper().size(); - - for (register label face=0; face Foam::fvMatrix::H1() const ); volScalarField& H1_ = tH1(); - // Loop over field components - /* - for (direction cmpt=0; cmpt& ptf = psi_.boundaryField()[patchI]; + + if (ptf.coupled() && ptf.size()) + { + addToInternalField + ( + lduAddr().patchAddr(patchI), + boundaryCoeffs_[patchI].component(0), + H1_ + ); + } + } + H1_.internalField() /= psi_.mesh().V(); H1_.correctBoundaryConditions(); diff --git a/tutorials/compressible/rhoSimplecFoam/squareBend/0/U b/tutorials/compressible/rhoSimplecFoam/squareBend/0/U index 70988607c8..d5bd0a74cf 100644 --- a/tutorials/compressible/rhoSimplecFoam/squareBend/0/U +++ b/tutorials/compressible/rhoSimplecFoam/squareBend/0/U @@ -28,7 +28,7 @@ boundaryField inlet { type flowRateInletVelocity; - flowRate constant 0.5; //0.75; + flowRate constant 0.5; value uniform (0 0 0); } outlet diff --git a/tutorials/compressible/rhoSimplecFoam/squareBend/system/decomposeParDict b/tutorials/compressible/rhoSimplecFoam/squareBend/system/decomposeParDict new file mode 100644 index 0000000000..9e844a62ba --- /dev/null +++ b/tutorials/compressible/rhoSimplecFoam/squareBend/system/decomposeParDict @@ -0,0 +1,45 @@ +/*--------------------------------*- C++ -*----------------------------------*\ +| ========= | | +| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | +| \\ / O peration | Version: dev | +| \\ / A nd | Web: www.OpenFOAM.org | +| \\/ M anipulation | | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class dictionary; + location "system"; + object decomposeParDict; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +numberOfSubdomains 8; + +method hierarchical; + +simpleCoeffs +{ + n (8 1 1); + delta 0.001; +} + +hierarchicalCoeffs +{ + n (4 2 1); + delta 0.001; + order xyz; +} + +manualCoeffs +{ + dataFile ""; +} + +distributed no; + +roots ( ); + + +// ************************************************************************* //