From 063936cfe00692d1a3c425f83992d8427a85b46b Mon Sep 17 00:00:00 2001 From: andy Date: Wed, 22 Oct 2008 16:13:37 +0100 Subject: [PATCH 01/16] bugfix: multiple edge gradings per block --- .../mesh/generation/blockMesh/blockPoints.C | 228 +++++++----------- 1 file changed, 86 insertions(+), 142 deletions(-) diff --git a/applications/utilities/mesh/generation/blockMesh/blockPoints.C b/applications/utilities/mesh/generation/blockMesh/blockPoints.C index 0cb024215d..2fecb21bf5 100644 --- a/applications/utilities/mesh/generation/blockMesh/blockPoints.C +++ b/applications/utilities/mesh/generation/blockMesh/blockPoints.C @@ -41,23 +41,23 @@ namespace Foam void block::blockPoints() { // set local variables for mesh specification - label ni = blockDef_.n().x(); - label nj = blockDef_.n().y(); - label nk = blockDef_.n().z(); + const label ni = blockDef_.n().x(); + const label nj = blockDef_.n().y(); + const label nk = blockDef_.n().z(); - point start = blockDef_.points()[blockDef_.blockShape()[0]]; - point xEnd = blockDef_.points()[blockDef_.blockShape()[1]]; - point xyEnd = blockDef_.points()[blockDef_.blockShape()[2]]; - point yEnd = blockDef_.points()[blockDef_.blockShape()[3]]; + const point p000 = blockDef_.points()[blockDef_.blockShape()[0]]; + const point p100 = blockDef_.points()[blockDef_.blockShape()[1]]; + const point p110 = blockDef_.points()[blockDef_.blockShape()[2]]; + const point p010 = blockDef_.points()[blockDef_.blockShape()[3]]; - point zEnd = blockDef_.points()[blockDef_.blockShape()[4]]; - point xzEnd = blockDef_.points()[blockDef_.blockShape()[5]]; - point xyzEnd = blockDef_.points()[blockDef_.blockShape()[6]]; - point yzEnd = blockDef_.points()[blockDef_.blockShape()[7]]; + const point p001 = blockDef_.points()[blockDef_.blockShape()[4]]; + const point p101 = blockDef_.points()[blockDef_.blockShape()[5]]; + const point p111 = blockDef_.points()[blockDef_.blockShape()[6]]; + const point p011 = blockDef_.points()[blockDef_.blockShape()[7]]; - // set reference to the list of edge point and weighting factors - const List >& edgePoints = blockDef_.blockEdgePoints(); - const scalarListList& edgeWeights = blockDef_.blockEdgeWeights(); + // list of edge point and weighting factors + const List >& p = blockDef_.blockEdgePoints(); + const scalarListList& w = blockDef_.blockEdgeWeights(); // generate vertices @@ -69,193 +69,138 @@ void block::blockPoints() { label vertexNo = vtxLabel(i, j, k); - vector edgex1 = start*(1.0 - edgeWeights[0][i]) - + xEnd*edgeWeights[0][i]; + // points on edges + vector edgex1 = p000 + (p100 - p000)*w[0][i]; + vector edgex2 = p010 + (p110 - p010)*w[1][i]; + vector edgex3 = p011 + (p111 - p011)*w[2][i]; + vector edgex4 = p001 + (p101 - p001)*w[3][i]; - vector edgex2 = yEnd*(1.0 - edgeWeights[1][i]) - + xyEnd*edgeWeights[1][i]; + vector edgey1 = p000 + (p010 - p000)*w[4][j]; + vector edgey2 = p100 + (p110 - p100)*w[5][j]; + vector edgey3 = p101 + (p111 - p101)*w[6][j]; + vector edgey4 = p001 + (p011 - p001)*w[7][j]; - vector edgex3 = yzEnd*(1.0 - edgeWeights[2][i]) - + xyzEnd*edgeWeights[2][i]; - - vector edgex4 = zEnd*(1.0 - edgeWeights[3][i]) - + xzEnd*edgeWeights[3][i]; - - - - vector edgey1 = start*(1.0 - edgeWeights[4][j]) - + yEnd*edgeWeights[4][j]; - - vector edgey2 = xEnd*(1.0 - edgeWeights[5][j]) - + xyEnd*edgeWeights[5][j]; - - vector edgey3 = xzEnd*(1.0 - edgeWeights[6][j]) - + xyzEnd*edgeWeights[6][j]; - - vector edgey4 = zEnd*(1.0 - edgeWeights[7][j]) - + yzEnd*edgeWeights[7][j]; - - - - vector edgez1 = start*(1.0 - edgeWeights[8][k]) - + zEnd*edgeWeights[8][k]; - - vector edgez2 = xEnd*(1.0 - edgeWeights[9][k]) - + xzEnd*edgeWeights[9][k]; - - vector edgez3 = xyEnd*(1.0 - edgeWeights[10][k]) - + xyzEnd*edgeWeights[10][k]; - - vector edgez4 = yEnd*(1.0 - edgeWeights[11][k]) - + yzEnd*edgeWeights[11][k]; + vector edgez1 = p000 + (p001 - p000)*w[8][k]; + vector edgez2 = p100 + (p101 - p100)*w[9][k]; + vector edgez3 = p110 + (p111 - p110)*w[10][k]; + vector edgez4 = p010 + (p011 - p010)*w[11][k]; // calculate the importance factors for all edges // x - direction scalar impx1 = ( - (1.0 - edgeWeights[0][i]) - *(1.0 - edgeWeights[4][j]) - *(1.0 - edgeWeights[8][k]) - + edgeWeights[0][i] - *(1.0 - edgeWeights[5][j]) - *(1.0 - edgeWeights[9][k]) + (1.0 - w[0][i])*(1.0 - w[4][j])*(1.0 - w[8][k]) + + w[0][i]*(1.0 - w[5][j])*(1.0 - w[9][k]) ); scalar impx2 = ( - (1.0 - edgeWeights[1][i]) - *edgeWeights[4][j] - *(1.0 - edgeWeights[11][k]) - + edgeWeights[1][i] - *edgeWeights[5][j] - *(1.0 - edgeWeights[10][k]) + (1.0 - w[1][i])*w[4][j]*(1.0 - w[11][k]) + + w[1][i]*w[5][j]*(1.0 - w[10][k]) ); - scalar impx3 = - ( - (1.0 - edgeWeights[2][i]) - *edgeWeights[7][j] - *edgeWeights[11][k] - + edgeWeights[2][i] - *edgeWeights[6][j] - *edgeWeights[10][k] - ); + scalar impx3 = + ( + (1.0 - w[2][i])*w[7][j]*w[11][k] + + w[2][i]*w[6][j]*w[10][k] + ); - scalar impx4 = - ( - (1.0 - edgeWeights[3][i]) - *(1.0 - edgeWeights[7][j]) - *edgeWeights[8][k] - + edgeWeights[3][i] - *(1.0 - edgeWeights[6][j]) - *edgeWeights[9][k] - ); + scalar impx4 = + ( + (1.0 - w[3][i])*(1.0 - w[7][j])*w[8][k] + + w[3][i]*(1.0 - w[6][j])*w[9][k] + ); + scalar magImpx = impx1 + impx2 + impx3 + impx4; + impx1 /= magImpx; + impx2 /= magImpx; + impx3 /= magImpx; + impx4 /= magImpx; // y - direction scalar impy1 = ( - (1.0 - edgeWeights[4][j]) - *(1.0 - edgeWeights[0][i]) - *(1.0 - edgeWeights[8][k]) - + edgeWeights[4][j] - *(1.0 - edgeWeights[1][i]) - *(1.0 - edgeWeights[11][k]) + (1.0 - w[4][j])*(1.0 - w[0][i])*(1.0 - w[8][k]) + + w[4][j]*(1.0 - w[1][i])*(1.0 - w[11][k]) ); scalar impy2 = ( - (1.0 - edgeWeights[5][j]) - *edgeWeights[0][i] - *(1.0 - edgeWeights[9][k]) - + edgeWeights[5][j] - *edgeWeights[1][i] - *(1.0 - edgeWeights[10][k]) + (1.0 - w[5][j])*w[0][i]*(1.0 - w[9][k]) + + w[5][j]*w[1][i]*(1.0 - w[10][k]) ); scalar impy3 = ( - (1.0 - edgeWeights[6][j]) - *edgeWeights[3][i] - *edgeWeights[9][k] - + edgeWeights[6][j] - *edgeWeights[2][i] - *edgeWeights[10][k] + (1.0 - w[6][j])*w[3][i]*w[9][k] + + w[6][j]*w[2][i]*w[10][k] ); scalar impy4 = ( - (1.0 - edgeWeights[7][j]) - *(1.0 - edgeWeights[3][i]) - *edgeWeights[8][k] - + edgeWeights[7][j] - *(1.0 - edgeWeights[2][i]) - *edgeWeights[11][k] + (1.0 - w[7][j])*(1.0 - w[3][i])*w[8][k] + + w[7][j]*(1.0 - w[2][i])*w[11][k] ); + scalar magImpy = impy1 + impy2 + impy3 + impy4; + impy1 /= magImpy; + impy2 /= magImpy; + impy3 /= magImpy; + impy4 /= magImpy; // z - direction scalar impz1 = ( - (1.0 - edgeWeights[8][k]) - *(1.0 - edgeWeights[0][i]) - *(1.0 - edgeWeights[4][j]) - + edgeWeights[8][k] - *(1.0 - edgeWeights[3][i]) - *(1.0 - edgeWeights[7][j]) + (1.0 - w[8][k])*(1.0 - w[0][i])*(1.0 - w[4][j]) + + w[8][k]*(1.0 - w[3][i])*(1.0 - w[7][j]) ); scalar impz2 = ( - (1.0 - edgeWeights[9][k]) - *edgeWeights[0][i] - *(1.0 - edgeWeights[5][j]) - + edgeWeights[9][k] - *edgeWeights[3][i] - *(1.0 - edgeWeights[6][j]) + (1.0 - w[9][k])*w[0][i]*(1.0 - w[5][j]) + + w[9][k]*w[3][i]*(1.0 - w[6][j]) ); scalar impz3 = ( - (1.0 - edgeWeights[10][k]) - *edgeWeights[1][i] - *edgeWeights[5][j] - + edgeWeights[10][k] - *edgeWeights[2][i] - *edgeWeights[6][j] + (1.0 - w[10][k])*w[1][i]*w[5][j] + + w[10][k]*w[2][i]*w[6][j] ); scalar impz4 = ( - (1.0 - edgeWeights[11][k]) - *(1.0 - edgeWeights[1][i]) - *edgeWeights[4][j] - + edgeWeights[11][k] - *(1.0 - edgeWeights[2][i]) - *edgeWeights[7][j] + (1.0 - w[11][k])*(1.0 - w[1][i])*w[4][j] + + w[11][k]*(1.0 - w[2][i])*w[7][j] ); + scalar magImpz = impz1 + impz2 + impz3 + impz4; + impz1 /= magImpz; + impz2 /= magImpz; + impz3 /= magImpz; + impz4 /= magImpz; + // calculate the correction vectors - vector corx1 = impx1*(edgePoints[0][i] - edgex1); - vector corx2 = impx2*(edgePoints[1][i] - edgex2); - vector corx3 = impx3*(edgePoints[2][i] - edgex3); - vector corx4 = impx4*(edgePoints[3][i] - edgex4); + vector corx1 = impx1*(p[0][i] - edgex1); + vector corx2 = impx2*(p[1][i] - edgex2); + vector corx3 = impx3*(p[2][i] - edgex3); + vector corx4 = impx4*(p[3][i] - edgex4); - vector cory1 = impy1*(edgePoints[4][j] - edgey1); - vector cory2 = impy2*(edgePoints[5][j] - edgey2); - vector cory3 = impy3*(edgePoints[6][j] - edgey3); - vector cory4 = impy4*(edgePoints[7][j] - edgey4); + vector cory1 = impy1*(p[4][j] - edgey1); + vector cory2 = impy2*(p[5][j] - edgey2); + vector cory3 = impy3*(p[6][j] - edgey3); + vector cory4 = impy4*(p[7][j] - edgey4); - vector corz1 = impz1*(edgePoints[8][k] - edgez1); - vector corz2 = impz2*(edgePoints[9][k] - edgez2); - vector corz3 = impz3*(edgePoints[10][k] - edgez3); - vector corz4 = impz4*(edgePoints[11][k] - edgez4); + vector corz1 = impz1*(p[8][k] - edgez1); + vector corz2 = impz2*(p[9][k] - edgez2); + vector corz3 = impz3*(p[10][k] - edgez3); + vector corz4 = impz4*(p[11][k] - edgez4); // multiply by the importance factor + // x - direction edgex1 *= impx1; edgex2 *= impx2; @@ -285,7 +230,6 @@ void block::blockPoints() vertices_[vertexNo] += corx1 + corx2 + corx3 + corx4; vertices_[vertexNo] += cory1 + cory2 + cory3 + cory4; vertices_[vertexNo] += corz1 + corz2 + corz3 + corz4; - } } } From 9ea84ab62e58a97f674be76d87932c2d995c3001 Mon Sep 17 00:00:00 2001 From: andy Date: Thu, 23 Oct 2008 16:59:27 +0100 Subject: [PATCH 02/16] added error msg when sampling patch is not found --- .../directMappedFixedValueFvPatchField.C | 10 ++++++++++ 1 file changed, 10 insertions(+) diff --git a/src/finiteVolume/fields/fvPatchFields/derived/directMappedFixedValue/directMappedFixedValueFvPatchField.C b/src/finiteVolume/fields/fvPatchFields/derived/directMappedFixedValue/directMappedFixedValueFvPatchField.C index 279e72a0c5..ab3af1717e 100644 --- a/src/finiteVolume/fields/fvPatchFields/derived/directMappedFixedValue/directMappedFixedValueFvPatchField.C +++ b/src/finiteVolume/fields/fvPatchFields/derived/directMappedFixedValue/directMappedFixedValueFvPatchField.C @@ -234,6 +234,16 @@ void directMappedFixedValueFvPatchField::updateCoeffs() ( mpp.samplePatch() ); + if (patchID < 0) + { + FatalErrorIn + ( + "void directMappedFixedValueFvPatchField::" + "updateCoeffs()" + )<< "Unable to find sample patch " << mpp.samplePatch() + << " for patch " << this->patch().name() << nl + << abort(FatalError); + } typedef GeometricField fieldType; const word& fieldName = this->dimensionedInternalField().name(); const fieldType& sendField = From 1fb786eb6af47737f0ab2bbee5670f3acfcf7bd1 Mon Sep 17 00:00:00 2001 From: mattijs Date: Thu, 30 Oct 2008 15:40:31 +0000 Subject: [PATCH 03/16] removed transfer; delegate to List instead --- src/OpenFOAM/fields/Fields/Field/Field.C | 14 -------------- src/OpenFOAM/fields/Fields/Field/Field.H | 6 ------ 2 files changed, 20 deletions(-) diff --git a/src/OpenFOAM/fields/Fields/Field/Field.C b/src/OpenFOAM/fields/Fields/Field/Field.C index 994f8ee630..c906bcbf55 100644 --- a/src/OpenFOAM/fields/Fields/Field/Field.C +++ b/src/OpenFOAM/fields/Fields/Field/Field.C @@ -575,20 +575,6 @@ void Field::replace } -template -void Field::transfer(Field& f) -{ - List::transfer(f); -} - - -template -void Field::transfer(List& lst) -{ - List::transfer(lst); -} - - template tmp > Field::T() const { diff --git a/src/OpenFOAM/fields/Fields/Field/Field.H b/src/OpenFOAM/fields/Fields/Field/Field.H index 8e168cb0ac..1b372925e3 100644 --- a/src/OpenFOAM/fields/Fields/Field/Field.H +++ b/src/OpenFOAM/fields/Fields/Field/Field.H @@ -297,12 +297,6 @@ public: //- Replace a component field of the field void replace(const direction, const cmptType&); - //- Transfer the contents of the argument Field into this Field - void transfer(Field&); - - //- Transfer the contents of the argument List into this Field - void transfer(List&); - //- Return the field transpose (only defined for second rank tensors) tmp > T() const; From ded4f0b2a4be7a2e1e7433360f637845b69587d2 Mon Sep 17 00:00:00 2001 From: andy Date: Thu, 30 Oct 2008 18:49:14 +0000 Subject: [PATCH 04/16] added option to set reference by point --- .../general/findRefCell/findRefCell.C | 59 ++++++++++++++++++- .../general/findRefCell/findRefCell.H | 3 +- .../fvMatrices/fvMatrix/fvMatrix.C | 2 +- 3 files changed, 59 insertions(+), 5 deletions(-) diff --git a/src/finiteVolume/cfdTools/general/findRefCell/findRefCell.C b/src/finiteVolume/cfdTools/general/findRefCell/findRefCell.C index 00ce11355a..321eeba6a8 100644 --- a/src/finiteVolume/cfdTools/general/findRefCell/findRefCell.C +++ b/src/finiteVolume/cfdTools/general/findRefCell/findRefCell.C @@ -33,15 +33,68 @@ void Foam::setRefCell const volScalarField& field, const dictionary& dict, label& refCelli, - scalar& refValue + scalar& refValue, + bool forceReference ) { - if (field.needReference()) + if (field.needReference() || forceReference) { word refCellName = field.name() + "RefCell"; + word refPointName = field.name() + "RefPoint"; + word refValueName = field.name() + "RefValue"; - refCelli = readLabel(dict.lookup(refCellName)); + if (dict.found(refCellName)) + { + if (Pstream::master()) + { + refCelli = readLabel(dict.lookup(refCellName)); + } + else + { + refCelli = -1; + } + } + else if (dict.found(refPointName)) + { + point refPointi(dict.lookup(refPointName)); + refCelli = field.mesh().findCell(refPointi); + label hasRef = (refCelli >= 0 ? 1 : 0); + label sumHasRef = returnReduce