Loading...
Searching...
No Matches
ITHACAassign.H
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ██╗████████╗██╗ ██╗ █████╗ ██████╗ █████╗ ███████╗██╗ ██╗
3 ██║╚══██╔══╝██║ ██║██╔══██╗██╔════╝██╔══██╗ ██╔════╝██║ ██║
4 ██║ ██║ ███████║███████║██║ ███████║█████╗█████╗ ██║ ██║
5 ██║ ██║ ██╔══██║██╔══██║██║ ██╔══██║╚════╝██╔══╝ ╚██╗ ██╔╝
6 ██║ ██║ ██║ ██║██║ ██║╚██████╗██║ ██║ ██║ ╚████╔╝
7 ╚═╝ ╚═╝ ╚═╝ ╚═╝╚═╝ ╚═╝ ╚═════╝╚═╝ ╚═╝ ╚═╝ ╚═══╝
8
9 * In real Time Highly Advanced Computational Applications for Finite Volumes
10 * Copyright (C) 2017 by the ITHACA-FV authors
11-------------------------------------------------------------------------------
12License
13 This file is part of ITHACA-FV
14 ITHACA-FV is free software: you can redistribute it and/or modify
15 it under the terms of the GNU Lesser General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18 ITHACA-FV is distributed in the hope that it will be useful,
19 but WITHOUT ANY WARRANTY; without even the implied warranty of
20 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
21 GNU Lesser General Public License for more details.
22 You should have received a copy of the GNU Lesser General Public License
23 along with ITHACA-FV. If not, see <http://www.gnu.org/licenses/>.
24Namespace
25 ITHACAutilites
26Description
27 Utilities to assign fields and BCs of OF fields
28SourceFiles
29 ITHACAassign.C
30\*---------------------------------------------------------------------------*/
31
36
37#ifndef ITHACAassign_H
38#define ITHACAassign_H
39
40#include "fvCFD.H"
41#include "IOmanip.H"
42#include "freestreamFvPatchField.H"
43#include <sys/stat.h>
44#include <unistd.h>
45#pragma GCC diagnostic push
46#pragma GCC diagnostic ignored "-Wold-style-cast"
47#include <Eigen/Eigen>
48#pragma GCC diagnostic pop
49#include <functional>
50#include "./colormod.H"
51#include "polyMeshTools.H"
52#include <chrono>
53#include "mixedFvPatchFields.H"
54#include "fvMeshSubset.H"
55using namespace std::placeholders;
56#include "Foam2Eigen.H"
57#include "ITHACAerror.H"
58#include "ITHACAparameters.H"
59
61namespace ITHACAutilities
62{
63
64//--------------------------------------------------------------------------
72template<typename Type>
73void assignIF(GeometricField<Type, fvPatchField, volMesh>& field,
74 Type value);
75
76//--------------------------------------------------------------------------
85template<typename Type>
86void assignIF(GeometricField<Type, fvPatchField, volMesh>& field,
87 Type& value, List<label>& indices);
88
89//--------------------------------------------------------------------------
98template<typename Type>
99void assignIF(GeometricField<Type, fvPatchField, volMesh>& field,
100 Type& value, label index);
101
102//--------------------------------------------------------------------------
109void assignONE(volScalarField& s, List<label>& L);
110
111
112//--------------------------------------------------------------------------
119void assignBC(GeometricField<scalar, fvPatchField, volMesh>& s, label BC_ind,
120 double value);
121
122//--------------------------------------------------------------------------
129void assignBC(GeometricField<scalar, fvPatchField, volMesh>& s, label BC_ind,
130 Eigen::MatrixXd valueVec);
131
132//--------------------------------------------------------------------------
139void assignBC(GeometricField<scalar, fvPatchField, volMesh>& s, label BC_ind,
140 List<double> valueList);
141//--------------------------------------------------------------------------
148
149void assignBC(GeometricField<vector, pointPatchField, pointMesh>& s,
150 label BC_ind,
151 List<double> valueList);
152
153//--------------------------------------------------------------------------
160void assignBC(GeometricField<vector, fvPatchField, volMesh>& s, label BC_ind,
161 vector value);
162
163//--------------------------------------------------------------------------
170void assignBC(GeometricField<vector, pointPatchField, pointMesh>& s,
171 label BC_ind,
172 vector value);
173
174//--------------------------------------------------------------------------
181void assignBC(GeometricField<tensor, fvPatchField, volMesh>& s, label BC_ind,
182 tensor value);
183
184//--------------------------------------------------------------------------
191void assignBC(GeometricField<vector, fvPatchField, volMesh>& s, label BC_ind,
192 Eigen::MatrixXd valueVec);
193
194//--------------------------------------------------------------------------
201void assignBC(GeometricField<tensor, fvPatchField, volMesh>& s, label BC_ind,
202 Eigen::MatrixXd valueVec);
203
204//--------------------------------------------------------------------------
211void assignBC(GeometricField<vector, pointPatchField, pointMesh>& s,
212 label BC_ind,
213 Eigen::MatrixXd valueVec);
214
215//--------------------------------------------------------------------------
222void assignBC(GeometricField<vector, fvPatchField, volMesh>& s, label BC_ind,
223 List<vector> valueList);
224
225//--------------------------------------------------------------------------
232void assignBC(GeometricField<vector, pointPatchField, pointMesh>& s,
233 label BC_ind,
234 List<vector> valueList);
235
236//--------------------------------------------------------------------------
243void assignBC(GeometricField<tensor, fvPatchField, volMesh>& s, label BC_ind,
244 List<tensor> valueList);
245
246//------------------------------------------------------------------------------
255void assignBC(GeometricField<scalar, fvsPatchField, surfaceMesh>& field,
256 label BC_ind, Eigen::MatrixXd value);
257
258void assignBC(GeometricField<vector, fvsPatchField, surfaceMesh>& field,
259 label BC_ind, Eigen::MatrixXd value);
260
261//------------------------------------------------------------------------------
270template<typename Type>
271void assignBC(GeometricField<Type, fvsPatchField, surfaceMesh>& field,
272 label BC_ind, List<Type>& value);
273
274
275
276//------------------------------------------------------------------------------
285template<typename Type>
286void assignBC(GeometricField<Type, fvsPatchField, surfaceMesh>& field,
287 label BC_ind, Type& value);
288
289
290//--------------------------------------------------------------------------
298template<typename Type>
299void changeNeumann2Dirichlet(GeometricField<Type, fvPatchField, volMesh>&
300 field,
301 Type& value);
302
303
304//--------------------------------------------------------------------------
313template<typename Type>
314void assignZeroDirichlet(GeometricField<Type, fvPatchField,
315 volMesh>& field);
316
317#pragma GCC diagnostic push
318#pragma GCC diagnostic ignored "-Wcomment"
319//------------------------------------------------------------------------------
340template<typename Type>
341void setBoxToValue(GeometricField<Type, fvPatchField, volMesh>& field,
342 Eigen::MatrixXd Box, Type value);
343#pragma GCC diagnostic pop
344
345//--------------------------------------------------------------------------
354template<class Type>
355void changeBCtype(GeometricField<Type, fvPatchField, volMesh>&
356 field, word BCtype, label BC_ind);
357
358
359//--------------------------------------------------------------------------
369template<typename Type>
370void setIndices2Value(labelList& ind2set, List<Type>& value2set,
371 labelList& movingIDS, List<Type>& originalList);
372
373//--------------------------------------------------------------------------
384template<typename Type>
385void assignMixedBC(GeometricField<Type, fvPatchField, volMesh>& field,
386 label BC_ind, Eigen::MatrixXd& value,
387 Eigen::MatrixXd& grad, Eigen::MatrixXd& valueFrac);
388
389
390//--------------------------------------------------------------------------
401template<typename Type>
402void assignMixedBC(GeometricField<Type, fvPatchField, volMesh>& field,
403 label BC_ind, List<Type>& value,
404 List<Type>& grad, List<scalar>& valueFrac);
405
406
407//------------------------------------------------------------------------------
412template<typename T>
413void setToZero(T& f1);
414
415}
416
417#endif
Header file of the Foam2Eigen class.
Header file of the ITHACAerror file.
Simple header and source file of the Color::Modifier class to change color to the output stream.
Namespace to implement some useful assign operation of OF fields.
void assignMixedBC(GeometricField< Type, fvPatchField, volMesh > &field, label BC_ind, List< Type > &value, List< Type > &grad, List< scalar > &valueFrac)
Assign value of a boundary condition of type "mixed".
void assignIF(GeometricField< Type, fvPatchField, volMesh > &s, Type value)
Assign internal field.
void changeNeumann2Dirichlet(GeometricField< Type, fvPatchField, volMesh > &field, Type &value)
Change all Neumann boundary conditions to Dirichlet boundary conditions.
void setToZero(T &f1)
Set a field of type vol[Scalar|Vector|Tensor]Field to 0.
void assignONE(volScalarField &s, List< label > &L)
Assign one to volScalarField.
void setIndices2Value(labelList &ind2set, List< Type > &value2set, labelList &movingIDS, List< Type > &originalList)
Sets some given Indices of a list of objects to given values.
void assignZeroDirichlet(GeometricField< Type, fvPatchField, volMesh > &field)
Assign zero internal field.
void assignBC(GeometricField< scalar, fvPatchField, volMesh > &s, label BC_ind, double value)
Assign uniform Boundary Condition to a volScalarField.
void changeBCtype(GeometricField< Type, fvPatchField, volMesh > &field, word BCtype, label BC_ind)
Change the boundary condition type for a GeometricField.
void setBoxToValue(GeometricField< Type, fvPatchField, volMesh > &field, Eigen::MatrixXd Box, Type value)
Set value of a volScalarField to a constant inside a given box.