Loading...
Searching...
No Matches
ITHACAerror.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-------------------------------------------------------------------------------
12
13License
14 This file is part of ITHACA-FV
15
16 ITHACA-FV is free software: you can redistribute it and/or modify
17 it under the terms of the GNU Lesser General Public License as published by
18 the Free Software Foundation, either version 3 of the License, or
19 (at your option) any later version.
20
21 ITHACA-FV is distributed in the hope that it will be useful,
22 but WITHOUT ANY WARRANTY; without even the implied warranty of
23 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
24 GNU Lesser General Public License for more details.
25
26 You should have received a copy of the GNU Lesser General Public License
27 along with ITHACA-FV. If not, see <http://www.gnu.org/licenses/>.
28\*---------------------------------------------------------------------------*/
29
30#ifndef ITHACAerror_H
31#define ITHACAerror_H
32
35
36
37#include "fvCFD.H"
38#include "IOmanip.H"
39#include "freestreamFvPatchField.H"
40#include <sys/stat.h>
41#include <unistd.h>
42#pragma GCC diagnostic push
43#pragma GCC diagnostic ignored "-Wold-style-cast"
44#include <Eigen/Eigen>
45#pragma GCC diagnostic pop
46#include <functional>
47#include "./colormod.H"
48#include "polyMeshTools.H"
49#include <chrono>
50#include "mixedFvPatchFields.H"
51#include "fvMeshSubset.H"
52using namespace std::placeholders;
53#include "Foam2Eigen.H"
54#include "ITHACAnorm.H"
55
56namespace ITHACAutilities
57{
58
59//------------------------------------------------------------------------------
72template<class Type, template<class> class PatchField, class GeoMesh>
73double errorFrobRel(GeometricField<Type, PatchField, GeoMesh>& field1,
74 GeometricField<Type, PatchField, GeoMesh>& field2, List<label>* labels = NULL);
75
76//--------------------------------------------------------------------------
87template<typename T>
88double errorL2Rel(GeometricField<T, fvPatchField, volMesh>& field1,
89 GeometricField<T, fvPatchField, volMesh>& field2, List<label>* labels = NULL);
90
91//--------------------------------------------------------------------------
102template<typename T>
103double errorLinfRel(GeometricField<T, fvPatchField, volMesh>& field1,
104 GeometricField<T, fvPatchField, volMesh>& field2, List<label>* labels = NULL);
105
106//--------------------------------------------------------------------------
118template<typename T>
119double errorL2Abs(GeometricField<T, fvPatchField, volMesh>& field1,
120 GeometricField<T, fvPatchField, volMesh>& field2, volScalarField& Volumes);
121
122//--------------------------------------------------------------------------
133template<typename T>
134double errorL2Abs(GeometricField<T, fvPatchField, volMesh>& field1,
135 GeometricField<T, fvPatchField, volMesh>& field2, List<label>* labels = NULL);
136
137//--------------------------------------------------------------------------
148template<typename T>
149Eigen::MatrixXd errorL2Rel(PtrList<GeometricField<T, fvPatchField, volMesh >>&
150 fields1,
151 PtrList<GeometricField<T, fvPatchField, volMesh >>& fields2,
152 List<label>* labels = NULL);
153
154//--------------------------------------------------------------------------
165template<class T, template<class> class PatchField, class GeoMesh>
166Eigen::MatrixXd errorFrobRel(PtrList<GeometricField<T, PatchField, GeoMesh >>&
167 fields1,
168 PtrList<GeometricField<T, PatchField, GeoMesh >>& fields2,
169 List<label>* labels = NULL);
170//--------------------------------------------------------------------------
171
183template<typename T>
184Eigen::MatrixXd errorL2Abs(
185 PtrList<GeometricField<T, fvPatchField, volMesh >>& fields1,
186 PtrList<GeometricField<T, fvPatchField, volMesh >>& fields2,
187 PtrList<volScalarField>& Volumes);
188
189//--------------------------------------------------------------------------
200template<typename T>
201Eigen::MatrixXd errorL2Abs(PtrList<GeometricField<T, fvPatchField, volMesh >>&
202 fields1,
203 PtrList<GeometricField<T, fvPatchField, volMesh >>& fields2,
204 List<label>* labels = NULL);
205
206}
207
208
209#endif
Header file of the Foam2Eigen class.
Header file of the ITHACAnorm 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.
double errorLinfRel(GeometricField< T, fvPatchField, volMesh > &field1, GeometricField< T, fvPatchField, volMesh > &field2, List< label > *labels)
Computes the relative error between two geometric Fields in Linf norm.
double errorL2Rel(GeometricField< T, fvPatchField, volMesh > &field1, GeometricField< T, fvPatchField, volMesh > &field2, List< label > *labels)
Computes the relative error between two geometric Fields in L2 norm.
double errorFrobRel(GeometricField< Type, PatchField, GeoMesh > &field1, GeometricField< Type, PatchField, GeoMesh > &field2, List< label > *labels)
Computes the relative error between two Fields in the Frobenius norm.
Definition ITHACAerror.C:41