steadyNS.H

Go to the documentation of this file.

/*---------------------------------------------------------------------------*\
     ██╗████████╗██╗  ██╗ █████╗  ██████╗ █████╗       ███████╗██╗   ██╗
     ██║╚══██╔══╝██║  ██║██╔══██╗██╔════╝██╔══██╗      ██╔════╝██║   ██║
     ██║   ██║   ███████║███████║██║     ███████║█████╗█████╗  ██║   ██║
     ██║   ██║   ██╔══██║██╔══██║██║     ██╔══██║╚════╝██╔══╝  ╚██╗ ██╔╝
     ██║   ██║   ██║  ██║██║  ██║╚██████╗██║  ██║      ██║      ╚████╔╝
     ╚═╝   ╚═╝   ╚═╝  ╚═╝╚═╝  ╚═╝ ╚═════╝╚═╝  ╚═╝      ╚═╝       ╚═══╝
 
 * In real Time Highly Advanced Computational Applications for Finite Volumes
 * Copyright (C) 2017 by the ITHACA-FV authors
-------------------------------------------------------------------------------
License
    This file is part of ITHACA-FV
    ITHACA-FV is free software: you can redistribute it and/or modify
    it under the terms of the GNU Lesser General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.
    ITHACA-FV is distributed in the hope that it will be useful,
    but WITHOUT ANY WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
    GNU Lesser General Public License for more details.
    You should have received a copy of the GNU Lesser General Public License
    along with ITHACA-FV. If not, see <http://www.gnu.org/licenses/>.
Class
    steadyNS
Description
    Reduction class for a Stationary Navier-Stokes problem.
SourceFiles
    steadyNS.C
\*---------------------------------------------------------------------------*/

 
#ifndef steadyNS_H
#define steadyNS_H
#include "fvCFD.H"
#include "singlePhaseTransportModel.H"
#include "turbulentTransportModel.H"
#include "simpleControl.H"
#include "pisoControl.H"
#include "fvOptions.H"
#include "reductionProblem.H"
#include "ITHACAstream.H"
#include "ITHACAparameters.H"
#if OPENFOAM >= 1812
#include "ITHACAforces18.H"
#else
#include "ITHACAforces.H"
#endif
#include "volFields.H"
#include <iostream>
#include "IPstream.H"
#include "OPstream.H"
#include "Modes.H"
 
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
/*---------------------------------------------------------------------------*\
                        Class SteadyNS Declaration
\*---------------------------------------------------------------------------*/


class steadyNS: public reductionProblem
{
 
 
    public:
        // Constructors
        steadyNS();

        steadyNS(int argc, char* argv[]);
        ~steadyNS() {};
 
        ITHACAparameters* para;
 
        // Member Functions
        PtrList<volScalarField> Pfield;

        PtrList<volVectorField> Ufield;

        PtrList<volVectorField> supfield;

        PtrList<surfaceScalarField> Phifield;

        volScalarModes Pmodes;

        volVectorModes Umodes;

        volVectorModes supmodes;

        surfaceScalarModes Phimodes;

        PtrList<volVectorField> liftfield;

        PtrList<volScalarField> liftfieldP;

        PtrList<volVectorField> Uomfield;

        volVectorModes L_U_SUPmodes;

        surfaceScalarModes L_PHImodes;

        scalar tolerance;

        scalar maxIter;

        label NUmodesOut;

        label NPmodesOut;

        label NSUPmodesOut;

        label NNutModesOut;

        label NUmodes;

        label NPmodes;

        label NSUPmodes;

        label NNutModes;


        Eigen::MatrixXd B_matrix;

        Eigen::MatrixXd M_matrix;

        Eigen::MatrixXd K_matrix;

        List <Eigen::MatrixXd> C_matrix;
        Eigen::Tensor<double, 3 > C_tensor;

        Eigen::MatrixXd P_matrix;

        Eigen::MatrixXd L_matrix;

        Eigen::MatrixXd L_D_matrix;

        Eigen::MatrixXd D_matrix;

        List <Eigen::MatrixXd> G_matrix;

        Eigen::Tensor<double, 3 > gTensor;

        Eigen::MatrixXd BC1_matrix;

        List <Eigen::MatrixXd> BC2_matrix;

        Eigen::Tensor<double, 3 > bc2Tensor;

        Eigen::MatrixXd BC3_matrix;

        Eigen::MatrixXd BC4_matrix;

        Eigen::MatrixXd W_matrix;

        Eigen::MatrixXd I_matrix;

