From 769dcffe0b4daedef7ddce9b89c7df1f826785fa Mon Sep 17 00:00:00 2001
From: Mark Olesen <Mark.Olesen@esi-group.com>
Date: Tue, 17 Oct 2023 15:00:36 +0200
Subject: [PATCH] ENH: support user-defined sampleFreq for noise models (#2999)

CONFIG: downgrade non-uniform time from error to warning

- can be a spurious error when the deltaT is very small

CONFIG: support keywords 'minFreq', 'maxFreq'

- these are the updated naming for 'fl' and 'fu' (still supported)
---
 .../noise/noiseModels/noiseModel/noiseModel.C | 62 +++++++++++++------
 .../noise/noiseModels/noiseModel/noiseModel.H | 17 +++--
 .../noise/noiseModels/pointNoise/pointNoise.C |  7 ++-
 .../noiseModels/surfaceNoise/surfaceNoise.C   | 18 ++++--
 4 files changed, 72 insertions(+), 32 deletions(-)

diff --git a/src/randomProcesses/noise/noiseModels/noiseModel/noiseModel.C b/src/randomProcesses/noise/noiseModels/noiseModel/noiseModel.C
index f72297c6e52..73cca53db73 100644
--- a/src/randomProcesses/noise/noiseModels/noiseModel/noiseModel.C
+++ b/src/randomProcesses/noise/noiseModels/noiseModel/noiseModel.C
@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2015-2022 OpenCFD Ltd.
+    Copyright (C) 2015-2023 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -185,29 +185,34 @@ Foam::scalar Foam::noiseModel::checkUniformTimeStep
     const scalarList& times
 ) const
 {
-    scalar deltaT = -1.0;
+    scalar deltaT = 1;
 
     if (times.size() > 1)
     {
         // Uniform time step
-        deltaT = (times.last() - times.first())/scalar(times.size() - 1);
+        deltaT = (times.back() - times.front())/scalar(times.size() - 1);
 
-        for (label i = 1; i < times.size(); ++i)
+        bool nonUniform = false;
+        for (label i = 1; i < times.size() && !nonUniform; ++i)
         {
-            scalar dT = times[i] - times[i-1];
+            const scalar dT = times[i] - times[i-1];
 
-            if (mag(dT/deltaT - 1) > 1e-8)
-            {
-                FatalErrorInFunction
-                    << "Unable to process data with a variable time step"
-                    << exit(FatalError);
-            }
+            nonUniform = (mag(dT/deltaT - 1) > 1e-8);
+        }
+
+        if (nonUniform)
+        {
+            WarningInFunction
+                << "Detected non-uniform time step:"
+                << " resulting FFT may be incorrect"
+                << " (or the deltaT is just very small): " << deltaT
+                << endl;
         }
     }
     else
     {
         FatalErrorInFunction
-            << "Unable to create FFT with a single value"
+            << "Unable to create FFT with 0 or 1 values"
             << exit(FatalError);
     }
 
@@ -617,6 +622,7 @@ Foam::noiseModel::noiseModel
     nSamples_(65536),
     fLower_(25),
     fUpper_(10000),
+    sampleFreq_(0),
     startTime_(0),
     windowModelPtr_(),
     SPLweighting_(weightingType::none),
@@ -654,10 +660,16 @@ bool Foam::noiseModel::read(const dictionary& dict)
         return false;
     }
 
+    // Re-assign defaults (to avoid stickiness)
+    fLower_ = 25;
+    fUpper_ = 10000;
+    sampleFreq_ = 0;
+
     dict.readIfPresent("rhoRef", rhoRef_);
     dict.readIfPresent("N", nSamples_);
-    dict.readIfPresent("fl", fLower_);
-    dict.readIfPresent("fu", fUpper_);
+    dict.readIfPresentCompat("minFreq", {{"fl", 2312}}, fLower_);
+    dict.readIfPresentCompat("maxFreq", {{"fu", 2312}}, fUpper_);
+    dict.readIfPresent("sampleFreq", sampleFreq_);
     dict.readIfPresent("startTime", startTime_);
     dict.readIfPresent("minPressure", minPressure_);
     dict.readIfPresent("maxPressure", maxPressure_);
@@ -666,29 +678,39 @@ bool Foam::noiseModel::read(const dictionary& dict)
     if (fLower_ < 0)
     {
         FatalIOErrorInFunction(dict)
-            << "fl: lower frequency bound must be greater than zero"
+            << "Lower frequency bound must be greater than zero"
             << exit(FatalIOError);
     }
 
     if (fUpper_ < 0)
     {
         FatalIOErrorInFunction(dict)
-            << "fu: upper frequency bound must be greater than zero"
+            << "Upper frequency bound must be greater than zero"
             << exit(FatalIOError);
     }
 
     if (fUpper_ < fLower_)
     {
         FatalIOErrorInFunction(dict)
-            << "fu: upper frequency bound must be greater than lower "
-            << "frequency bound (fl)"
+            << "Upper frequency bound (" << fUpper_
+            << ") must be greater than lower frequency bound ("
+            << fLower_ << ")"
             << exit(FatalIOError);
-
     }
 
     Info<< "    Frequency bounds:" << nl
         << "        lower: " << fLower_ << nl
