From ccb38822ea566887eb5e783c8020feccbf25c161 Mon Sep 17 00:00:00 2001 From: mattijs Date: Mon, 8 Apr 2013 16:29:45 +0100 Subject: [PATCH 01/16] BUG: foamToEnsight: using 3 digit extension --- .../postProcessing/dataConversion/foamToEnsight/ensightMesh.C | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/applications/utilities/postProcessing/dataConversion/foamToEnsight/ensightMesh.C b/applications/utilities/postProcessing/dataConversion/foamToEnsight/ensightMesh.C index df6fe4da41..c90e2319e8 100644 --- a/applications/utilities/postProcessing/dataConversion/foamToEnsight/ensightMesh.C +++ b/applications/utilities/postProcessing/dataConversion/foamToEnsight/ensightMesh.C @@ -1026,7 +1026,7 @@ void Foam::ensightMesh::write if (timeIndex == 0) { - timeFile += "000."; + timeFile += "0000."; } else if (meshMoving) { From 6f2e475ba7cf2f2345f32519b2a124d513305704 Mon Sep 17 00:00:00 2001 From: Sergio Ferraris Date: Tue, 9 Apr 2013 09:45:38 +0100 Subject: [PATCH 02/16] BUG: Adding entry to the table for solidChemistryModels and correcting chemistryThermo entry in chemistryProperties for pyrolisis --- .../basicSolidChemistryModelNew.C | 8 +++++--- .../constant/panelRegion/chemistryProperties | 2 +- 2 files changed, 6 insertions(+), 4 deletions(-) diff --git a/src/thermophysicalModels/solidChemistryModel/basicSolidChemistryModel/basicSolidChemistryModelNew.C b/src/thermophysicalModels/solidChemistryModel/basicSolidChemistryModel/basicSolidChemistryModelNew.C index adf318299f..523fff0305 100644 --- a/src/thermophysicalModels/solidChemistryModel/basicSolidChemistryModel/basicSolidChemistryModelNew.C +++ b/src/thermophysicalModels/solidChemistryModel/basicSolidChemistryModel/basicSolidChemistryModelNew.C @@ -53,11 +53,12 @@ New Info<< "Selecting chemistry type " << chemistryTypeDict << endl; - const int nCmpt = 12; + const int nCmpt = 13; const char* cmptNames[nCmpt] = { "chemistrySolver", "chemistryThermo", + "baseChemistry", "transport", "thermo", "equationOfState", @@ -107,8 +108,9 @@ New word chemistryTypeName ( word(chemistryTypeDict.lookup("chemistrySolver")) + '<' - + word(chemistryTypeDict.lookup("chemistryThermo")) + ',' - + solidThermoTypeName + ',' + gasThermoTypeName + ">" + + word(chemistryTypeDict.lookup("chemistryThermo")) + '<' + + typeName + ',' + + solidThermoTypeName + ',' + gasThermoTypeName + ">>" ); Info<< "chemistryTypeName " << chemistryTypeName << endl; diff --git a/tutorials/combustion/fireFoam/les/oppositeBurningPanels/constant/panelRegion/chemistryProperties b/tutorials/combustion/fireFoam/les/oppositeBurningPanels/constant/panelRegion/chemistryProperties index 9f8ecdd196..dd30353e69 100644 --- a/tutorials/combustion/fireFoam/les/oppositeBurningPanels/constant/panelRegion/chemistryProperties +++ b/tutorials/combustion/fireFoam/les/oppositeBurningPanels/constant/panelRegion/chemistryProperties @@ -18,7 +18,7 @@ FoamFile chemistryType { chemistrySolver ode; - chemistryThermo pyrolysis; + chemistryThermo pyrolysisChemistryModel; } chemistry on; From b5dc0921678b5a9de413615daec449bd0c8b4a9e Mon Sep 17 00:00:00 2001 From: Henry Date: Tue, 9 Apr 2013 10:49:47 +0100 Subject: [PATCH 03/16] Minor formatting --- .../schemes/FitData/FitData.C | 18 +++++++++--------- 1 file changed, 9 insertions(+), 9 deletions(-) diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/FitData/FitData.C b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/FitData/FitData.C index 5b227e12eb..c1d322ff4e 100644 --- a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/FitData/FitData.C +++ b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/FitData/FitData.C @@ -45,11 +45,11 @@ Foam::FitData::FitData linearCorrection_(linearCorrection), linearLimitFactor_(linearLimitFactor), centralWeight_(centralWeight), -# ifdef SPHERICAL_GEOMETRY + #ifdef SPHERICAL_GEOMETRY dim_(2), -# else + #else dim_(mesh.nGeometricD()), -# endif + #endif minSize_(Polynomial::nTerms(dim_)) { // Check input @@ -79,7 +79,7 @@ void Foam::FitData::findFaceDirs idir = mesh.faceAreas()[facei]; idir /= mag(idir); -# ifndef SPHERICAL_GEOMETRY + #ifndef SPHERICAL_GEOMETRY if (mesh.nGeometricD() <= 2) // find the normal direction { if (mesh.geometricD()[0] == -1) @@ -100,10 +100,10 @@ void Foam::FitData::findFaceDirs const face& f = mesh.faces()[facei]; kdir = mesh.points()[f[0]] - mesh.faceCentres()[facei]; } -# else + #else // Spherical geometry so kdir is the radial direction kdir = mesh.faceCentres()[facei]; -# endif + #endif if (mesh.nGeometricD() == 3) { @@ -170,11 +170,11 @@ void Foam::FitData::calcFit d.x() = (p - p0)&idir; d.y() = (p - p0)&jdir; -# ifndef SPHERICAL_GEOMETRY + #ifndef SPHERICAL_GEOMETRY d.z() = (p - p0)&kdir; -# else + #else d.z() = mag(p) - mag(p0); -# endif + #endif if (ip == 0) { From 8629dca1f5d3d403a69c6e0876c1a6beb22bd23e Mon Sep 17 00:00:00 2001 From: Henry Date: Tue, 9 Apr 2013 10:50:32 +0100 Subject: [PATCH 04/16] P1 radiation model: only evaluate Qr on non-coupled patches --- .../radiationModels/radiationModel/P1/P1.C | 10 +++++++--- 1 file changed, 7 insertions(+), 3 deletions(-) diff --git a/src/thermophysicalModels/radiationModels/radiationModel/P1/P1.C b/src/thermophysicalModels/radiationModels/radiationModel/P1/P1.C index a07142897b..cf7a2329d4 100644 --- a/src/thermophysicalModels/radiationModels/radiationModel/P1/P1.C +++ b/src/thermophysicalModels/radiationModels/radiationModel/P1/P1.C @@ -241,10 +241,14 @@ void Foam::radiation::P1::calculate() ); // Calculate radiative heat flux on boundaries. - forAll(mesh_.boundaryMesh(), patchI) + forAll(mesh_.boundaryMesh(), patchi) { - Qr_.boundaryField()[patchI] = - -gamma.boundaryField()[patchI]*G_.boundaryField()[patchI].snGrad(); + if (!G_.boundaryField()[patchi].coupled()) + { + Qr_.boundaryField()[patchi] = + -gamma.boundaryField()[patchi] + *G_.boundaryField()[patchi].snGrad(); + } } } From 928cef3003b7fa173a304c65af1a472c217fa692 Mon Sep 17 00:00:00 2001 From: Henry Date: Tue, 9 Apr 2013 11:40:23 +0100 Subject: [PATCH 05/16] incompressibleTwoPhaseMixture: Correct the phase name lookup in read --- .../incompressibleTwoPhaseMixture.C | 10 ++++++++-- 1 file changed, 8 insertions(+), 2 deletions(-) diff --git a/src/transportModels/incompressible/incompressibleTwoPhaseMixture/incompressibleTwoPhaseMixture.C b/src/transportModels/incompressible/incompressibleTwoPhaseMixture/incompressibleTwoPhaseMixture.C index 98046c9b76..397f9566b1 100644 --- a/src/transportModels/incompressible/incompressibleTwoPhaseMixture/incompressibleTwoPhaseMixture.C +++ b/src/transportModels/incompressible/incompressibleTwoPhaseMixture/incompressibleTwoPhaseMixture.C @@ -175,8 +175,14 @@ bool Foam::incompressibleTwoPhaseMixture::read() { if ( - nuModel1_().read(subDict(phase1Name_)) - && nuModel2_().read(subDict(phase2Name_)) + nuModel1_().read + ( + subDict(phase1Name_ == "1" ? "phase1": phase1Name_) + ) + && nuModel2_().read + ( + subDict(phase2Name_ == "2" ? "phase2": phase2Name_) + ) ) { nuModel1_->viscosityProperties().lookup("rho") >> rho1_; From b15cf81217be3d6ca203a6bca9df41fc94914d32 Mon Sep 17 00:00:00 2001 From: mattijs Date: Tue, 9 Apr 2013 12:09:55 +0100 Subject: [PATCH 06/16] ENH: autoSnapDriverFeature: protect for zero sized faces --- .../autoHexMesh/autoHexMeshDriver/autoSnapDriverFeature.C | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoSnapDriverFeature.C b/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoSnapDriverFeature.C index cecee33cef..37e61cbbb7 100644 --- a/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoSnapDriverFeature.C +++ b/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoSnapDriverFeature.C @@ -767,7 +767,7 @@ void Foam::autoSnapDriver::binFeatureFace ) const { // What to do with very far attraction? For now just ignore the face - if (magSqr(faceDisp) < sqr(snapDist)) + if (magSqr(faceDisp) < sqr(snapDist) && mag(faceSurfaceNormal) > VSMALL) { const point pt = fc + faceDisp; From 876d6b52d598acbdc152481cc022304d281e471b Mon Sep 17 00:00:00 2001 From: mattijs Date: Tue, 9 Apr 2013 12:10:58 +0100 Subject: [PATCH 07/16] BUG: incompressibleTwoPhaseMixture: re-reading --- .../incompressibleTwoPhaseMixture.C | 10 ++++++++-- 1 file changed, 8 insertions(+), 2 deletions(-) diff --git a/src/transportModels/incompressible/incompressibleTwoPhaseMixture/incompressibleTwoPhaseMixture.C b/src/transportModels/incompressible/incompressibleTwoPhaseMixture/incompressibleTwoPhaseMixture.C index 98046c9b76..397f9566b1 100644 --- a/src/transportModels/incompressible/incompressibleTwoPhaseMixture/incompressibleTwoPhaseMixture.C +++ b/src/transportModels/incompressible/incompressibleTwoPhaseMixture/incompressibleTwoPhaseMixture.C @@ -175,8 +175,14 @@ bool Foam::incompressibleTwoPhaseMixture::read() { if ( - nuModel1_().read(subDict(phase1Name_)) - && nuModel2_().read(subDict(phase2Name_)) + nuModel1_().read + ( + subDict(phase1Name_ == "1" ? "phase1": phase1Name_) + ) + && nuModel2_().read + ( + subDict(phase2Name_ == "2" ? "phase2": phase2Name_) + ) ) { nuModel1_->viscosityProperties().lookup("rho") >> rho1_; From 142d9333c76334023383e871d433363f501f2b0a Mon Sep 17 00:00:00 2001 From: Henry Date: Tue, 9 Apr 2013 14:28:54 +0100 Subject: [PATCH 08/16] CentredFitSnGrad: New fit-based snGrad --- src/finiteVolume/Make/files | 1 + .../CentredFitSnGrad/CentredFitSnGradData.C | 272 ++++++++++++++++++ .../CentredFitSnGrad/CentredFitSnGradData.H | 126 ++++++++ .../CentredFitSnGrad/CentredFitSnGradScheme.H | 182 ++++++++++++ .../quadraticFitSnGrad/quadraticFitSnGrads.C | 51 ++++ .../schemes/FitData/FitData.H | 32 ++- 6 files changed, 658 insertions(+), 6 deletions(-) create mode 100644 src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradData.C create mode 100644 src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradData.H create mode 100644 src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradScheme.H create mode 100644 src/finiteVolume/finiteVolume/snGradSchemes/quadraticFitSnGrad/quadraticFitSnGrads.C diff --git a/src/finiteVolume/Make/files b/src/finiteVolume/Make/files index 1ce1f78949..2a650f784f 100644 --- a/src/finiteVolume/Make/files +++ b/src/finiteVolume/Make/files @@ -356,6 +356,7 @@ $(snGradSchemes)/faceCorrectedSnGrad/faceCorrectedSnGrads.C $(snGradSchemes)/limitedSnGrad/limitedSnGrads.C $(snGradSchemes)/uncorrectedSnGrad/uncorrectedSnGrads.C $(snGradSchemes)/orthogonalSnGrad/orthogonalSnGrads.C +$(snGradSchemes)/quadraticFitSnGrad/quadraticFitSnGrads.C convectionSchemes = finiteVolume/convectionSchemes $(convectionSchemes)/convectionScheme/convectionSchemes.C diff --git a/src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradData.C b/src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradData.C new file mode 100644 index 0000000000..135deac4bf --- /dev/null +++ b/src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradData.C @@ -0,0 +1,272 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2013 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 . + +\*---------------------------------------------------------------------------*/ + +#include "CentredFitSnGradData.H" +#include "surfaceFields.H" +#include "volFields.H" +#include "SVD.H" +#include "syncTools.H" +#include "extendedCentredCellToFaceStencil.H" + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +template +Foam::CentredFitSnGradData::CentredFitSnGradData +( + const fvMesh& mesh, + const extendedCentredCellToFaceStencil& stencil, + const scalar linearLimitFactor, + const scalar centralWeight +) +: + FitData + < + CentredFitSnGradData, + extendedCentredCellToFaceStencil, + Polynomial + > + ( + mesh, stencil, true, linearLimitFactor, centralWeight + ), + coeffs_(mesh.nFaces()) +{ + if (debug) + { + Info<< "Contructing CentredFitSnGradData" << endl; + } + + calcFit(); + + if (debug) + { + Info<< "CentredFitSnGradData::CentredFitSnGradData() :" + << "Finished constructing polynomialFit data" + << endl; + } +} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +template +void Foam::CentredFitSnGradData::calcFit +( + scalarList& coeffsi, + const List& C, + const scalar wLin, + const scalar deltaCoeff, + const label facei +) +{ + vector idir(1,0,0); + vector jdir(0,1,0); + vector kdir(0,0,1); + this->findFaceDirs(idir, jdir, kdir, facei); + + // Setup the point weights + scalarList wts(C.size(), scalar(1)); + wts[0] = this->centralWeight(); + wts[1] = this->centralWeight(); + + // Reference point + point p0 = this->mesh().faceCentres()[facei]; + + // p0 -> p vector in the face-local coordinate system + vector d; + + // Local coordinate scaling + scalar scale = 1; + + // Matrix of the polynomial components + scalarRectangularMatrix B(C.size(), this->minSize(), scalar(0)); + + forAll(C, ip) + { + const point& p = C[ip]; + const vector p0p = p - p0; + + d.x() = p0p & idir; + d.y() = p0p & jdir; + d.z() = p0p & kdir; + + if (ip == 0) + { + scale = cmptMax(cmptMag((d))); + } + + // Scale the radius vector + d /= scale; + + Polynomial::addCoeffs(B[ip], d, wts[ip], this->dim()); + } + + // Additional weighting for constant and linear terms + for (label i = 0; i < B.n(); i++) + { + B[i][0] *= wts[0]; + B[i][1] *= wts[0]; + } + + // Set the fit + label stencilSize = C.size(); + coeffsi.setSize(stencilSize); + + bool goodFit = false; + for (int iIt = 0; iIt < 8 && !goodFit; iIt++) + { + SVD svd(B, SMALL); + + for (label i=0; ilinearLimitFactor()*wLin) + && (mag(wts[0]*wts[1]*svd.VSinvUt()[0][1] - (1 - wLin) + ) < this->linearLimitFactor()*(1 - wLin)) + && coeffsi[0] < 0 && coeffsi[1] > 0 + && mag(coeffsi[0] + deltaCoeff) < 0.5*deltaCoeff + && mag(coeffsi[1] - deltaCoeff) < 0.5*deltaCoeff; + + if (!goodFit) + { + // (not good fit so increase weight in the centre and weight + // for constant and linear terms) + + WarningIn + ( + "CentredFitSnGradData::calcFit" + "(const List& C, const label facei" + ) << "Cannot fit face " << facei << " iteration " << iIt + << " with sum of weights " << sum(coeffsi) << nl + << " Weights " << coeffsi << nl + << " Linear weights " << wLin << " " << 1 - wLin << nl + << " deltaCoeff " << deltaCoeff << nl + << " sing vals " << svd.S() << nl + << "Components of goodFit:\n" + << " wts[0]*wts[0]*svd.VSinvUt()[0][0] = " + << wts[0]*wts[0]*svd.VSinvUt()[0][0] << nl + << " wts[0]*wts[1]*svd.VSinvUt()[0][1] = " + << wts[0]*wts[1]*svd.VSinvUt()[0][1] + << " dim = " << this->dim() << endl; + + wts[0] *= 10; + wts[1] *= 10; + + for (label j = 0; j < B.m(); j++) + { + B[0][j] *= 10; + B[1][j] *= 10; + } + + for (label i = 0; i < B.n(); i++) + { + B[i][0] *= 10; + B[i][1] *= 10; + } + } + } + + if (goodFit) + { + // Remove the uncorrected coefficients + coeffsi[0] += deltaCoeff; + coeffsi[1] -= deltaCoeff; + } + else + { + WarningIn + ( + "CentredFitSnGradData::calcFit(..)" + ) << "Could not fit face " << facei + << " Coefficients = " << coeffsi + << ", reverting to uncorrected." << endl; + + coeffsi = 0; + } +} + + +template +void Foam::CentredFitSnGradData::calcFit() +{ + const fvMesh& mesh = this->mesh(); + + // Get the cell/face centres in stencil order. + // Centred face stencils no good for triangles or tets. + // Need bigger stencils + List > stencilPoints(mesh.nFaces()); + this->stencil().collectData(mesh.C(), stencilPoints); + + // find the fit coefficients for every face in the mesh + + const surfaceScalarField& w = mesh.surfaceInterpolation::weights(); + const surfaceScalarField& dC = mesh.deltaCoeffs(); + + for (label facei = 0; facei < mesh.nInternalFaces(); facei++) + { + calcFit + ( + coeffs_[facei], + stencilPoints[facei], + w[facei], + dC[facei], + facei + ); + } + + const surfaceScalarField::GeometricBoundaryField& bw = w.boundaryField(); + const surfaceScalarField::GeometricBoundaryField& bdC = dC.boundaryField(); + + forAll(bw, patchi) + { + const fvsPatchScalarField& pw = bw[patchi]; + const fvsPatchScalarField& pdC = bdC[patchi]; + + if (pw.coupled()) + { + label facei = pw.patch().start(); + + forAll(pw, i) + { + calcFit + ( + coeffs_[facei], + stencilPoints[facei], + pw[i], + pdC[i], + facei + ); + facei++; + } + } + } +} + + +// ************************************************************************* // diff --git a/src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradData.H b/src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradData.H new file mode 100644 index 0000000000..c025a17dd5 --- /dev/null +++ b/src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradData.H @@ -0,0 +1,126 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2013 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 . + +Class + Foam::CentredFitSnGradData + +Description + Data for the quadratic fit correction interpolation scheme + +SourceFiles + CentredFitSnGradData.C + +\*---------------------------------------------------------------------------*/ + +#ifndef CentredFitSnGradData_H +#define CentredFitSnGradData_H + +#include "FitData.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +class extendedCentredCellToFaceStencil; + +/*---------------------------------------------------------------------------*\ + Class CentredFitSnGradData Declaration +\*---------------------------------------------------------------------------*/ + +template +class CentredFitSnGradData +: + public FitData + < + CentredFitSnGradData, + extendedCentredCellToFaceStencil, + Polynomial + > +{ + // Private data + + //- For each cell in the mesh store the values which multiply the + // values of the stencil to obtain the gradient for each direction + List coeffs_; + +public: + + TypeName("CentredFitSnGradData"); + + + // Constructors + + //- Construct from components + CentredFitSnGradData + ( + const fvMesh& mesh, + const extendedCentredCellToFaceStencil& stencil, + const scalar linearLimitFactor, + const scalar centralWeight + ); + + + //- Destructor + virtual ~CentredFitSnGradData() + {} + + + // Member functions + + //- Return reference to fit coefficients + const List& coeffs() const + { + return coeffs_; + } + + //- Calculate the fit for the specified face and set the coefficients + void calcFit + ( + scalarList& coeffsi, // coefficients to be set + const List&, // Stencil points + const scalar wLin, // Weight for linear approximation (weights + // nearest neighbours) + const scalar deltaCoeff, // uncorrected delta coefficient + const label faci // Current face index + ); + + void calcFit(); +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#ifdef NoRepository +# include "CentredFitSnGradData.C" +#endif + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradScheme.H b/src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradScheme.H new file mode 100644 index 0000000000..ca307811e6 --- /dev/null +++ b/src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradScheme.H @@ -0,0 +1,182 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2013 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 . + +Class + Foam::CentredFitSnGradScheme + +Description + Centred fit snGrad scheme which applies an explicit correction to snGrad + +\*---------------------------------------------------------------------------*/ + +#ifndef CentredFitSnGradScheme_H +#define CentredFitSnGradScheme_H + +#include "CentredFitSnGradData.H" +#include "snGradScheme.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +class fvMesh; + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace fv +{ +/*---------------------------------------------------------------------------*\ + Class CentredFitSnGradSnGradScheme Declaration +\*---------------------------------------------------------------------------*/ + +template +class CentredFitSnGradScheme +: + public snGradScheme +{ + // Private Data + + //- Factor the fit is allowed to deviate from linear. + // This limits the amount of high-order correction and increases + // stability on bad meshes + const scalar linearLimitFactor_; + + //- Weights for central stencil + const scalar centralWeight_; + + + // Private Member Functions + + //- Disallow default bitwise copy construct + CentredFitSnGradScheme(const CentredFitSnGradScheme&); + + //- Disallow default bitwise assignment + void operator=(const CentredFitSnGradScheme&); + + +public: + + //- Runtime type information + TypeName("CentredFitSnGradScheme"); + + + // Constructors + + //- Construct from mesh and Istream + CentredFitSnGradScheme(const fvMesh& mesh, Istream& is) + : + snGradScheme(mesh), + linearLimitFactor_(readScalar(is)), + centralWeight_(1000) + {} + + + // Member Functions + + //- Return the interpolation weighting factors for the given field + virtual tmp deltaCoeffs + ( + const GeometricField& + ) const + { + return this->mesh().deltaCoeffs(); + } + + //- Return true if this scheme uses an explicit correction + virtual bool corrected() const + { + return true; + } + + //- Return the explicit correction to the face-interpolate + virtual tmp > + correction + ( + const GeometricField& vf + ) const + { + const fvMesh& mesh = this->mesh(); + + const extendedCentredCellToFaceStencil& stencil = Stencil::New + ( + mesh + ); + + const CentredFitSnGradData& cfd = + CentredFitSnGradData::New + ( + mesh, + stencil, + linearLimitFactor_, + centralWeight_ + ); + + const List& f = cfd.coeffs(); + + tmp > sft + = stencil.weightedSum(vf, f); + + sft().dimensions() /= dimLength; + + return sft; + } +}; + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace fv + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +// Add the patch constructor functions to the hash tables + +#define makeCentredFitSnGradTypeScheme(SS, POLYNOMIAL, STENCIL, TYPE) \ + \ +typedef CentredFitSnGradScheme \ + CentredFitSnGradScheme##TYPE##POLYNOMIAL##STENCIL##_; \ +defineTemplateTypeNameAndDebugWithName \ + (CentredFitSnGradScheme##TYPE##POLYNOMIAL##STENCIL##_, #SS, 0); \ + \ +snGradScheme::addMeshConstructorToTable \ + > \ + add##SS##STENCIL##TYPE##MeshConstructorToTable_; + +#define makeCentredFitSnGradScheme(SS, POLYNOMIAL, STENCIL) \ + \ +makeCentredFitSnGradTypeScheme(SS,POLYNOMIAL,STENCIL,scalar) \ +makeCentredFitSnGradTypeScheme(SS,POLYNOMIAL,STENCIL,vector) \ +makeCentredFitSnGradTypeScheme(SS,POLYNOMIAL,STENCIL,sphericalTensor) \ +makeCentredFitSnGradTypeScheme(SS,POLYNOMIAL,STENCIL,symmTensor) \ +makeCentredFitSnGradTypeScheme(SS,POLYNOMIAL,STENCIL,tensor) + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/finiteVolume/finiteVolume/snGradSchemes/quadraticFitSnGrad/quadraticFitSnGrads.C b/src/finiteVolume/finiteVolume/snGradSchemes/quadraticFitSnGrad/quadraticFitSnGrads.C new file mode 100644 index 0000000000..63e4ea98fa --- /dev/null +++ b/src/finiteVolume/finiteVolume/snGradSchemes/quadraticFitSnGrad/quadraticFitSnGrads.C @@ -0,0 +1,51 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2013 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 . + +\*---------------------------------------------------------------------------*/ + +#include "CentredFitSnGradScheme.H" +#include "quadraticFitPolynomial.H" +#include "centredCFCCellToFaceStencilObject.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + defineTemplateTypeNameAndDebug + ( + CentredFitSnGradData, + 0 + ); + + namespace fv + { + makeCentredFitSnGradScheme + ( + quadraticFit, + quadraticFitPolynomial, + centredCFCCellToFaceStencilObject + ); + } +} + +// ************************************************************************* // diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/FitData/FitData.H b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/FitData/FitData.H index cee373b6c4..5faeed6eb1 100644 --- a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/FitData/FitData.H +++ b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/FitData/FitData.H @@ -29,7 +29,7 @@ Description The linearCorrection_ determines whether the fit is for a corrected linear scheme (first two coefficients are corrections for owner and neighbour) or a pure upwind scheme (first coefficient is correction for - owner ; weight on face taken as 1). + owner; weight on face taken as 1). SourceFiles FitData.C @@ -80,7 +80,7 @@ class FitData const label minSize_; - // Private Member Functions +protected: //- Find the normal direction (i) and j and k directions for face faci void findFaceDirs @@ -93,9 +93,6 @@ class FitData public: - //TypeName("FitData"); - - // Constructors //- Construct from components @@ -122,6 +119,30 @@ public: return stencil_; } + //- Factor the fit is allowed to deviate from the base scheme + scalar linearLimitFactor() const + { + return linearLimitFactor_; + } + + //- Return weight for central stencil + scalar centralWeight() const + { + return centralWeight_; + } + + //- Dimensionality of the geometry + label dim() const + { + return dim_; + } + + //- Minimum stencil size + label minSize() const + { + return minSize_; + } + bool linearCorrection() const { return linearCorrection_; @@ -140,7 +161,6 @@ public: //- Calculate the fit for all the faces virtual void calcFit() = 0; - //- Recalculate weights (but not stencil) when the mesh moves bool movePoints(); }; From da535a449a7c069171eaf466bea28710e1ec1879 Mon Sep 17 00:00:00 2001 From: Henry Date: Tue, 9 Apr 2013 22:28:40 +0100 Subject: [PATCH 09/16] linearFitSnGrad: Added new snGrad scheme which handles non-orthogonality with a compact molecule --- src/finiteVolume/Make/files | 1 + .../linearFitSnGrad/linearFitSnGrads.C | 51 +++++++++++++++++++ 2 files changed, 52 insertions(+) create mode 100644 src/finiteVolume/finiteVolume/snGradSchemes/linearFitSnGrad/linearFitSnGrads.C diff --git a/src/finiteVolume/Make/files b/src/finiteVolume/Make/files index 2a650f784f..dd37136533 100644 --- a/src/finiteVolume/Make/files +++ b/src/finiteVolume/Make/files @@ -357,6 +357,7 @@ $(snGradSchemes)/limitedSnGrad/limitedSnGrads.C $(snGradSchemes)/uncorrectedSnGrad/uncorrectedSnGrads.C $(snGradSchemes)/orthogonalSnGrad/orthogonalSnGrads.C $(snGradSchemes)/quadraticFitSnGrad/quadraticFitSnGrads.C +$(snGradSchemes)/linearFitSnGrad/linearFitSnGrads.C convectionSchemes = finiteVolume/convectionSchemes $(convectionSchemes)/convectionScheme/convectionSchemes.C diff --git a/src/finiteVolume/finiteVolume/snGradSchemes/linearFitSnGrad/linearFitSnGrads.C b/src/finiteVolume/finiteVolume/snGradSchemes/linearFitSnGrad/linearFitSnGrads.C new file mode 100644 index 0000000000..4335cc8e69 --- /dev/null +++ b/src/finiteVolume/finiteVolume/snGradSchemes/linearFitSnGrad/linearFitSnGrads.C @@ -0,0 +1,51 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2013 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 . + +\*---------------------------------------------------------------------------*/ + +#include "CentredFitSnGradScheme.H" +#include "linearFitPolynomial.H" +#include "centredFECCellToFaceStencilObject.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + defineTemplateTypeNameAndDebug + ( + CentredFitSnGradData, + 0 + ); + + namespace fv + { + makeCentredFitSnGradScheme + ( + linearFit, + linearFitPolynomial, + centredFECCellToFaceStencilObject + ); + } +} + +// ************************************************************************* // From eecbe11085d87204b69651e650d2d01a63535e4d Mon Sep 17 00:00:00 2001 From: Henry Date: Wed, 10 Apr 2013 12:15:47 +0100 Subject: [PATCH 10/16] CentredFitSnGrad: Change deltaCoeffs -> nonOrthDeltaCoeffs for very non-orthogonal meshes --- .../CentredFitSnGrad/CentredFitSnGradData.C | 4 +--- .../CentredFitSnGrad/CentredFitSnGradData.H | 3 ++- .../CentredFitSnGrad/CentredFitSnGradScheme.H | 23 ++++++++++--------- 3 files changed, 15 insertions(+), 15 deletions(-) diff --git a/src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradData.C b/src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradData.C index 135deac4bf..36a363ea1f 100644 --- a/src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradData.C +++ b/src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradData.C @@ -25,9 +25,7 @@ License #include "CentredFitSnGradData.H" #include "surfaceFields.H" -#include "volFields.H" #include "SVD.H" -#include "syncTools.H" #include "extendedCentredCellToFaceStencil.H" // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // @@ -226,7 +224,7 @@ void Foam::CentredFitSnGradData::calcFit() // find the fit coefficients for every face in the mesh const surfaceScalarField& w = mesh.surfaceInterpolation::weights(); - const surfaceScalarField& dC = mesh.deltaCoeffs(); + const surfaceScalarField& dC = mesh.nonOrthDeltaCoeffs(); for (label facei = 0; facei < mesh.nInternalFaces(); facei++) { diff --git a/src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradData.H b/src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradData.H index c025a17dd5..be834684a7 100644 --- a/src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradData.H +++ b/src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradData.H @@ -25,7 +25,7 @@ Class Foam::CentredFitSnGradData Description - Data for the quadratic fit correction interpolation scheme + Data for centred fit snGrad schemes SourceFiles CentredFitSnGradData.C @@ -64,6 +64,7 @@ class CentredFitSnGradData // values of the stencil to obtain the gradient for each direction List coeffs_; + public: TypeName("CentredFitSnGradData"); diff --git a/src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradScheme.H b/src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradScheme.H index ca307811e6..d92478c317 100644 --- a/src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradScheme.H +++ b/src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradScheme.H @@ -100,7 +100,7 @@ public: const GeometricField& ) const { - return this->mesh().deltaCoeffs(); + return this->mesh().nonOrthDeltaCoeffs(); } //- Return true if this scheme uses an explicit correction @@ -124,18 +124,18 @@ public: ); const CentredFitSnGradData& cfd = - CentredFitSnGradData::New - ( - mesh, - stencil, - linearLimitFactor_, - centralWeight_ - ); - - const List& f = cfd.coeffs(); + CentredFitSnGradData::New + ( + mesh, + stencil, + linearLimitFactor_, + centralWeight_ + ); tmp > sft - = stencil.weightedSum(vf, f); + ( + stencil.weightedSum(vf, cfd.coeffs()) + ); sft().dimensions() /= dimLength; @@ -143,6 +143,7 @@ public: } }; + // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // } // End namespace fv From 7e376c96851307ea72ca758748cd7cacc9032a9e Mon Sep 17 00:00:00 2001 From: laurence Date: Wed, 10 Apr 2013 16:35:39 +0100 Subject: [PATCH 11/16] BUG: extendedFeatureEdgeMesh: Swap feature point algorithm for a faster one The population of featurePointFeatureEdges was slow. Added a List called edgeMap that contains the featureEdge index for every edge in a patch (returns -1 if it is not a feature edge). --- .../extendedFeatureEdgeMeshTemplates.C | 39 ++++++++----------- 1 file changed, 16 insertions(+), 23 deletions(-) diff --git a/src/edgeMesh/extendedFeatureEdgeMesh/extendedFeatureEdgeMeshTemplates.C b/src/edgeMesh/extendedFeatureEdgeMesh/extendedFeatureEdgeMeshTemplates.C index 63cef912b1..ecc049995b 100644 --- a/src/edgeMesh/extendedFeatureEdgeMesh/extendedFeatureEdgeMeshTemplates.C +++ b/src/edgeMesh/extendedFeatureEdgeMesh/extendedFeatureEdgeMeshTemplates.C @@ -65,8 +65,7 @@ void Foam::extendedFeatureEdgeMesh::sortPointsAndEdges labelListList featurePointFeatureEdges(nFeatPts); forAll(featurePointFeatureEdges, pI) { - featurePointFeatureEdges[pI] = - labelList(pointEdges[featurePoints[pI]].size(), -1); + featurePointFeatureEdges[pI] = pointEdges[featurePoints[pI]]; } // Mapping between old and new indices, there is entry in the map for each @@ -74,6 +73,10 @@ void Foam::extendedFeatureEdgeMesh::sortPointsAndEdges // >= 0 corresponds to the index labelList pointMap(sFeatLocalPts.size(), -1); + // Mapping between surface edge index and its feature edge index. -1 if it + // is not a feature edge + labelList edgeMap(sFeatEds.size(), -1); + // Noting when the normal of a face has been used so not to duplicate labelList faceMap(surf.size(), -1); @@ -98,6 +101,8 @@ void Foam::extendedFeatureEdgeMesh::sortPointsAndEdges { label sFEI = featureEdges[i]; + edgeMap[sFEI] = i; + const edge& fE(sFeatEds[sFEI]); // Check to see if the points have been already used @@ -156,43 +161,31 @@ void Foam::extendedFeatureEdgeMesh::sortPointsAndEdges { regionEdges.append(i); } - - forAll(featurePointFeatureEdges, pI) - { - const labelList& fpfEdges = pointEdges[featurePoints[pI]]; - - labelList& fpfe = featurePointFeatureEdges[pI]; - - forAll(fpfEdges, eI) - { - if (sFEI == fpfEdges[eI]) - { - fpfe[eI] = i; - } - } - } } + // Populate feature point feature edges + DynamicList