Loading...
Searching...
No Matches
SteadyNSSimple.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/>.
24Class
25 steadyNS
26Description
27 Reduction class for a Stationary Navier-Stokes problem.
28SourceFiles
29 steadyNS.C
30\*---------------------------------------------------------------------------*/
31
36
37#ifndef SteadyNSSimple_H
38#define SteadyNSSimple_H
39#include "fvCFD.H"
40#include "singlePhaseTransportModel.H"
41#include "turbulentTransportModel.H"
42#include "simpleControl.H"
43#include "pisoControl.H"
44#include "fvOptions.H"
45#include "steadyNS.H"
46#include "reductionProblem.H"
47#include "ITHACAstream.H"
48#include "ITHACAparameters.H"
49#include "ITHACAforces.H"
50#include "volFields.H"
51#include <iostream>
52#include <datatable.h>
53#include <bspline.h>
54#include <bsplinebuilder.h>
55#include <rbfspline.h>
56#include <spline.h>
57
58
59// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
60
61/*---------------------------------------------------------------------------*\
62 Class SteadyNS Declaration
63\*---------------------------------------------------------------------------*/
64
66
70{
71
72 public:
73 // Constructors
76
78 SteadyNSSimple(int argc, char* argv[]);
80
82 PtrList<volScalarField> nutFields;
83
86
88 std::vector<SPLINTER::DataTable*> samples;
89
91 std::vector<SPLINTER::RBFSpline*> rbfSplines;
92
94 Eigen::MatrixXd coeffL2;
95
98
101
103 label saver;
104
106 label folderN;
107
108 // Variables
110 fvVectorMatrix* Ueqn_global;
112 fvScalarMatrix* Peqn_global;
113
114 //Functions
120 void getTurbRBF(label NNutModes);
126 void truthSolve2(List<scalar> mu_now, word Folder = "./ITHACAoutput/Offline/");
127};
128
129#endif
130
131
132
133
134
135
136
137
138
139
Header file of the ITHACAstream class, it contains the implementation of several methods for input ou...
Implementation of a parametrized full order steady NS problem and preparation of the the reduced ma...
std::vector< SPLINTER::DataTable * > samples
Create a samples for interpolation.
std::vector< SPLINTER::RBFSpline * > rbfSplines
Create a RBF splines for interpolation.
SteadyNSSimple()
Null constructor.
label folderN
Counter to save intermediate steps in the correct folder (for turbulent case only)
void getTurbRBF(label NNutModes)
Function to calculate RBF weights for turbulence.
void truthSolve2(List< scalar > mu_now, word Folder="./ITHACAoutput/Offline/")
Offline solver for the whole Navier-Stokes problem.
volScalarModes nutModes
List of POD modes for eddy viscosity.
bool middleExport
Export also intermediate fields.
Eigen::MatrixXd coeffL2
The matrix of L2 projection coefficients for the eddy viscosity.
fvScalarMatrix * Peqn_global
Initialization for the full pressure linear system.
label saver
Counter to check if the middleStep has been reached or not (for turbulent case only)
fvVectorMatrix * Ueqn_global
Initialization for the full velocity linear system.
label middleStep
Distancing between intermediate steps (for turbulent case only)
PtrList< volScalarField > nutFields
List of snapshots for the solution for eddy viscosity.
Implementation of a parametrized full order steady NS problem and preparation of the the reduced ma...
Definition steadyNS.H:70
label NNutModes
Number of nut modes used for the projection.
Definition steadyNS.H:149
Header file of the reductionProblem class.
Header file of the steadyNS class.