diff --git a/src/OpenFOAM/graph/graph.C b/src/OpenFOAM/graph/graph.C index a7698a20f450fb4190188e5c7fda59475fb0a96d..5eb1599d6e4658a6c76a4068dd6fa7632ca9b86f 100644 --- a/src/OpenFOAM/graph/graph.C +++ b/src/OpenFOAM/graph/graph.C @@ -28,6 +28,7 @@ License #include "IOmanip.H" #include "Pair.H" #include "OSspecific.H" +#include "SubField.H" // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // @@ -160,6 +161,45 @@ Foam::scalarField& Foam::graph::y() } +void Foam::graph::setXRange(const scalar x0, const scalar x1) +{ + if (x1 < x0) + { + FatalErrorInFunction + << "When setting limits, x1 must be greater than x0" << nl + << " x0: " << x0 << nl + << " x1: " << x1 << nl + << abort(FatalError); + } + + label i0 = 0; + label i1 = 0; + + forAll(x_, i) + { + if (x_[i] < x0) + { + i0 = i + 1; + } + if (x_[i] < x1) + { + i1 = i; + } + } + + label nX = i1 - i0 + 1; + scalarField xNew(SubField<scalar>(x_, nX, i0)); + x_.transfer(xNew); + + forAllIter(HashPtrTable<curve>, *this, iter) + { + curve* c = iter(); + scalarField cNew(SubField<scalar>(*c, nX, i0)); + c->transfer(cNew); + } +} + + Foam::autoPtr<Foam::graph::writer> Foam::graph::writer::New ( const word& graphFormat diff --git a/src/OpenFOAM/graph/graph.H b/src/OpenFOAM/graph/graph.H index b3ff00c37d5938a826bc3cc889d88622fc78ba8f..d77d468d34cec2c52d0aa1b3e61e0f454c55e2c1 100644 --- a/src/OpenFOAM/graph/graph.H +++ b/src/OpenFOAM/graph/graph.H @@ -178,6 +178,9 @@ public: scalarField& y(); + // Limit the data for all axes based on x-limits + void setXRange(const scalar x0, const scalar x1); + // Write diff --git a/src/functionObjects/forces/forces/forces.C b/src/functionObjects/forces/forces/forces.C index 97de38581f28c37ba7fa76e93f8996b639b82cf2..669fa94d8506a6015b6c9e0761acc84f2ba5f880 100644 --- a/src/functionObjects/forces/forces/forces.C +++ b/src/functionObjects/forces/forces/forces.C @@ -353,7 +353,7 @@ Foam::functionObjects::forces::devRhoReff() const else if (foundObject<dictionary>("transportProperties")) { const dictionary& transportProperties = - lookupObject<dictionary>("transportProperties"); + lookupObject<dictionary>("transportProperties"); dimensionedScalar nu(transportProperties.lookup("nu")); diff --git a/src/meshTools/coordinateSystems/cylindricalCS.H b/src/meshTools/coordinateSystems/cylindricalCS.H index 16b46a72456c8f19d15d66483ef7c030a1edcff7..4504231caeec791f1ce012f02afbd1ce1372f667 100644 --- a/src/meshTools/coordinateSystems/cylindricalCS.H +++ b/src/meshTools/coordinateSystems/cylindricalCS.H @@ -142,7 +142,7 @@ public: // Member Functions //- Are angles in degrees? - bool inDegrees() const; + bool inDegrees() const; //- Non-const access to inDegrees bool& inDegrees(); diff --git a/src/randomProcesses/noise/noiseModels/noiseModel/noiseModel.C b/src/randomProcesses/noise/noiseModels/noiseModel/noiseModel.C index ae80887fe89f8d15fff01d275f050dd86d837e3b..b98f00c8bda13ec2afc036a9bb0e2ae7c9a24849 100644 --- a/src/randomProcesses/noise/noiseModels/noiseModel/noiseModel.C +++ b/src/randomProcesses/noise/noiseModels/noiseModel/noiseModel.C @@ -35,6 +35,26 @@ namespace Foam // * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * // +void Foam::noiseModel::readWriteOption +( + const dictionary& dict, + const word& lookup, + bool& option +) const +{ + dict.readIfPresent(lookup, option); + + if (option) + { + Info<< " " << lookup << ": " << "yes" << endl; + } + else + { + Info<< " " << lookup << ": " << "no" << endl; + } +} + + Foam::scalar Foam::noiseModel::checkUniformTimeStep ( const scalarList& times @@ -92,6 +112,22 @@ Foam::label Foam::noiseModel::findStartTimeIndex } +Foam::fileName Foam::noiseModel::baseFileDir(const label dataseti) const +{ + fileName baseDir("$FOAM_CASE"); + word datasetName("input" + Foam::name(dataseti)); + baseDir = + baseDir.expand() + /"postProcessing" + /"noise" + /outputPrefix_ + /type() + /datasetName; + + return baseDir; +} + + // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // Foam::noiseModel::noiseModel(const dictionary& dict, const bool readFields) @@ -101,9 +137,16 @@ Foam::noiseModel::noiseModel(const dictionary& dict, const bool readFields) nSamples_(65536), fLower_(25), fUpper_(10000), + customBounds_(false), startTime_(0), windowModelPtr_(), - graphFormat_("raw") + graphFormat_("raw"), + outputPrefix_(), + writePrmsf_(true), + writeSPL_(true), + writePSD_(true), + writePSDf_(true), + writeOctaves_(true) { if (readFields) { @@ -124,10 +167,18 @@ bool Foam::noiseModel::read(const dictionary& dict) { dict.readIfPresent("rhoRef", rhoRef_); dict.readIfPresent("N", nSamples_); - dict.readIfPresent("fl", fLower_); - dict.readIfPresent("fu", fUpper_); + customBounds_ = false; + if (dict.readIfPresent("fl", fLower_)) + { + customBounds_ = true; + } + if (dict.readIfPresent("fu", fUpper_)) + { + customBounds_ = true; + } dict.readIfPresent("startTime", startTime_); dict.readIfPresent("graphFormat", graphFormat_); + dict.readIfPresent("outputPrefix", outputPrefix_); // Check number of samples - must be a power of 2 for our FFT bool powerOf2 = ((nSamples_ != 0) && !(nSamples_ & (nSamples_ - 1))); @@ -164,6 +215,14 @@ bool Foam::noiseModel::read(const dictionary& dict) } + Info<< " Write options:" << endl; + dictionary optDict(dict.subOrEmptyDict("writeOptions")); + readWriteOption(optDict, "writePrmsf", writePrmsf_); + readWriteOption(optDict, "writeSPL", writeSPL_); + readWriteOption(optDict, "writePSD", writePSD_); + readWriteOption(optDict, "writeOctaves", writeOctaves_); + + windowModelPtr_ = windowModel::New(dict, nSamples_); return true; diff --git a/src/randomProcesses/noise/noiseModels/noiseModel/noiseModel.H b/src/randomProcesses/noise/noiseModels/noiseModel/noiseModel.H index 89c41edb84c927554f505235f19fe83668d1ced5..7b625096cea8b742e2a28b5b6d5539c60ea0e9a7 100644 --- a/src/randomProcesses/noise/noiseModels/noiseModel/noiseModel.H +++ b/src/randomProcesses/noise/noiseModels/noiseModel/noiseModel.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2015-2016 OpenCFD Ltd. + \\ / A nd | Copyright (C) 2015-2017 OpenCFD Ltd. \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -35,17 +35,35 @@ Description fl 25; fu 25; startTime 0; + + outputPrefix "test1"; + + // Optional write options dictionary + writeOptions + { + writePrmsf no; + writeSPL yes; + writePSD yes; + writePSDf no; + writeOctaves yes; + } \endverbatim where \table - 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 - startTime | Start time | no | 0 - graphFormat | Graph format | no | raw + 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 + startTime | Start time | no | 0 + outputPrefix | Prefix applied to output files| no | '' + graphFormat | Graph format | no | raw + writePrmsf | Write Prmsf data | no | yes + writeSPL | Write SPL data | no | yes + writePSD | Write PSD data | no | yes + writePSDf | Write PSDf data | no | yes + writeOctaves | Write octaves data | no | yes \endtable Note @@ -108,6 +126,9 @@ protected: //- Upper frequency limit, default = 10kHz scalar fUpper_; + //- Flag to indicate that custom frequenct bounds are being used + bool customBounds_; + //- Start time, default = 0s scalar startTime_; @@ -118,8 +139,37 @@ protected: word graphFormat_; + // Write options + + //- Output file prefix, default = '' + fileName outputPrefix_; + + //- Write Prmsf; default = yes + bool writePrmsf_; + + //- Write SPL; default = yes + bool writeSPL_; + + //- Write PSD; default = yes + bool writePSD_; + + //- Write PSDf; default = yes + bool writePSDf_; + + //- Write writeOctaves; default = yes + bool writeOctaves_; + + // Protected Member Functions + //- Helper function to read write options and provide info feedback + void readWriteOption + ( + const dictionary& dict, + const word& lookup, + bool& option + ) const; + //- Check and return uniform time step scalar checkUniformTimeStep ( @@ -133,6 +183,9 @@ protected: const scalar startTime ) const; + //- Return the base output directory + fileName baseFileDir(const label dataseti) const; + public: diff --git a/src/randomProcesses/noise/noiseModels/pointNoise/pointNoise.C b/src/randomProcesses/noise/noiseModels/pointNoise/pointNoise.C index 88b245b62c4773cef598789803d9f1265750297c..083bcdf31bf54ae2ab153c46440d07312945998e 100644 --- a/src/randomProcesses/noise/noiseModels/pointNoise/pointNoise.C +++ b/src/randomProcesses/noise/noiseModels/pointNoise/pointNoise.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2015-2016 OpenCFD Ltd. + \\ / A nd | Copyright (C) 2015-2017 OpenCFD Ltd. \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -66,7 +66,11 @@ void pointNoise::filterTimeData } -void pointNoise::processData(const Function1Types::CSV<scalar>& data) +void pointNoise::processData +( + const label dataseti, + const Function1Types::CSV<scalar>& data +) { Info<< "Reading data file " << data.fName() << endl; @@ -86,7 +90,7 @@ void pointNoise::processData(const Function1Types::CSV<scalar>& data) windowModelPtr_->validate(t.size()); const windowModel& win = windowModelPtr_(); const scalar deltaf = 1.0/(deltaT*win.nSamples()); - fileName outDir(fileName("postProcessing")/"noise"/typeName/fNameBase); + fileName outDir(baseFileDir(dataseti)/fNameBase); // Create the fft noiseFFT nfft(deltaT, p); @@ -97,13 +101,27 @@ void pointNoise::processData(const Function1Types::CSV<scalar>& data) // RMS pressure [Pa] graph Prmsf(nfft.RMSmeanPf(win)); - Info<< " Creating graph for " << Prmsf.title() << endl; - Prmsf.write(outDir, graph::wordify(Prmsf.title()), graphFormat_); + if (customBounds_) + { + Prmsf.setXRange(fLower_, fUpper_); + } + if (writePrmsf_) + { + Info<< " Creating graph for " << Prmsf.title() << endl; + Prmsf.write(outDir, graph::wordify(Prmsf.title()), graphFormat_); + } // PSD [Pa^2/Hz] graph PSDf(nfft.PSDf(win)); - Info<< " Creating graph for " << PSDf.title() << endl; - PSDf.write(outDir, graph::wordify(PSDf.title()), graphFormat_); + if (customBounds_) + { + PSDf.setXRange(fLower_, fUpper_); + } + if (writePSDf_) + { + Info<< " Creating graph for " << PSDf.title() << endl; + PSDf.write(outDir, graph::wordify(PSDf.title()), graphFormat_); + } // PSD [dB/Hz] graph PSDg @@ -114,8 +132,12 @@ void pointNoise::processData(const Function1Types::CSV<scalar>& data) Prmsf.x(), noiseFFT::PSD(PSDf.y()) ); - Info<< " Creating graph for " << PSDg.title() << endl; - PSDg.write(outDir, graph::wordify(PSDg.title()), graphFormat_); + + if (writePSD_) + { + Info<< " Creating graph for " << PSDg.title() << endl; + PSDg.write(outDir, graph::wordify(PSDg.title()), graphFormat_); + } // SPL [dB] graph SPLg @@ -126,79 +148,58 @@ void pointNoise::processData(const Function1Types::CSV<scalar>& data) Prmsf.x(), noiseFFT::SPL(PSDf.y()*deltaf) ); - Info<< " Creating graph for " << SPLg.title() << endl; - SPLg.write(outDir, graph::wordify(SPLg.title()), graphFormat_); - - labelList octave13BandIDs; - scalarField octave13FreqCentre; - noiseFFT::octaveBandInfo - ( - Prmsf.x(), - fLower_, - fUpper_, - 3, - octave13BandIDs, - octave13FreqCentre - ); - - // 1/3 octave data - // --------------- - - // PSD [Pa^2/Hz] - graph PSD13f(nfft.octaves(PSDf, octave13BandIDs, false)); - - // Integrated PSD = P(rms)^2 [Pa^2] - graph Prms13f2(nfft.octaves(PSDf, octave13BandIDs, true)); - - graph PSD13g - ( - "PSD13_dB_Hz(fm)", - "fm [Hz]", - "PSD(fm) [dB_Hz]", - octave13FreqCentre, - noiseFFT::PSD(PSD13f.y()) - ); - Info<< " Creating graph for " << PSD13g.title() << endl; - PSD13g.write(outDir, graph::wordify(PSD13g.title()), graphFormat_); - - graph SPL13g - ( - "SPL13_dB(fm)", - "fm [Hz]", - "SPL(fm) [dB]", - octave13FreqCentre, - noiseFFT::SPL(Prms13f2.y()) - ); - Info<< " Creating graph for " << SPL13g.title() << endl; - SPL13g.write(outDir, graph::wordify(SPL13g.title()), graphFormat_); -} - - -// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // - -void pointNoise::calculate() -{ - // Point data only handled by master - if (!Pstream::master()) + if (writeSPL_) { - return; + Info<< " Creating graph for " << SPLg.title() << endl; + SPLg.write(outDir, graph::wordify(SPLg.title()), graphFormat_); } - - if (inputFileNames_.size()) - { - forAll(inputFileNames_, i) - { - const fileName fName = inputFileNames_[i].expand(); - Function1Types::CSV<scalar> data("pressure", dict_, "Data", fName); - processData(data); - } - } - else + if (writeOctaves_) { - Function1Types::CSV<scalar> data("pressure", dict_, "Data"); - processData(data); + labelList octave13BandIDs; + scalarField octave13FreqCentre; + noiseFFT::octaveBandInfo + ( + Prmsf.x(), + fLower_, + fUpper_, + 3, + octave13BandIDs, + octave13FreqCentre + ); + + + // 1/3 octave data + // --------------- + + // PSD [Pa^2/Hz] + graph PSD13f(nfft.octaves(PSDf, octave13BandIDs, false)); + + // Integrated PSD = P(rms)^2 [Pa^2] + graph Prms13f2(nfft.octaves(PSDf, octave13BandIDs, true)); + + graph PSD13g + ( + "PSD13_dB_Hz(fm)", + "fm [Hz]", + "PSD(fm) [dB_Hz]", + octave13FreqCentre, + noiseFFT::PSD(PSD13f.y()) + ); + Info<< " Creating graph for " << PSD13g.title() << endl; + PSD13g.write(outDir, graph::wordify(PSD13g.title()), graphFormat_); + + graph SPL13g + ( + "SPL13_dB(fm)", + "fm [Hz]", + "SPL(fm) [dB]", + octave13FreqCentre, + noiseFFT::SPL(Prms13f2.y()) + ); + Info<< " Creating graph for " << SPL13g.title() << endl; + SPL13g.write(outDir, graph::wordify(SPL13g.title()), graphFormat_); } } @@ -222,11 +223,45 @@ pointNoise::~pointNoise() {} +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +void pointNoise::calculate() +{ + // Point data only handled by master + if (!Pstream::master()) + { + return; + } + + + forAll(inputFileNames_, filei) + { + fileName fName = inputFileNames_[filei]; + fName.expand(); + if (!fName.isAbsolute()) + { + fName = "$FOAM_CASE"/fName; + } + fName.expand(); + Function1Types::CSV<scalar> data("pressure", dict_, "Data", fName); + processData(filei, data); + } +} + + bool pointNoise::read(const dictionary& dict) { if (noiseModel::read(dict)) { - dict.readIfPresent("inputFiles", inputFileNames_); + if (!dict.readIfPresent("inputFiles", inputFileNames_)) + { + inputFileNames_.setSize(1); + + // Note: unsafe lookup of file name from pressureData sub-dict! + dict.subDict("pressureData").lookup("fileName") + >> inputFileNames_[0]; + } + return true; } diff --git a/src/randomProcesses/noise/noiseModels/pointNoise/pointNoise.H b/src/randomProcesses/noise/noiseModels/pointNoise/pointNoise.H index 374678487f5089f83ec31bcdd7009a5372cf77de..693418939a8bacea33782a1726705dc9aa5504df 100644 --- a/src/randomProcesses/noise/noiseModels/pointNoise/pointNoise.H +++ b/src/randomProcesses/noise/noiseModels/pointNoise/pointNoise.H @@ -117,7 +117,11 @@ protected: ) const; //- Process the CSV data - void processData(const Function1Types::CSV<scalar>& data); + void processData + ( + const label dataseti, + const Function1Types::CSV<scalar>& data + ); public: diff --git a/src/randomProcesses/noise/noiseModels/surfaceNoise/surfaceNoise.C b/src/randomProcesses/noise/noiseModels/surfaceNoise/surfaceNoise.C index 051f016969b456a30aacf80c15e9b4e436ad4206..7b9e6180b95fd27d1b2f069a75e3c9f330b8fab7 100644 --- a/src/randomProcesses/noise/noiseModels/surfaceNoise/surfaceNoise.C +++ b/src/randomProcesses/noise/noiseModels/surfaceNoise/surfaceNoise.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2015-2016 OpenCFD Ltd. + \\ / A nd | Copyright (C) 2015-2017 OpenCFD Ltd. \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -226,20 +226,18 @@ void surfaceNoise::readSurfaceData Foam::scalar surfaceNoise::writeSurfaceData ( + const fileName& outDirBase, const word& fName, - const word& groupName, const word& title, const scalar freq, const scalarField& data, - const labelList& procFaceOffset + const labelList& procFaceOffset, + const bool writeSurface ) const { Info<< " processing " << title << " for frequency " << freq << endl; - fileName outDir - ( - fileName("postProcessing")/"noise"/groupName/Foam::name(freq) - ); + const fileName outDir(outDirBase/Foam::name(freq)); if (Pstream::parRun()) { @@ -279,20 +277,23 @@ Foam::scalar surfaceNoise::writeSurfaceData } } - // could also have meshedSurface implement meshedSurf - fileName outFileName = writerPtr_->write - ( - outDir, - fName, - meshedSurfRef + // Could also have meshedSurface implement meshedSurf + if (writeSurface) + { + fileName outFileName = writerPtr_->write ( - surf.points(), - surf.surfFaces() - ), - title, - allData, - false - ); + outDir, + fName, + meshedSurfRef + ( + surf.points(), + surf.surfFaces() + ), + title, + allData, + false + ); + } // TO BE VERIFIED: area-averaged values // areaAverage = sum(allData*surf.magSf())/sum(surf.magSf()); @@ -306,20 +307,23 @@ Foam::scalar surfaceNoise::writeSurfaceData { const meshedSurface& surf = readerPtr_->geometry(); - // could also have meshedSurface implement meshedSurf - writerPtr_->write - ( - outDir, - fName, - meshedSurfRef + // Could also have meshedSurface implement meshedSurf + if (writeSurface) + { + writerPtr_->write ( - surf.points(), - surf.surfFaces() - ), - title, - data, - false - ); + outDir, + fName, + meshedSurfRef + ( + surf.points(), + surf.surfFaces() + ), + title, + data, + false + ); + } // TO BE VERIFIED: area-averaged values // return sum(data*surf.magSf())/sum(surf.magSf()); @@ -459,9 +463,15 @@ bool surfaceNoise::read(const dictionary& dict) void surfaceNoise::calculate() { - forAll(inputFileNames_, i) + forAll(inputFileNames_, filei) { - fileName fName = inputFileNames_[i]; + fileName fName = inputFileNames_[filei]; + fName.expand(); + + if (!fName.isAbsolute()) + { + fName = "$FOAM_CASE"/fName; + } initialise(fName.expand()); @@ -494,10 +504,16 @@ void surfaceNoise::calculate() Info<< "Creating noise FFTs" << endl; + const scalarField freq1(noiseFFT::frequencies(nSamples_, deltaT_)); + + // Reset desired frequency range if outside actual frequency range + fLower_ = min(fLower_, max(freq1)); + fUpper_ = min(fUpper_, max(freq1)); + // Storage for FFT data const label nLocalFace = pData.size(); - const scalarField freq1(noiseFFT::frequencies(nSamples_, deltaT_)); - const label nFFT = freq1.size()/fftWriteInterval_; + const label nFFT = ceil(freq1.size()/scalar(fftWriteInterval_)); + List<scalarField> surfPrmsf(nFFT); List<scalarField> surfPSDf(nFFT); forAll(surfPrmsf, freqI) @@ -553,7 +569,7 @@ void surfaceNoise::calculate() // Store the frequency results in slot for face of surface forAll(surfPrmsf, i) { - label freqI = (i + 1)*fftWriteInterval_ - 1; + label freqI = i*fftWriteInterval_; surfPrmsf[i][faceI] = Prmsf.y()[freqI]; surfPSDf[i][faceI] = PSDf.y()[freqI]; } @@ -575,60 +591,90 @@ void surfaceNoise::calculate() const word fNameBase = fName.nameLessExt(); // Output directory for graphs - fileName outDir - ( - fileName("postProcessing")/"noise"/typeName/fNameBase - ); + fileName outDirBase(baseFileDir(filei)/fNameBase); const scalar deltaf = 1.0/(deltaT_*win.nSamples()); - Info<< "Writing fft surface data" << endl; + Info<< "Writing fft surface data"; + if (fftWriteInterval_ == 1) { - scalarField PrmsfAve(surfPrmsf.size(), 0); - scalarField PSDfAve(surfPrmsf.size(), 0); - scalarField fOut(surfPrmsf.size(), 0); + Info<< endl; + } + else + { + Info<< " at every " << fftWriteInterval_ << " frequency points" + << endl; + } - forAll(surfPrmsf, i) - { - label freqI = (i + 1)*fftWriteInterval_ - 1; - fOut[i] = freq1[freqI]; - const word gName = "fft"; - PrmsfAve[i] = writeSurfaceData - ( - fNameBase, - gName, - "Prmsf", - freq1[freqI], - surfPrmsf[i], - procFaceOffset - ); + { + fileName outDir(outDirBase/"fft"); - PSDfAve[i] = writeSurfaceData - ( - fNameBase, - gName, - "PSDf", - freq1[freqI], - surfPSDf[i], - procFaceOffset - ); - writeSurfaceData - ( - fNameBase, - gName, - "PSD", - freq1[freqI], - noiseFFT::PSD(surfPSDf[i]), - procFaceOffset - ); - writeSurfaceData - ( - fNameBase, - gName, - "SPL", - freq1[freqI], - noiseFFT::SPL(surfPSDf[i]*deltaf), - procFaceOffset - ); + // Determine frequency range of interest + // Note: freqencies have fixed interval, and are in the range + // 0 to fftWriteInterval_*(n-1)*deltaf + label f0 = ceil(fLower_/deltaf/scalar(fftWriteInterval_)); + label f1 = floor(fUpper_/deltaf/scalar(fftWriteInterval_)); + label nFreq = f1 - f0; + + scalarField PrmsfAve(nFreq, 0); + scalarField PSDfAve(nFreq, 0); + scalarField fOut(nFreq, 0); + + if (nFreq == 0) + { + WarningInFunction + << "No surface data available using a fftWriteInterval of " + << fftWriteInterval_ << endl; + } + else + { + forAll(fOut, i) + { + label freqI = (i + f0)*fftWriteInterval_; + fOut[i] = freq1[freqI]; + + + PrmsfAve[i] = writeSurfaceData + ( + outDir, + fNameBase, + "Prmsf", + freq1[freqI], + surfPrmsf[i + f0], + procFaceOffset, + writePrmsf_ + ); + + PSDfAve[i] = writeSurfaceData + ( + outDir, + fNameBase, + "PSDf", + freq1[freqI], + surfPSDf[i + f0], + procFaceOffset, + writePSDf_ + ); + writeSurfaceData + ( + outDir, + fNameBase, + "PSD", + freq1[freqI], + noiseFFT::PSD(surfPSDf[i + f0]), + procFaceOffset, + writePSD_ + ); + writeSurfaceData + ( + outDir, + fNameBase, + "SPL", + freq1[freqI], + noiseFFT::SPL(surfPSDf[i + f0]*deltaf), + procFaceOffset, + writeSPL_ + ); + } } graph Prmsfg @@ -675,38 +721,42 @@ void surfaceNoise::calculate() Info<< "Writing one-third octave surface data" << endl; { + fileName outDir(outDirBase/"oneThirdOctave"); + scalarField PSDfAve(surfPSD13f.size(), 0); scalarField Prms13f2Ave(surfPSD13f.size(), 0); forAll(surfPSD13f, i) { - const word gName = "oneThirdOctave"; PSDfAve[i] = writeSurfaceData ( + outDir, fNameBase, - gName, "PSD13f", octave13FreqCentre[i], surfPSD13f[i], - procFaceOffset + procFaceOffset, + writeOctaves_ ); writeSurfaceData ( + outDir, fNameBase, - gName, "PSD13", octave13FreqCentre[i], noiseFFT::PSD(surfPSD13f[i]), - procFaceOffset + procFaceOffset, + writeOctaves_ ); writeSurfaceData ( + outDir, fNameBase, - gName, "SPL13", octave13FreqCentre[i], noiseFFT::SPL(surfPrms13f2[i]), - procFaceOffset + procFaceOffset, + writeOctaves_ ); Prms13f2Ave[i] = diff --git a/src/randomProcesses/noise/noiseModels/surfaceNoise/surfaceNoise.H b/src/randomProcesses/noise/noiseModels/surfaceNoise/surfaceNoise.H index 165a1cefd77963b882d10c4d44285eaa79c6258e..6450881d76f4006fc300983312d4de123126552c 100644 --- a/src/randomProcesses/noise/noiseModels/surfaceNoise/surfaceNoise.H +++ b/src/randomProcesses/noise/noiseModels/surfaceNoise/surfaceNoise.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2015-2016 OpenCFD Ltd. + \\ / A nd | Copyright (C) 2015-2017 OpenCFD Ltd. \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -71,7 +71,23 @@ Description { collateTimes true; } + + // Write Prmsf; default = yes + writePrmsf no; + + // Write SPL; default = yes + writeSPL yes; + + // Write PSD; default = yes + writePSD yes; + + // Write PSDf; default = yes + writePSDf no; + + // Write writeOctaves; default = yes + writeOctaves yes; } + \endverbatim SourceFiles @@ -167,12 +183,13 @@ protected: // Returns the area average value scalar writeSurfaceData ( + const fileName& outDirBase, const word& fName, - const word& groupName, const word& title, const scalar freq, const scalarField& data, - const labelList& procFaceOffset + const labelList& procFaceOffset, + const bool writeSurface ) const; //- Calculate the area average value