-        << "        upper: " << fUpper_ << endl;
+        << "        upper: " << fUpper_ << nl
+        << "       sample: ";
+
+    if (sampleFreq_ > 0)
+    {
+        Info<< sampleFreq_ << nl;
+    }
+    else
+    {
+        Info<< "auto" << nl;
+    }
 
     weightingTypeNames_.readIfPresent("SPLweighting", dict, SPLweighting_);
 
diff --git a/src/randomProcesses/noise/noiseModels/noiseModel/noiseModel.H b/src/randomProcesses/noise/noiseModels/noiseModel/noiseModel.H
index c3f1267e434..14c8d6db443 100644
--- a/src/randomProcesses/noise/noiseModels/noiseModel/noiseModel.H
+++ b/src/randomProcesses/noise/noiseModels/noiseModel/noiseModel.H
@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2015-2022 OpenCFD Ltd.
+    Copyright (C) 2015-2023 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -34,8 +34,8 @@ Description
     \verbatim
     rhoRef          1;
     N               4096;
-    fl              25;
-    fu              10000;
+    minFreq         25;
+    maxFreq         10000;
     startTime       0;
 
     outputPrefix    "test1";
@@ -58,8 +58,9 @@ Description
         Property     | Description                   | Required  | Default value
         rhoRef       | Reference density             | no        | 1
         N            | Number of samples in sampling window | no | 65536 (2^16)
-        fl           | Lower frequency bounds        | no        | 25
-        fu           | Upper frequency bounds        | no        | 10000
+        minFreq      | Lower frequency bounds (fl)   | no        | 25
+        maxFreq      | Upper frequency bounds (fu)   | no        | 10000
+        sampleFreq   | Sample frequency              | no        | 0
         startTime    | Start time                    | no        | 0
         outputPrefix | Prefix applied to output files| no        | ''
         SPLweighting | Weighting: dBA, dBB, dBC, DBD | no        | none
@@ -71,6 +72,9 @@ Description
         writeOctaves | Write octaves data            | no        | yes
     \endtable
 
+    If provided, the sampleFreq is used to define the deltaT between samples.
+    Otherwise the deltaT is determined from the time samples themselves.
+
 SourceFiles
     noiseModel.C
 
@@ -158,6 +162,9 @@ protected:
         //- Upper frequency limit, default = 10kHz
         scalar fUpper_;
 
+        //- Prescribed sample frequency
+        scalar sampleFreq_;
+
         //- Start time, default = 0s
         scalar startTime_;
 
diff --git a/src/randomProcesses/noise/noiseModels/pointNoise/pointNoise.C b/src/randomProcesses/noise/noiseModels/pointNoise/pointNoise.C
index 88b741efb9a..462abc3ba9f 100644
--- a/src/randomProcesses/noise/noiseModels/pointNoise/pointNoise.C
+++ b/src/randomProcesses/noise/noiseModels/pointNoise/pointNoise.C
@@ -101,7 +101,12 @@ void pointNoise::processData
 
     Info<< "Creating noise FFT" << endl;
 
-    const scalar deltaT = checkUniformTimeStep(t);
+    const scalar deltaT =
+    (
+        sampleFreq_ > 0
+      ? (1.0/sampleFreq_)
+      : checkUniformTimeStep(t)
+    );
 
     // Apply conversions
     p *= rhoRef_;
diff --git a/src/randomProcesses/noise/noiseModels/surfaceNoise/surfaceNoise.C b/src/randomProcesses/noise/noiseModels/surfaceNoise/surfaceNoise.C
index 806b815d817..fa4473ee4a6 100644
--- a/src/randomProcesses/noise/noiseModels/surfaceNoise/surfaceNoise.C
+++ b/src/randomProcesses/noise/noiseModels/surfaceNoise/surfaceNoise.C
@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2015-2022 OpenCFD Ltd.
+    Copyright (C) 2015-2023 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -96,15 +96,21 @@ void surfaceNoise::initialise(const fileName& fName)
     // Note: all processors should call the windowing validate function
     label nRequiredTimes = windowModelPtr_->validate(nAvailableTimes);
 
-    if (Pstream::master())
+    if (UPstream::master())
     {
         // Restrict times
-        times_.setSize(nRequiredTimes);
-        forAll(times_, timeI)
+        times_.resize_nocopy(nRequiredTimes);
+        forAll(times_, timei)
         {
-            times_[timeI] = allTimes[timeI + startTimeIndex_].value();
+            times_[timei] = allTimes[timei + startTimeIndex_].value();
         }
-        deltaT_ = checkUniformTimeStep(times_);
+
+        deltaT_ =
+        (
+            sampleFreq_ > 0
+          ? (1.0/sampleFreq_)
+          : checkUniformTimeStep(times_)
+        );
 
         // Read the surface geometry
         // Note: hard-coded to read mesh from first time index
-- 
GitLab