applications/utilities/preProcessing/setFields/setFields.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 1991-2010 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
  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     Selects a cell set through a dictionary.
  26
  27 \*---------------------------------------------------------------------------*/
  28
  29 #include "argList.H"
  30 #include "timeSelector.H"
  31 #include "Time.H"
  32 #include "fvMesh.H"
  33 #include "topoSetSource.H"
  34 #include "cellSet.H"
  35 #include "volFields.H"
  36
  37 using namespace Foam;
  38
  39 template<class GeoField>
  40 void setFieldType
  41 (
  42     const fvMesh& mesh,
  43     const labelList& selectedCells,
  44     Istream& fieldValueStream
  45 )
  46 {
  47     word fieldName(fieldValueStream);
  48
  49     IOobject fieldHeader
  50     (
  51         fieldName,
  52         mesh.time().timeName(),
  53         mesh,
  54         IOobject::MUST_READ
  55     );
  56
  57     // Check field exists
  58     if (fieldHeader.headerOk())
  59     {
  60         Info<< "    Setting " << fieldHeader.headerClassName()
  61             << " " << fieldName << endl;
  62
  63         GeoField field(fieldHeader, mesh);
  64
  65         typename GeoField::value_type value
  66         (
  67             static_cast<const typename GeoField::value_type&>
  68             (
  69                 pTraits<typename GeoField::value_type>(fieldValueStream)
  70             )
  71         );
  72
  73         if (selectedCells.size() == field.size())
  74         {
  75             field.internalField() = value;
  76         }
  77         else
  78         {
  79             forAll(selectedCells, celli)
  80             {
  81                 field[selectedCells[celli]] = value;
  82             }
  83         }
  84
  85         forAll(field.boundaryField(), patchi)
  86         {
  87             field.boundaryField()[patchi] =
  88                 field.boundaryField()[patchi].patchInternalField();
  89         }
  90
  91         field.write();
  92     }
  93     else
  94     {
  95         WarningIn
  96         (
  97             "void setFieldType"
  98             "(const fvMesh& mesh, const labelList& selectedCells,"
  99             "Istream& fieldValueStream)"
 100         ) << "Field " << fieldName << " not found" << endl;
 101     }
 102 }
 103
 104
 105 class setField
 106 {
 107
 108 public:
 109
 110     setField()
 111     {}
 112
 113     autoPtr<setField> clone() const
 114     {
 115         return autoPtr<setField>(new setField());
 116     }
 117
 118     class iNew
 119     {
 120         const fvMesh& mesh_;
 121         const labelList& selectedCells_;
 122
 123     public:
 124
 125         iNew(const fvMesh& mesh, const labelList& selectedCells)
 126         :
 127             mesh_(mesh),
 128             selectedCells_(selectedCells)
 129         {}
 130
 131         autoPtr<setField> operator()(Istream& fieldValues) const
 132         {
 133             word fieldType(fieldValues);
 134
 135             if (fieldType == "volScalarFieldValue")
 136             {
 137                 setFieldType<volScalarField>
 138                     (mesh_, selectedCells_, fieldValues);
 139             }
 140             else if (fieldType == "volVectorFieldValue")
 141             {
 142                 setFieldType<volVectorField>
 143                     (mesh_, selectedCells_, fieldValues);
 144             }
 145             else if (fieldType == "volSphericalTensorFieldValue")
 146             {
 147                 setFieldType<volSphericalTensorField>
 148                     (mesh_, selectedCells_, fieldValues);
 149             }
 150             else if (fieldType == "volSymmTensorFieldValue")
 151             {
 152                 setFieldType<volSymmTensorField>
 153                     (mesh_, selectedCells_, fieldValues);
 154             }
 155             else if (fieldType == "volTensorFieldValue")
 156             {
 157                 setFieldType<volTensorField>
 158                     (mesh_, selectedCells_, fieldValues);
 159             }
 160             else
 161             {
 162                 WarningIn("setField::iNew::operator()(Istream& is)")
 163                     << "field type " << fieldType << " not currently supported"
 164                     << endl;
 165             }
 166
 167             return autoPtr<setField>(new setField());
 168         }
 169     };
 170 };
 171
 172
 173 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 174
 175 int main(int argc, char *argv[])
 176 {
 177     timeSelector::addOptions();
 178
 179 #   include "setRootCase.H"
 180 #   include "createTime.H"
 181
 182     // Get times list
 183     instantList timeDirs = timeSelector::select0(runTime, args);
 184
 185 #   include "createMesh.H"
 186
 187     Info<< "Reading setFieldsDict\n" << endl;
 188
 189     IOdictionary setFieldsDict
 190     (
 191         IOobject
 192         (
 193             "setFieldsDict",
 194             runTime.system(),
 195             mesh,
 196             IOobject::MUST_READ,
 197             IOobject::NO_WRITE
 198         )
 199     );
 200
 201     if (setFieldsDict.found("defaultFieldValues"))
 202     {
 203         Info<< "Setting field default values" << endl;
 204         PtrList<setField> defaultFieldValues
 205         (
 206             setFieldsDict.lookup("defaultFieldValues"),
 207             setField::iNew(mesh, labelList(mesh.nCells()))
 208         );
 209         Info<< endl;
 210     }
 211
 212
 213     Info<< "Setting field region values" << endl;
 214
 215     PtrList<entry> regions(setFieldsDict.lookup("regions"));
 216
 217     forAll(regions, regionI)
 218     {
 219         const entry& region = regions[regionI];
 220
 221         autoPtr<topoSetSource> cellSelector =
 222             topoSetSource::New(region.keyword(), mesh, region.dict());
 223
 224         cellSet selectedCellSet
 225         (
 226             mesh,
 227             "cellSet",
 228             mesh.nCells()/10+1  // Reasonable size estimate.
 229         );
 230
 231         cellSelector->applyToSet
 232         (
 233             topoSetSource::NEW,
 234             selectedCellSet
 235         );
 236
 237         PtrList<setField> fieldValues
 238         (
 239             region.dict().lookup("fieldValues"),
 240             setField::iNew(mesh, selectedCellSet.toc())
 241         );
 242     }
 243
 244     Info<< "\nEnd" << endl;
 245
 246     return 0;
 247 }
 248
 249
 250 // ************************************************************************* //