src/lagrangian/molecularDynamics/molecularMeasurements/correlationFunction/correlationFunction.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 2008-2009 OpenCFD Ltd.
   6      \\/     M anipulation  |
   7 -------------------------------------------------------------------------------
   8 License
   9     This file is part of OpenFOAM.
  10
  11     OpenFOAM is free software; you can redistribute it and/or modify it
  12     under the terms of the GNU General Public License as published by the
  13     Free Software Foundation; either version 2 of the License, or (at your
  14     option) any later version.
  15
  16     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
  17     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
  18     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
  19     for more details.
  20
  21     You should have received a copy of the GNU General Public License
  22     along with OpenFOAM; if not, write to the Free Software Foundation,
  23     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
  24
  25 \*---------------------------------------------------------------------------*/
  26
  27 #include "correlationFunction.H"
  28
  29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
  30
  31 template<class Type>
  32 const char* const
  33     Foam::correlationFunction<Type>::typeName("correlationFunction");
  34
  35
  36 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
  37
  38 template<class Type>
  39 void Foam::correlationFunction<Type>::setTimesAndSizes
  40 (
  41     const label tZeroBufferSize
  42 )
  43 {
  44     sampleSteps_  = ceil(sampleInterval_/mesh_.time().deltaT().value());
  45
  46     sampleInterval_ = sampleSteps_*mesh_.time().deltaT().value();
  47
  48     label bufferLength(ceil(duration_/sampleInterval_));
  49
  50     duration_ = bufferLength*sampleInterval_;
  51
  52     label bufferingInterval(ceil(averagingInterval_/sampleInterval_));
  53
  54     averagingInterval_ = bufferingInterval*sampleInterval_;
  55
  56     label nBuffers(ceil(duration_/averagingInterval_));
  57
  58     this->setSizes
  59     (
  60         nBuffers,
  61         bufferLength,
  62         bufferingInterval
  63     );
  64
  65     tZeroBuffers_ =
  66         Field< Field<Type> >
  67         (
  68             nBuffers,
  69             Field<Type>
  70             (
  71                 tZeroBufferSize,
  72                 pTraits<Type>::zero
  73             )
  74         );
  75 }
  76
  77
  78 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
  79
  80 template<class Type>
  81 Foam::correlationFunction<Type>::correlationFunction
  82 (
  83     const polyMesh& mesh,
  84     const dictionary& cfDict,
  85     const label tZeroBufferSize
  86 )
  87 :
  88     bufferedAccumulator<scalar>(),
  89     mesh_(mesh)
  90 {
  91     duration_ = readScalar(cfDict.lookup("duration"));
  92
  93     sampleInterval_ = readScalar(cfDict.lookup("sampleInterval"));
  94
  95     averagingInterval_ = readScalar(cfDict.lookup("averagingInterval"));
  96
  97     setTimesAndSizes(tZeroBufferSize);
  98 }
  99
 100
 101 template<class Type>
 102 Foam::correlationFunction<Type>::correlationFunction
 103 (
 104     const polyMesh& mesh,
 105     const label tZeroBufferSize,
 106     const scalar duration,
 107     const scalar sampleInterval,
 108     const scalar averagingInterval
 109 )
 110 :
 111     bufferedAccumulator<scalar>(),
 112     mesh_(mesh),
 113     duration_(duration),
 114     sampleInterval_(sampleInterval),
 115     averagingInterval_(averagingInterval)
 116 {
 117     setTimesAndSizes(tZeroBufferSize);
 118 }
 119
 120
 121 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 122
 123 template<class Type>
 124 Foam::correlationFunction<Type>::~correlationFunction()
 125 {}
 126
 127
 128 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 129
 130 template<class Type>
 131 void Foam::correlationFunction<Type>::calculateCorrelationFunction
 132 (
 133     const Field<Type>& currentValues
 134 )
 135 {
 136     if (measurandFieldSize() != currentValues.size())
 137     {
 138         FatalErrorIn("correlationFunction<Type>::calculateCorrelationFunction")
 139             << "Trying to supply a Field of length"
 140             << currentValues.size()
 141             << " to calculate the correlation function. "
 142             << "Expecting a Field of length "
 143             << measurandFieldSize() << nl
 144             << abort(FatalError);
 145     }
 146
 147     List<scalar> cFSums(nBuffers(),0.0);
 148
 149     forAll(tZeroBuffers_, tZB)
 150     {
 151         scalar& cFSum = cFSums[tZB];
 152
 153         const Field<Type>& tZeroBuffer = tZeroBuffers_[tZB];
 154
 155         forAll(currentValues, cV)
 156         {
 157             const Type& tZeroBufferValue = tZeroBuffer[cV];
 158
 159             const Type& currentValue = currentValues[cV];
 160
 161             forAll(currentValue, component)
 162             {
 163                 cFSum +=
 164                 (
 165                     tZeroBufferValue[component]*currentValue[component]
 166                 );
 167             }
 168         }
 169
 170         cFSum /= (measurandFieldSize()*currentValues[0].size());
 171     }
 172
 173     label bufferToRefill = addToBuffers(cFSums);
 174
 175     if (bufferToRefill != -1)
 176     {
 177         tZeroBuffers_[bufferToRefill] = currentValues;
 178     }
 179 }
 180
 181
 182 template<class Type>
 183 void Foam::correlationFunction<Type>::calculateCorrelationFunction
 184 (
 185     const Type& currentValue
 186 )
 187 {
 188     if( measurandFieldSize() != 1)
 189     {
 190         FatalErrorIn("correlationFunction<Type>::calculateCorrelationFunction")
 191             << "Trying to supply a single value to calculate the correlation "
 192             << "function.  Expecting a Field of length "
 193             << measurandFieldSize()
 194             << abort(FatalError);
 195     }
 196
 197     calculateCorrelationFunction(Field<Type>(1, currentValue));
 198 }
 199
 200
 201 template<class Type>
 202 Foam::scalar Foam::correlationFunction<Type>::integral() const
 203 {
 204     Field<scalar> averageCF(averaged());
 205
 206     scalar cFIntegral = 0.0;
 207
 208     for (label v = 0; v < averageCF.size() - 1; v++)
 209     {
 210         cFIntegral +=
 211             0.5
 212            *sampleInterval_
 213            *(averageCF[v+1] + averageCF[v]);
 214     }
 215
 216     return cFIntegral;
 217 }
 218
 219
 220 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 221
 222 #   include "correlationFunctionIO.C"
 223
 224 // ************************************************************************* //