View Issue Details
ID | Project | Category | View Status | Date Submitted | Last Update |
---|---|---|---|---|---|
0002561 | OpenFOAM | Bug | public | 2017-05-24 09:52 | 2022-05-18 11:14 |
Reporter | niklas.wikstrom | Assigned To | will | ||
Priority | normal | Severity | crash | Reproducibility | always |
Status | resolved | Resolution | fixed | ||
Platform | GNU/Linux | OS | Other | OS Version | (please specify) |
Product Version | dev | ||||
Fixed in Version | dev | ||||
Summary | 0002561: cyclicACMI makes rhoCentral(DyM)Foam crash due to surface interpolation issue | ||||
Description | rhoCentralFoam crashes (sigFpe) immediately in time loop on a cases containing cyclicACMI patches. The cause of the sigFpe seems to be that "interpolate(rho,neg)" produces a zero rho_neg field, which is then used in the denominator to calculate U_neg (around line 88, rhoCentralFoam.C). Thank you. | ||||
Steps To Reproduce | Unpack attached small case. run included makeMesh.sh script run rhoCentralDyMFoam (or rhoCentralFoam) | ||||
Additional Information | Case runs as expected in OpenFOAM 3.x. Case runs in OpenFOAM 4.x/dev using "linear" interpolation. Case runs in OpenFOAM 4.x/dev using e.g. sonicFoam (For completeness: Tested on Fedora-25;gcc-6.3.1 and Centos-7.3;gcc-4.8.5) | ||||
Tags | ACMI, cyclicACMI | ||||
|
|
|
I should perhaps mention that the same issue appears in version 4.x |
|
The issue is in cyclicACMIFvPatchField<Type>::patchNeighbourField(). Uncoupled faces assume a patchNeighbour value of 0. Attached version assumes a patchNeighbour value of the blockage patch (for uncoupled faces). Needs a bit of testing. |
|
cyclicACMIFvPatchField.C (8,252 bytes)
/*---------------------------------------------------------------------------*\ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | \\ / A nd | Copyright (C) 2013-2017 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 <http://www.gnu.org/licenses/>. \*---------------------------------------------------------------------------*/ #include "cyclicACMIFvPatchField.H" #include "transformField.H" // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // template<class Type> Foam::cyclicACMIFvPatchField<Type>::cyclicACMIFvPatchField ( const fvPatch& p, const DimensionedField<Type, volMesh>& iF ) : cyclicACMILduInterfaceField(), coupledFvPatchField<Type>(p, iF), cyclicACMIPatch_(refCast<const cyclicACMIFvPatch>(p)) {} template<class Type> Foam::cyclicACMIFvPatchField<Type>::cyclicACMIFvPatchField ( const fvPatch& p, const DimensionedField<Type, volMesh>& iF, const dictionary& dict ) : cyclicACMILduInterfaceField(), coupledFvPatchField<Type>(p, iF, dict, dict.found("value")), cyclicACMIPatch_(refCast<const cyclicACMIFvPatch>(p)) { if (!isA<cyclicACMIFvPatch>(p)) { FatalIOErrorInFunction ( dict ) << " patch type '" << p.type() << "' not constraint type '" << typeName << "'" << "\n for patch " << p.name() << " of field " << this->internalField().name() << " in file " << this->internalField().objectPath() << exit(FatalIOError); } if (!dict.found("value") && this->coupled()) { this->evaluate(Pstream::commsTypes::blocking); } } template<class Type> Foam::cyclicACMIFvPatchField<Type>::cyclicACMIFvPatchField ( const cyclicACMIFvPatchField<Type>& ptf, const fvPatch& p, const DimensionedField<Type, volMesh>& iF, const fvPatchFieldMapper& mapper ) : cyclicACMILduInterfaceField(), coupledFvPatchField<Type>(ptf, p, iF, mapper), cyclicACMIPatch_(refCast<const cyclicACMIFvPatch>(p)) { if (!isA<cyclicACMIFvPatch>(this->patch())) { FatalErrorInFunction << "' not constraint type '" << typeName << "'" << "\n for patch " << p.name() << " of field " << this->internalField().name() << " in file " << this->internalField().objectPath() << exit(FatalIOError); } } template<class Type> Foam::cyclicACMIFvPatchField<Type>::cyclicACMIFvPatchField ( const cyclicACMIFvPatchField<Type>& ptf ) : cyclicACMILduInterfaceField(), coupledFvPatchField<Type>(ptf), cyclicACMIPatch_(ptf.cyclicACMIPatch_) {} template<class Type> Foam::cyclicACMIFvPatchField<Type>::cyclicACMIFvPatchField ( const cyclicACMIFvPatchField<Type>& ptf, const DimensionedField<Type, volMesh>& iF ) : cyclicACMILduInterfaceField(), coupledFvPatchField<Type>(ptf, iF), cyclicACMIPatch_(ptf.cyclicACMIPatch_) {} // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // template<class Type> bool Foam::cyclicACMIFvPatchField<Type>::coupled() const { return cyclicACMIPatch_.coupled(); } template<class Type> Foam::tmp<Foam::Field<Type>> Foam::cyclicACMIFvPatchField<Type>::patchNeighbourField() const { const Field<Type>& iField = this->primitiveField(); const cyclicACMIPolyPatch& cpp = cyclicACMIPatch_.cyclicACMIPatch(); tmp<Field<Type>> tpnf ( cyclicACMIPatch_.interpolate ( Field<Type> ( iField, cpp.neighbPatch().faceCells() ) ) ); if (doTransform()) { tpnf.ref() = transform(forwardT(), tpnf()); } return cpp.mask()*tpnf + (1.0-cpp.mask())*nonOverlapPatchField(); //return tpnf; } template<class Type> const Foam::cyclicACMIFvPatchField<Type>& Foam::cyclicACMIFvPatchField<Type>::neighbourPatchField() const { const GeometricField<Type, fvPatchField, volMesh>& fld = static_cast<const GeometricField<Type, fvPatchField, volMesh>&> ( this->primitiveField() ); return refCast<const cyclicACMIFvPatchField<Type>> ( fld.boundaryField()[cyclicACMIPatch_.neighbPatchID()] ); } template<class Type> const Foam::fvPatchField<Type>& Foam::cyclicACMIFvPatchField<Type>::nonOverlapPatchField() const { const GeometricField<Type, fvPatchField, volMesh>& fld = static_cast<const GeometricField<Type, fvPatchField, volMesh>&> ( this->primitiveField() ); return fld.boundaryField()[cyclicACMIPatch_.nonOverlapPatchID()]; } template<class Type> void Foam::cyclicACMIFvPatchField<Type>::updateInterfaceMatrix ( scalarField& result, const scalarField& psiInternal, const scalarField& coeffs, const direction cmpt, const Pstream::commsTypes ) const { const cyclicACMIPolyPatch& cpp = cyclicACMIPatch_.cyclicACMIPatch(); // note: only applying coupled contribution const labelUList& nbrFaceCellsCoupled = cpp.neighbPatch().faceCells(); scalarField pnf(psiInternal, nbrFaceCellsCoupled); // Transform according to the transformation tensors transformCoupleField(pnf, cmpt); const labelUList& faceCells = cyclicACMIPatch_.faceCells(); pnf = cyclicACMIPatch_.interpolate(pnf); forAll(faceCells, elemI) { result[faceCells[elemI]] -= coeffs[elemI]*pnf[elemI]; } } template<class Type> void Foam::cyclicACMIFvPatchField<Type>::updateInterfaceMatrix ( Field<Type>& result, const Field<Type>& psiInternal, const scalarField& coeffs, const Pstream::commsTypes ) const { const cyclicACMIPolyPatch& cpp = cyclicACMIPatch_.cyclicACMIPatch(); // note: only applying coupled contribution const labelUList& nbrFaceCellsCoupled = cpp.neighbPatch().faceCells(); Field<Type> pnf(psiInternal, nbrFaceCellsCoupled); // Transform according to the transformation tensors transformCoupleField(pnf); const labelUList& faceCells = cyclicACMIPatch_.faceCells(); pnf = cyclicACMIPatch_.interpolate(pnf); forAll(faceCells, elemI) { result[faceCells[elemI]] -= coeffs[elemI]*pnf[elemI]; } } template<class Type> void Foam::cyclicACMIFvPatchField<Type>::manipulateMatrix ( fvMatrix<Type>& matrix ) { const scalarField& mask = cyclicACMIPatch_.cyclicACMIPatch().mask(); // nothing to be done by the AMI, but re-direct to non-overlap patch // with non-overlap patch weights const fvPatchField<Type>& npf = nonOverlapPatchField(); const_cast<fvPatchField<Type>&>(npf).manipulateMatrix(matrix, 1.0 - mask); } template<class Type> void Foam::cyclicACMIFvPatchField<Type>::updateCoeffs() { // Update non-overlap patch - some will implement updateCoeffs, and // others will implement evaluate // Pass in (1 - mask) to give non-overlap patch the chance to do // manipulation of non-face based data const scalarField& mask = cyclicACMIPatch_.cyclicACMIPatch().mask(); const fvPatchField<Type>& npf = nonOverlapPatchField(); const_cast<fvPatchField<Type>&>(npf).updateWeightedCoeffs(1.0 - mask); } template<class Type> void Foam::cyclicACMIFvPatchField<Type>::write(Ostream& os) const { fvPatchField<Type>::write(os); this->writeEntry("value", os); } // ************************************************************************* // |
|
Above fix adjusts only the patchNeighbour value so will unfortunately make e.g. Gauss gradient incorrect on the partially overlapping bits since the coupled part should only include coupled contributions and the blockage part only the 'wall' contributions. We probably have to have the value field calculated separately, instead of being calculated from internal and neighbour field. |
|
Numerical issues with ACMI-based patch couplings have been resolved by the new Non-Conformal Coupled (NCC) development. See the following commit for explanation, and instructions for how to use NCC and convert cases from AMI. https://github.com/OpenFOAM/OpenFOAM-dev/commit/569fa31d09f98e29d1aaf84d40bb16043f104ec6 |
Date Modified | Username | Field | Change |
---|---|---|---|
2017-05-24 09:52 | niklas.wikstrom | New Issue | |
2017-05-24 09:52 | niklas.wikstrom | File Added: slidingBlockTestCase.tgz | |
2017-05-24 09:52 | niklas.wikstrom | Tag Attached: ACMI | |
2017-05-24 09:52 | niklas.wikstrom | Tag Attached: cyclicACMI | |
2017-05-24 10:04 | niklas.wikstrom | Note Added: 0008185 | |
2017-06-23 16:30 | MattijsJ | Note Added: 0008238 | |
2017-06-23 16:30 | MattijsJ | File Added: cyclicACMIFvPatchField.C | |
2017-06-30 10:42 | MattijsJ | Note Added: 0008294 | |
2022-05-18 11:14 | will | Assigned To | => will |
2022-05-18 11:14 | will | Status | new => resolved |
2022-05-18 11:14 | will | Resolution | open => fixed |
2022-05-18 11:14 | will | Fixed in Version | => dev |
2022-05-18 11:14 | will | Note Added: 0012593 |