View Issue Details

IDProjectCategoryView StatusLast Update
0002511OpenFOAMBugpublic2017-03-22 10:21
Reporterhk318i Assigned Tohenry  
PrioritylowSeveritytextReproducibilityalways
Status resolvedResolutionfixed 
PlatformGNU/LinuxOSUbuntuOS Version14.04
Product Versiondev 
Summary0002511: "Calculate nut" comment
DescriptionWhile I am reading your commit today of applyBoundaryLayer.C, I noticed this comment (Line 133)
    // Calculate nut - reference nut is calculated by the turbulence model
    // on its construction

As I recall this is not valid since OpenFOAM-dev commit 82ccde3269a33a695ca50db326a35a6866d9b6b8
https://github.com/OpenFOAM/OpenFOAM-dev/commit/82ccde3269a33a695ca50db326a35a6866d9b6b8

I'm not sure if nut values used in the calculations are correct or not, unfortunately I do not use this utility and I do not have a test case now.

Best Regards,
Hassan
TagsNo tags attached.

Activities

henry

2017-03-21 17:05

manager   ~0007961

I ran a range of tests of this utility and the results were fine; if you believe there is a problem you will need to provide a test-case which reproduces it.

hk318i

2017-03-22 09:28

reporter   ~0007962

Hello,

I run a quick rough comparison between OF-dev and OF-2.3, there is a slight difference which can be spotted using "writenut" option. Since nut is not calculated in the constructor, its value is based on the initial conditions. The same behaviour as OF-2.3 could be achieved by calling "validate()" function after selecting the model like in solvers.
Generally, it is not a problem for nut because nut will be corrected at the beginning of any solver. But k and epsilon values in this utility are based on nut values. I cannot judge the importance of this minor difference because I do not use this utility. What grabbed my attention is the comment which is not up to date with the code.

Best Regards,
Hassan

henry

2017-03-22 09:45

manager   ~0007963

It is not clear to me what you would like changed. If you could provide a patch that would help.

hk318i

2017-03-22 09:59

reporter  

validate.patch (808 bytes)   
diff --git a/applications/utilities/preProcessing/applyBoundaryLayer/applyBoundaryLayer.C b/applications/utilities/preProcessing/applyBoundaryLayer/applyBoundaryLayer.C
index f2584cb93..a3cf0f33d 100644
--- a/applications/utilities/preProcessing/applyBoundaryLayer/applyBoundaryLayer.C
+++ b/applications/utilities/preProcessing/applyBoundaryLayer/applyBoundaryLayer.C
@@ -131,7 +131,7 @@ int main(int argc, char *argv[])
     if (isA<incompressible::RASModel>(turbulence()))
     {
         // Calculate nut - reference nut is calculated by the turbulence model
-        // on its construction
+        turbulence->validate();
         tmp<volScalarField> tnut = turbulence->nut();
         volScalarField& nut = const_cast<volScalarField&>(tnut());
         volScalarField S(mag(dev(symm(fvc::grad(U)))));
validate.patch (808 bytes)   

hk318i

2017-03-22 10:01

reporter   ~0007964

I attached a patch here, I changed the comment and added a call for validate() to recalculate nut.

henry

2017-03-22 10:21

manager   ~0007966

OK, I see the problem, thanks for the report and patch.

Resolved by commit b3e4c547e8fe39200ff47d28c8b128a9610720c9

Issue History

Date Modified Username Field Change
2017-03-21 16:59 hk318i New Issue
2017-03-21 17:05 henry Note Added: 0007961
2017-03-22 09:28 hk318i Note Added: 0007962
2017-03-22 09:45 henry Note Added: 0007963
2017-03-22 09:59 hk318i File Added: validate.patch
2017-03-22 10:01 hk318i Note Added: 0007964
2017-03-22 10:21 henry Assigned To => henry
2017-03-22 10:21 henry Status new => resolved
2017-03-22 10:21 henry Resolution open => fixed
2017-03-22 10:21 henry Note Added: 0007966