        Eigen::MatrixXd DF_matrix;

        Eigen::MatrixXd KF_matrix;

        //

        Eigen::MatrixXd tauMatrix;

        Eigen::MatrixXd nMatrix;

        List <Eigen::MatrixXd> bcVelVec;

        List <Eigen::MatrixXd> bcVelMat;

        Eigen::MatrixXd BP_matrix;

        List<Eigen::MatrixXd> RD_matrix;

        List<Eigen::MatrixXd> RC_matrix;

        List<Eigen::MatrixXd> SD_matrix;

        List<Eigen::MatrixXd> SC_matrix;

        Eigen::Tensor<double, 3 > Cf_tensor;

        Eigen::Tensor<double, 3 > Ci_tensor;

        List<Eigen::MatrixXd> LinSysDiv;

        List<Eigen::MatrixXd> LinSysDiff;

        List<Eigen::MatrixXd> LinSysConv;
 
 
        // Other Variables
        bool supex;
 
        // Dummy variables to transform simplefoam into a class
        autoPtr<volScalarField> _p;

        autoPtr<volVectorField> _U;

        autoPtr<volScalarField> _p0;

        autoPtr<volVectorField> _U0;

        autoPtr<volVectorField> Uinl;

        autoPtr<dimensionedScalar> dt_dummy;

        autoPtr<dimensionedScalar> nu_dummy;

        mutable autoPtr<fvMesh> _mesh;

        autoPtr<simpleControl> _simple;

        autoPtr<fv::options> _fvOptions;

        autoPtr<Time> _runTime;

        autoPtr<surfaceScalarField> _phi;

        autoPtr<surfaceScalarField> _phi0;

        autoPtr<incompressible::turbulenceModel> turbulence;

        autoPtr<singlePhaseTransportModel> _laminarTransport;

        autoPtr<IOMRFZoneList> _MRF;

        label pRefCell;

        scalar pRefValue;

        scalar cumulativeContErr = 0;

        word bcMethod;

        word fluxMethod;
 
        // Functions
 
        //--------------------------------------------------------------------------
        void truthSolve(List<scalar> mu_now);

        void solvesupremizer(word type = "snapshots");

        void liftSolve();
 
        // Wrapped Proj. Methods;
 
        //--------------------------------------------------------------------------
        void projectPPE(fileName folder, label NUmodes, label NPmodes,
                        label NSUPmodes = 0);
 
        //--------------------------------------------------------------------------
        void projectSUP(fileName folder, label NUmodes, label NPmodes, label NSUPmodes);
 
        //--------------------------------------------------------------------------
        void discretizeThenProject(fileName folder, label NUmodes, label NPmodes,
                                   label NSUPmodes = 0);
 
        //--------------------------------------------------------------------------
        //  Projection Methods Momentum Equation
        Eigen::MatrixXd diffusive_term(label NUmodes, label NPmodes, label NSUPmodes);
 
        //--------------------------------------------------------------------------
        //  Projection Methods Momentum Equation
        Eigen::MatrixXd diffusive_term_sym(label NUmodes, label NPmodes,
                                           label NSUPmodes);
 
        //--------------------------------------------------------------------------
        Eigen::MatrixXd pressure_gradient_term(label NUmodes, label NPmodes,
                                               label NSUPmodes);
 
        //--------------------------------------------------------------------------
        List < Eigen::MatrixXd >  convective_term(label NUmodes, label NPmodes,
                label NSUPmodes);
 
        //--------------------------------------------------------------------------
        Eigen::MatrixXd  mass_term(label NUmodes, label NPmodes, label NSUPmodes);
 
        // Projection Methods Continuity Equation
 
        //--------------------------------------------------------------------------
        Eigen::MatrixXd divergence_term(label NUmodes, label NPmodes, label NSUPmodes);
 
        //--------------------------------------------------------------------------
        List < Eigen::MatrixXd > div_momentum(label NUmodes, label NPmodes);

        Eigen::Tensor<double, 3 > divMomentum(label NUmodes, label NPmodes);
        
        //--------------------------------------------------------------------------
        Eigen::Tensor<double, 3 > divMomentum_cache(label NUmodes, label NPmodes);
 
        //--------------------------------------------------------------------------
        Eigen::MatrixXd laplacian_pressure(label NPmodes);
 
        //--------------------------------------------------------------------------
        Eigen::MatrixXd convective_background(label NUmodes, volVectorField vls);
 
