View Issue Details
ID | Project | Category | View Status | Date Submitted | Last Update |
---|---|---|---|---|---|
0002634 | OpenFOAM | Bug | public | 2017-07-26 23:21 | 2017-07-27 16:57 |
Reporter | Shorty | Assigned To | will | ||
Priority | normal | Severity | minor | Reproducibility | always |
Status | resolved | Resolution | fixed | ||
Platform | GNU/Linux | OS | Ubuntu | OS Version | 16.04 |
Product Version | dev | ||||
Summary | 0002634: limitVelocity function | ||||
Description | Hi all, I guess the limitVelocity function has some bug. As far as I understand the function, it limit the magnitude of the velocity field to the one we specify with the max variable, right? However, the way it is implemented does not provide what I expect. I expect the following: - If the magnitude of the velocity of a cell is higher the one we specified, change the velocity in a way to get the magnitude value we specified (max). Right now it is as follow: If the magSqr of the velocity is larger than the sqr of max, we modify the velocity (reduce it) by the ratio of sqr(max)/magSqr(Ui). If we do so, the magnitude of the new vector is much smaller than the one we specified with the max keyword (example below). In addition, the actual implementation reduces the magnitude of the vector even more, the higher its magnitude is. But reducing to the one we specified will never happen. Again: If I got it right, we should do it as follow. If the magnitude of the velocity is larger than the one we specified, reduce the velocity in a way that the magnitude of the new vector is equal to the one we specified with the max keyword, right? If this is correct, then we can use the simple patch I attached. Otherwise I would be happy to understand the actual imeplementation because even with a discussion with my colleague today, we did not get the point. The example below demonstrate the above text: | ||||
Steps To Reproduce | Example given: U = (2 1 4) // of the cell we check out max = 2; mag(U) = 4.58258 --> bigger than max so reduce it. Actual calculation will give us a new velocity vector as: ========================================================= sqrMax = sqr(2) = 4 magSqrUi = magSqr(U) = 21 Unew = U * sqrMax / magSqrUi = 4 / 21 = (0.190476 0.0952381 0.380952) However, the magnitude of the new vector is: mag(Unew) = 0.436436 (I expect it to be 2) Patched version calculates the new velocity vector as: ========================================================= max = 2 magUi = mag(U) = 4.58258 Unew = U * max / magUi = (0.872872 0.436436 1.74574) Now we get the magnitude of the new vector exact to be 2 (as we specified with the max keyword) mag(Unew) = 2 | ||||
Tags | limitVelocity | ||||
|
limitVelocityPatch.C (3,254 bytes)
/*---------------------------------------------------------------------------*\ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | \\ / A nd | Copyright (C) 2016 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 "limitVelocity.H" #include "volFields.H" #include "addToRunTimeSelectionTable.H" // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // namespace Foam { namespace fv { defineTypeNameAndDebug(limitVelocity, 0); addToRunTimeSelectionTable ( option, limitVelocity, dictionary ); } } // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // Foam::fv::limitVelocity::limitVelocity ( const word& name, const word& modelType, const dictionary& dict, const fvMesh& mesh ) : cellSetOption(name, modelType, dict, mesh), UName_(coeffs_.lookupOrDefault<word>("U", "U")), max_(readScalar(coeffs_.lookup("max"))) { fieldNames_.setSize(1, UName_); applied_.setSize(1, false); } // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // bool Foam::fv::limitVelocity::read(const dictionary& dict) { if (cellSetOption::read(dict)) { coeffs_.lookup("max") >> max_; return true; } else { return false; } } void Foam::fv::limitVelocity::correct(volVectorField& U) { const scalar maxU = max_; vectorField& Uif = U.primitiveFieldRef(); forAll(cells_, i) { const label celli = cells_[i]; const scalar magUi = mag(Uif[celli]); if (magUi > maxU) { Uif[celli] *= maxU/magUi; } } // handle boundaries in the case of 'all' if (selectionMode_ == smAll) { volVectorField::Boundary& Ubf = U.boundaryFieldRef(); forAll(Ubf, patchi) { fvPatchVectorField& Up = Ubf[patchi]; if (!Up.fixesValue()) { forAll(Up, facei) { const scalar magUi = mag(Up[facei]); if (magUi > maxU) { Up[facei] *= maxU/magUi; } } } } } } // ************************************************************************* // |
|
Thanks for the report Fixed in dev by commit 16b559c1093019a4a85aa8c810c9d837dfa8fa4c Fixed in 5.x by commit 3ece3e9e81ef05c8336bb7fa5e9a8e0bd78aa6b2 |
|
Just one remark to the patch. Taking the sqrt() of the prefaktor is correct but we could reduce the call of functions as given in my patch :) Well I know, this will not influence the speed of the code too much and is tweak in the micro level :) By the way, you are welcomed Will. Nice to see the patch already there. |
|
Your patch calls mag (which contains sqrt) on every velocity value. My fix only calls sqrt for the velocity values that are larger than the maximum. Sqrt is the expensive bit, so mine will be significantly quicker. |
|
Oh yes you are right :) Thanks Will for the note. Then everything is fine. |
Date Modified | Username | Field | Change |
---|---|---|---|
2017-07-26 23:21 | Shorty | New Issue | |
2017-07-26 23:21 | Shorty | File Added: limitVelocityPatch.C | |
2017-07-26 23:21 | Shorty | Tag Attached: limitVelocity | |
2017-07-27 15:21 | will | Assigned To | => will |
2017-07-27 15:21 | will | Status | new => resolved |
2017-07-27 15:21 | will | Resolution | open => fixed |
2017-07-27 15:21 | will | Note Added: 0008463 | |
2017-07-27 16:04 | Shorty | Status | resolved => feedback |
2017-07-27 16:04 | Shorty | Resolution | fixed => reopened |
2017-07-27 16:04 | Shorty | Note Added: 0008465 | |
2017-07-27 16:11 | will | Status | feedback => resolved |
2017-07-27 16:11 | will | Resolution | reopened => fixed |
2017-07-27 16:11 | will | Note Added: 0008466 | |
2017-07-27 16:54 | Shorty | Status | resolved => feedback |
2017-07-27 16:54 | Shorty | Resolution | fixed => reopened |
2017-07-27 16:54 | Shorty | Note Added: 0008468 | |
2017-07-27 16:57 | will | Status | feedback => resolved |
2017-07-27 16:57 | will | Resolution | reopened => fixed |