applications/utilities/postProcessing/miscellaneous/temporalInterpolate/temporalInterpolate.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
   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
  13     the Free Software Foundation, either version 3 of the License, or
  14     (at your 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, see <http://www.gnu.org/licenses/>.
  23
  24 Description
  25     Interpolate fields between time-steps e.g. for animation.
  26
  27 \*---------------------------------------------------------------------------*/
  28
  29 #include "argList.H"
  30 #include "timeSelector.H"
  31
  32 #include "fvMesh.H"
  33 #include "Time.H"
  34 #include "volMesh.H"
  35 #include "surfaceMesh.H"
  36 #include "volFields.H"
  37 #include "surfaceFields.H"
  38 #include "pointFields.H"
  39 #include "ReadFields.H"
  40
  41 using namespace Foam;
  42
  43 class fieldInterpolator
  44 {
  45     Time& runTime_;
  46     const fvMesh& mesh_;
  47     const IOobjectList& objects_;
  48     const HashSet<word>& selectedFields_;
  49     instant ti_;
  50     instant ti1_;
  51     int divisions_;
  52
  53 public:
  54
  55     fieldInterpolator
  56     (
  57         Time& runTime,
  58         const fvMesh& mesh,
  59         const IOobjectList& objects,
  60         const HashSet<word>& selectedFields,
  61         const instant& ti,
  62         const instant& ti1,
  63         int divisions
  64     )
  65     :
  66         runTime_(runTime),
  67         mesh_(mesh),
  68         objects_(objects),
  69         selectedFields_(selectedFields),
  70         ti_(ti),
  71         ti1_(ti1),
  72         divisions_(divisions)
  73     {}
  74
  75     template<class GeoFieldType>
  76     void interpolate();
  77 };
  78
  79
  80 template<class GeoFieldType>
  81 void fieldInterpolator::interpolate()
  82 {
  83     const word& fieldClassName = GeoFieldType::typeName;
  84
  85     IOobjectList fields = objects_.lookupClass(fieldClassName);
  86
  87     if (fields.size())
  88     {
  89         Info<< "    " << fieldClassName << "s:";
  90
  91         forAllConstIter(IOobjectList, fields, fieldIter)
  92         {
  93             if
  94             (
  95                 selectedFields_.empty()
  96              || selectedFields_.found(fieldIter()->name())
  97             )
  98             {
  99                 Info<< " " << fieldIter()->name() << '(';
 100
 101                 GeoFieldType fieldi
 102                 (
 103                     IOobject
 104                     (
 105                         fieldIter()->name(),
 106                         ti_.name(),
 107                         fieldIter()->db(),
 108                         IOobject::MUST_READ,
 109                         IOobject::NO_WRITE,
 110                         false
 111                     ),
 112                     mesh_
 113                 );
 114
 115                 GeoFieldType fieldi1
 116                 (
 117                     IOobject
 118                     (
 119                         fieldIter()->name(),
 120                         ti1_.name(),
 121                         fieldIter()->db(),
 122                         IOobject::MUST_READ,
 123                         IOobject::NO_WRITE,
 124                         false
 125                     ),
 126                     mesh_
 127                 );
 128
 129                 scalar deltaT = (ti1_.value() - ti_.value())/(divisions_ + 1);
 130
 131                 for (int j=0; j<divisions_; j++)
 132                 {
 133                     instant timej = instant(ti_.value() + (j + 1)*deltaT);
 134
 135                     runTime_.setTime(timej.name(), 0);
 136
 137                     Info<< timej.name();
 138
 139                     if (j < divisions_-1)
 140                     {
 141                         Info<< " ";
 142                     }
 143
 144                     scalar lambda = scalar(j + 1)/scalar(divisions_ + 1);
 145
 146                     GeoFieldType fieldj
 147                     (
 148                         IOobject
 149                         (
 150                             fieldIter()->name(),
 151                             timej.name(),
 152                             fieldIter()->db(),
 153                             IOobject::NO_READ,
 154                             IOobject::NO_WRITE,
 155                             false
 156                         ),
 157                         (1.0 - lambda)*fieldi + lambda*fieldi1
 158                     );
 159
 160                     fieldj.write();
 161                 }
 162
 163                 Info<< ')';
 164             }
 165         }
 166
 167         Info<< endl;
 168     }
 169 }
 170
 171
 172 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 173
 174 // Main program:
 175 int main(int argc, char *argv[])
 176 {
 177     timeSelector::addOptions();
 178     argList::addOption
 179     (
 180         "fields",
 181         "list",
 182         "specify a list of fields to be interpolated. Eg, '(U T p)' - "
 183         "regular expressions not currently supported"
 184     );
 185     argList::addOption
 186     (
 187         "divisions",
 188         "integer",
 189         "specify number of temporal sub-divisions to create (default = 1)."
 190     );
 191
 192     #include "setRootCase.H"
 193     #include "createTime.H"
 194     runTime.functionObjects().off();
 195
 196     HashSet<word> selectedFields;
 197     if (args.optionFound("fields"))
 198     {
 199         args.optionLookup("fields")() >> selectedFields;
 200     }
 201
 202     int divisions = 1;
 203     if (args.optionFound("divisions"))
 204     {
 205         args.optionLookup("divisions")() >> divisions;
 206     }
 207
 208     instantList timeDirs = timeSelector::select0(runTime, args);
 209
 210     #include "createMesh.H"
 211
 212     Info<< "Interpolating fields for times:" << endl;
 213
 214     for (label timei = 0; timei < timeDirs.size() - 1; timei++)
 215     {
 216         runTime.setTime(timeDirs[timei], timei);
 217
 218         // Read objects in time directory
 219         IOobjectList objects(mesh, runTime.timeName());
 220
 221         fieldInterpolator interpolator
 222         (
 223             runTime,
 224             mesh,
 225             objects,
 226             selectedFields,
 227             timeDirs[timei],
 228             timeDirs[timei+1],
 229             divisions
 230         );
 231
 232         // Interpolate vol fields
 233         interpolator.interpolate<volScalarField>();
 234         interpolator.interpolate<volVectorField>();
 235         interpolator.interpolate<volSphericalTensorField>();
 236         interpolator.interpolate<volSymmTensorField>();
 237         interpolator.interpolate<volTensorField>();
 238
 239         // Interpolate surface fields
 240         interpolator.interpolate<surfaceScalarField>();
 241         interpolator.interpolate<surfaceVectorField>();
 242         interpolator.interpolate<surfaceSphericalTensorField>();
 243         interpolator.interpolate<surfaceSymmTensorField>();
 244         interpolator.interpolate<surfaceTensorField>();
 245     }
 246
 247     Info<< "End\n" << endl;
 248
 249     return 0;
 250 }
 251
 252
 253 // ************************************************************************* //