        //--------------------------------------------------------------------------
        Eigen::MatrixXd divergent_convective_background(label NPmodes, label NUmodes,
                volVectorField vls);
 
        //--------------------------------------------------------------------------
        Eigen::MatrixXd pressure_BC1(label NPmodes, label NUmodes);
 
        //--------------------------------------------------------------------------
        List < Eigen::MatrixXd > pressure_BC2(label NPmodes, label NUmodes);
 
        //--------------------------------------------------------------------------
        Eigen::Tensor<double, 3 > pressureBC2(label NPmodes, label NUmodes);
 
        //--------------------------------------------------------------------------
        Eigen::MatrixXd pressure_BC3(label NPmodes, label NUmodes);
 
        //--------------------------------------------------------------------------
        Eigen::MatrixXd pressure_BC4(label NPmodes, label NUmodes);
 
        //--------------------------------------------------------------------------
        List< Eigen::MatrixXd > bcVelocityVec(label NUmodes, label NSUPmodes);
 
        //--------------------------------------------------------------------------
        List< Eigen::MatrixXd > bcVelocityMat(label NUmodes, label NSUPmodes);
 
        //--------------------------------------------------------------------------
        //  Projection Methods Momentum Equation
        Eigen::MatrixXd diffusive_term_flux_method(label NUmodes, label NPmodes,
                label NSUPmodes);
 
        //--------------------------------------------------------------------------
        List<Eigen::MatrixXd> boundary_vector_diffusion(label NUmodes, label NPmodes,
                label NSUPmodes);
 
        //--------------------------------------------------------------------------
        List<Eigen::MatrixXd> boundary_vector_convection(label NUmodes, label NPmodes,
                label NSUPmodes);
 
        //--------------------------------------------------------------------------
        Eigen::Tensor<double, 3 >  convective_term_flux_tens(label NUmodes,
                label NPmodes,
                label NSUPmodes);
 
        //--------------------------------------------------------------------------
        List < Eigen::MatrixXd > pressure_gradient_term_linsys_div(label NPmodes);
 
        //--------------------------------------------------------------------------
        List < Eigen::MatrixXd > pressure_gradient_term_linsys_diff(label NPmodes);
 
        //--------------------------------------------------------------------------
        List < Eigen::MatrixXd > pressure_gradient_term_linsys_conv(label NPmodes);
 
        //--------------------------------------------------------------------------
        //  Projection Methods Flux Equation
        Eigen::MatrixXd diffusive_term_consistent(label NUmodes, label NPmodes,
                label NSUPmodes);
 
        //--------------------------------------------------------------------------
        List < Eigen::MatrixXd > boundary_vector_diffusion_consistent(label NUmodes,
                label NSUPmodes);
 
        //--------------------------------------------------------------------------
        List < Eigen::MatrixXd > boundary_vector_convection_consistent(label NUmodes,
                label NSUPmodes);
 
        //--------------------------------------------------------------------------
        Eigen::MatrixXd mass_matrix_newtime_consistent(label NUmodes, label NPmodes,
                label NSUPmodes);
 
        //--------------------------------------------------------------------------
        Eigen::MatrixXd mass_matrix_oldtime_consistent(label NUmodes, label NPmodes,
                label NSUPmodes);
 
        //--------------------------------------------------------------------------
        Eigen::MatrixXd pressure_gradient_term_consistent(label NUmodes, label NPmodes,
                label NSUPmodes);
        //--------------------------------------------------------------------------
        Eigen::Tensor<double, 3 >  convective_term_consistent_tens(label NUmodes,
                label NPmodes,
                label NSUPmodes);
 
        //--------------------------------------------------------------------------
        void change_viscosity(double mu);
 
        //--------------------------------------------------------------------------
        void forcesMatrices(label NUmodes, label NPmodes, label NSUPmodes);
 
        //--------------------------------------------------------------------------
        void forcesMatrices(label nModes);
 
 
        //--------------------------------------------------------------------------
        void reconstructLiftAndDrag(const Eigen::MatrixXd& velCoeffs,
                                    const Eigen::MatrixXd& pressureCoeffs, fileName folder);
 
        //--------------------------------------------------------------------------
        Eigen::Tensor<double, 3 > convective_term_tens(label NUmodes,
                label NPmodes,
                label NSUPmodes);
        
        //--------------------------------------------------------------------------
        Eigen::Tensor<double, 3 > convective_term_tens_cache(label NUmodes,
                label NPmodes,
                label NSUPmodes);

        void restart();
 
 
};
 
#endif