From f4427f5bdb40cd1416f65c6417570e301e8cce8a Mon Sep 17 00:00:00 2001 From: bigfooted Date: Thu, 15 Feb 2024 12:27:57 +0100 Subject: [PATCH] remove duplicated setvorticity_strainmag --- SU2_CFD/include/variables/CIncNSVariable.hpp | 5 - SU2_CFD/include/variables/CNSVariable.hpp | 5 - .../include/variables/CPBIncNSVariable.hpp | 21 +- SU2_CFD/include/variables/CVariable.hpp | 29 +- SU2_CFD/src/numerics/heat.cpp | 322 ------------------ SU2_CFD/src/solvers/CPBIncNSSolver.cpp | 15 +- SU2_CFD/src/variables/CIncNSVariable.cpp | 52 --- SU2_CFD/src/variables/CNSVariable.cpp | 53 --- SU2_CFD/src/variables/CPBIncNSVariable.cpp | 61 +--- 9 files changed, 35 insertions(+), 528 deletions(-) delete mode 100644 SU2_CFD/src/numerics/heat.cpp diff --git a/SU2_CFD/include/variables/CIncNSVariable.hpp b/SU2_CFD/include/variables/CIncNSVariable.hpp index 48aa8ea0d3e..a5b42a512cc 100644 --- a/SU2_CFD/include/variables/CIncNSVariable.hpp +++ b/SU2_CFD/include/variables/CIncNSVariable.hpp @@ -62,11 +62,6 @@ class CIncNSVariable final : public CIncEulerVariable { Primitive(iPoint, indices.LaminarViscosity()) = laminarViscosity; } - /*! - * \brief Set the vorticity value. - */ - bool SetVorticity_StrainMag() override; - /*! * \overload * \param[in] eddy_visc - Value of the eddy viscosity. diff --git a/SU2_CFD/include/variables/CNSVariable.hpp b/SU2_CFD/include/variables/CNSVariable.hpp index f55b728c322..d23b7a16eb4 100644 --- a/SU2_CFD/include/variables/CNSVariable.hpp +++ b/SU2_CFD/include/variables/CNSVariable.hpp @@ -79,11 +79,6 @@ class CNSVariable final : public CEulerVariable { Primitive(iPoint, indices.CpTotal()) = val_Cp; } - /*! - * \brief Set the vorticity value. - */ - bool SetVorticity_StrainMag() override; - /*! * \overload * \param[in] eddy_visc - Value of the eddy viscosity. diff --git a/SU2_CFD/include/variables/CPBIncNSVariable.hpp b/SU2_CFD/include/variables/CPBIncNSVariable.hpp index 4c42bc732cf..26c5d849513 100644 --- a/SU2_CFD/include/variables/CPBIncNSVariable.hpp +++ b/SU2_CFD/include/variables/CPBIncNSVariable.hpp @@ -50,12 +50,12 @@ class CPBIncNSVariable final : public CPBIncEulerVariable { VectorType DES_LengthScale; public: - + /*! * \brief Constructor of the class. */ CPBIncNSVariable(void); - + /*! * \param[in] val_pressure - value of the pressure. * \param[in] val_velocity - Value of the flow velocity (initialization value). @@ -64,20 +64,20 @@ class CPBIncNSVariable final : public CPBIncEulerVariable { * \param[in] val_nvar - Number of variables of the problem. * \param[in] config - Definition of the particular problem. */ - CPBIncNSVariable(su2double val_density, su2double val_pressure, su2double *val_velocity, + CPBIncNSVariable(su2double val_density, su2double val_pressure, su2double *val_velocity, unsigned long npoint, unsigned short val_nDim, unsigned short val_nvar, CConfig *config); - + /*! * \brief Destructor of the class. */ virtual ~CPBIncNSVariable() = default; - + /*! * \brief Set all the primitive variables for incompressible flows */ bool SetPrimVar(unsigned long iPoint, su2double Density_Inf, su2double Viscosity_Inf, su2double eddy_visc, su2double turb_ke, CConfig *config); using CVariable::SetPrimVar; - + /*! * \brief Set the laminar viscosity. */ @@ -85,11 +85,6 @@ class CPBIncNSVariable final : public CPBIncEulerVariable { Primitive(iPoint,nDim+2) = laminarViscosity; } - /*! - * \brief Set the vorticity value. - */ - bool SetVorticity_StrainMag() override; - /*! * \overload * \param[in] eddy_visc - Value of the eddy viscosity. @@ -134,11 +129,11 @@ class CPBIncNSVariable final : public CPBIncEulerVariable { * \return Value of the DES length Scale. */ inline su2double GetDES_LengthScale(unsigned long iPoint) const override { return DES_LengthScale(iPoint); } - + /*! * \brief Get the thermal conductivity of the flow. * \return Value of the laminar viscosity of the flow. */ inline su2double GetThermalConductivity(unsigned long iPoint) const override { return 0.0; } - + }; diff --git a/SU2_CFD/include/variables/CVariable.hpp b/SU2_CFD/include/variables/CVariable.hpp index 5bc19e6d871..ceaca8d1cf8 100644 --- a/SU2_CFD/include/variables/CVariable.hpp +++ b/SU2_CFD/include/variables/CVariable.hpp @@ -1615,11 +1615,6 @@ class CVariable { */ inline virtual void SetSpecificHeatCv(unsigned long iPoint, su2double Cv) {} - /*! - * \brief A virtual member. - */ - inline virtual bool SetVorticity_StrainMag() { return false; } - /*! * \brief A virtual member. */ @@ -2438,32 +2433,32 @@ class CVariable { inline virtual su2double GetMassFlux(unsigned long iPoint) const { return 0.0; } inline virtual void SetMassFlux(unsigned long iPoint, su2double val_MassFlux) {} - + inline virtual void AddMassFlux(unsigned long iPoint, su2double val_MassFlux) {} - + inline virtual void SubtractMassFlux(unsigned long iPoint, su2double val_MassFlux) {} - + inline virtual void Set_Mom_Coeff(unsigned long iPoint, su2double *val_Mom_Coeff) {} - + inline virtual su2double Get_Mom_Coeff(unsigned long iPoint, unsigned short val_Var) { return 0.0; } - + inline virtual su2double Get_Mom_Coeff_nb(unsigned long iPoint, unsigned short val_Var) { return 0.0; } - + inline virtual void Set_Mom_Coeff_nb(unsigned long iPoint, su2double *val_Mom_Coeff) {} - + inline virtual void Add_Mom_Coeff_nb(unsigned long iPoint, su2double val_coeff_nb, unsigned short val_Var) {} inline virtual void Add_Mom_Coeff(unsigned long iPoint, su2double val_coeff, unsigned short val_Var) {} - + inline virtual void Set_Mom_Coeff_nbZero(unsigned long iPoint) {} - + inline virtual void Set_Mom_CoeffZero(unsigned long iPoint) {} - + inline virtual void Set_Mom_Coeff(unsigned long iPoint, unsigned short val_Var, su2double val_Mom_Coeff) {} - + inline virtual void Set_Mom_Coeff_nb(unsigned long iPoint, unsigned short val_Var, su2double val_Mom_Coeff) {} inline virtual su2double GetPoisson_Coeff(unsigned long iPoint) { return 0.0; } - + inline virtual void SetPoisson_Coeff(unsigned long iPoint, su2double val_Poisson_Coeff) {} }; diff --git a/SU2_CFD/src/numerics/heat.cpp b/SU2_CFD/src/numerics/heat.cpp deleted file mode 100644 index d9e44f54ee5..00000000000 --- a/SU2_CFD/src/numerics/heat.cpp +++ /dev/null @@ -1,322 +0,0 @@ -/*! - * \file heat.cpp - * \brief Implementation of numerics classes for heat transfer. - * \author F. Palacios, T. Economon - * \version 7.0.7 "Blackbird" - * - * SU2 Project Website: https://su2code.github.io - * - * The SU2 Project is maintained by the SU2 Foundation - * (http://su2foundation.org) - * - * Copyright 2012-2020, SU2 Contributors (cf. AUTHORS.md) - * - * SU2 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 2.1 of the License, or (at your option) any later version. - * - * SU2 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 SU2. If not, see . - */ - -#include "../../include/numerics/heat.hpp" - -CCentSca_Heat::CCentSca_Heat(unsigned short val_nDim, unsigned short val_nVar, CConfig *config) : - CNumerics(val_nDim, val_nVar, config) { - - implicit = (config->GetKind_TimeIntScheme_Turb() == EULER_IMPLICIT); - /* A grid is defined as dynamic if there's rigid grid movement or grid deformation AND the problem is time domain */ - dynamic_grid = config->GetDynamic_Grid(); - - MeanVelocity = new su2double [nDim]; - - Laminar_Viscosity_i = config->GetViscosity_FreeStreamND(); - Laminar_Viscosity_j = config->GetViscosity_FreeStreamND(); - - Param_Kappa_4 = config->GetKappa_4th_Heat(); - Diff_Lapl = new su2double [nVar]; -} - -CCentSca_Heat::~CCentSca_Heat(void) { - - delete [] MeanVelocity; - delete [] Diff_Lapl; - -} - -void CCentSca_Heat::ComputeResidual(su2double *val_residual, su2double **val_Jacobian_i, - su2double **val_Jacobian_j, CConfig *config) { - - AD::StartPreacc(); - AD::SetPreaccIn(V_i, nDim+3); AD::SetPreaccIn(V_j, nDim+3); - AD::SetPreaccIn(Temp_i); AD::SetPreaccIn(Temp_j); - AD::SetPreaccIn(Und_Lapl_i, nVar); AD::SetPreaccIn(Und_Lapl_j, nVar); - AD::SetPreaccIn(Normal, nDim); - if (dynamic_grid) { - AD::SetPreaccIn(GridVel_i, nDim); AD::SetPreaccIn(GridVel_j, nDim); - } - - /*--- Primitive variables at point i and j ---*/ - - Pressure_i = V_i[0]; Pressure_j = V_j[0]; - DensityInc_i = V_i[nDim+2]; DensityInc_j = V_j[nDim+2]; - BetaInc2_i = V_i[nDim+3]; BetaInc2_j = V_j[nDim+3]; - - /*--- Projected velocities at the current edge ---*/ - - ProjVelocity = 0.0; ProjVelocity_i = 0.0; ProjVelocity_j = 0.0; Area = 0.0; - for (iDim = 0; iDim < nDim; iDim++) { - ProjVelocity_i += V_i[iDim+1]*Normal[iDim]; - ProjVelocity_j += V_j[iDim+1]*Normal[iDim]; - Area += Normal[iDim]*Normal[iDim]; - } - Area = sqrt(Area); - - /*--- Computing the second order centered scheme part ---*/ - - ProjVelocity = 0.5*(ProjVelocity_i+ProjVelocity_j); - - val_residual[0] = 0.5*(Temp_i + Temp_j)*ProjVelocity; - - if (implicit) { - val_Jacobian_i[0][0] = 0.5*ProjVelocity; - val_Jacobian_j[0][0] = 0.5*ProjVelocity; - } - - /*--- Adding artificial dissipation to stabilize the centered scheme ---*/ - - Diff_Lapl[0] = Und_Lapl_i[0]-Und_Lapl_j[0]; - - SoundSpeed_i = sqrt(ProjVelocity_i*ProjVelocity_i + (BetaInc2_i/DensityInc_i)*Area*Area); - SoundSpeed_j = sqrt(ProjVelocity_j*ProjVelocity_j + (BetaInc2_j/DensityInc_j)*Area*Area); - - Local_Lambda_i = fabs(ProjVelocity_i)+SoundSpeed_i; - Local_Lambda_j = fabs(ProjVelocity_j)+SoundSpeed_j; - MeanLambda = 0.5*(Local_Lambda_i+Local_Lambda_j); - - val_residual[0] += -Param_Kappa_4*Diff_Lapl[0]*MeanLambda; - - if (implicit) { - cte_0 = Param_Kappa_4*su2double(Neighbor_i+1)*MeanLambda; - cte_1 = Param_Kappa_4*su2double(Neighbor_j+1)*MeanLambda; - - val_Jacobian_i[0][0] += cte_0; - val_Jacobian_j[0][0] -= cte_1; - } - - AD::SetPreaccOut(val_residual[0]); - AD::EndPreacc(); - -} - -CUpwSca_Heat::CUpwSca_Heat(unsigned short val_nDim, unsigned short val_nVar, CConfig *config) : - CNumerics(val_nDim, val_nVar, config) { - - implicit = (config->GetKind_TimeIntScheme_Turb() == EULER_IMPLICIT); - /* A grid is defined as dynamic if there's rigid grid movement or grid deformation AND the problem is time domain */ - dynamic_grid = config->GetDynamic_Grid(); - - Velocity_i = new su2double [nDim]; - Velocity_j = new su2double [nDim]; - - Laminar_Viscosity_i = config->GetViscosity_FreeStreamND(); - Laminar_Viscosity_j = config->GetViscosity_FreeStreamND(); -} - -CUpwSca_Heat::~CUpwSca_Heat(void) { - - delete [] Velocity_i; - delete [] Velocity_j; - -} - -void CUpwSca_Heat::ComputeResidual(su2double *val_residual, su2double **val_Jacobian_i, - su2double **val_Jacobian_j, CConfig *config) { - - q_ij = 0.0; - - AD::StartPreacc(); - AD::SetPreaccIn(V_i, nDim+1); AD::SetPreaccIn(V_j, nDim+1); - AD::SetPreaccIn(Temp_i); AD::SetPreaccIn(Temp_j); - AD::SetPreaccIn(Normal, nDim); - if (dynamic_grid) { - AD::SetPreaccIn(GridVel_i, nDim); AD::SetPreaccIn(GridVel_j, nDim); - } - - for (iDim = 0; iDim < nDim; iDim++) { - Velocity_i[iDim] = V_i[iDim+1]; - Velocity_j[iDim] = V_j[iDim+1]; - q_ij += 0.5*(Velocity_i[iDim]+Velocity_j[iDim])*Normal[iDim]; - } - - a0 = 0.5*(q_ij+fabs(q_ij)); - a1 = 0.5*(q_ij-fabs(q_ij)); - val_residual[0] = a0*Temp_i+a1*Temp_j; - - if (implicit) { - val_Jacobian_i[0][0] = a0; - val_Jacobian_j[0][0] = a1; - } - - AD::SetPreaccOut(val_residual[0]); - AD::EndPreacc(); - -} - -CAvgGrad_Heat::CAvgGrad_Heat(unsigned short val_nDim, unsigned short val_nVar, CConfig *config) : - CNumerics(val_nDim, val_nVar, config) { - - implicit = (config->GetKind_TimeIntScheme_Heat() == EULER_IMPLICIT); - - Edge_Vector = new su2double [nDim]; - Proj_Mean_GradHeatVar_Normal = new su2double [nVar]; - Proj_Mean_GradHeatVar_Corrected = new su2double [nVar]; - Mean_GradHeatVar = new su2double* [nVar]; - for (iVar = 0; iVar < nVar; iVar++) - Mean_GradHeatVar[iVar] = new su2double [nDim]; - -} - -CAvgGrad_Heat::~CAvgGrad_Heat(void) { - - delete [] Edge_Vector; - delete [] Proj_Mean_GradHeatVar_Normal; - delete [] Proj_Mean_GradHeatVar_Corrected; - for (iVar = 0; iVar < nVar; iVar++) - delete [] Mean_GradHeatVar[iVar]; - delete [] Mean_GradHeatVar; - -} - -void CAvgGrad_Heat::ComputeResidual(su2double *val_residual, su2double **Jacobian_i, - su2double **Jacobian_j, CConfig *config) { - - AD::StartPreacc(); - AD::SetPreaccIn(Coord_i, nDim); AD::SetPreaccIn(Coord_j, nDim); - AD::SetPreaccIn(Normal, nDim); - AD::SetPreaccIn(Temp_i); AD::SetPreaccIn(Temp_j); - AD::SetPreaccIn(ConsVar_Grad_i[0],nDim); AD::SetPreaccIn(ConsVar_Grad_j[0],nDim); - AD::SetPreaccIn(Thermal_Diffusivity_i); AD::SetPreaccIn(Thermal_Conductivity_j); - - Thermal_Diffusivity_Mean = 0.5*(Thermal_Diffusivity_i + Thermal_Diffusivity_j); - - /*--- Compute vector going from iPoint to jPoint ---*/ - - dist_ij_2 = 0; proj_vector_ij = 0; - for (iDim = 0; iDim < nDim; iDim++) { - Edge_Vector[iDim] = Coord_j[iDim]-Coord_i[iDim]; - dist_ij_2 += Edge_Vector[iDim]*Edge_Vector[iDim]; - proj_vector_ij += Edge_Vector[iDim]*Normal[iDim]; - } - if (dist_ij_2 == 0.0) proj_vector_ij = 0.0; - else proj_vector_ij = proj_vector_ij/dist_ij_2; - - /*--- Mean gradient approximation. Projection of the mean gradient in the direction of the edge ---*/ - for (iVar = 0; iVar < nVar; iVar++) { - Proj_Mean_GradHeatVar_Normal[iVar] = 0.0; - Proj_Mean_GradHeatVar_Corrected[iVar] = 0.0; - for (iDim = 0; iDim < nDim; iDim++) { - Mean_GradHeatVar[iVar][iDim] = 0.5*(ConsVar_Grad_i[iVar][iDim] + ConsVar_Grad_j[iVar][iDim]); - Proj_Mean_GradHeatVar_Normal[iVar] += Mean_GradHeatVar[iVar][iDim]*Normal[iDim]; - } - Proj_Mean_GradHeatVar_Corrected[iVar] = Proj_Mean_GradHeatVar_Normal[iVar]; - } - - val_residual[0] = Thermal_Diffusivity_Mean*Proj_Mean_GradHeatVar_Corrected[0]; - - /*--- For Jacobians -> Use of TSL approx. to compute derivatives of the gradients ---*/ - if (implicit) { - Jacobian_i[0][0] = -Thermal_Diffusivity_Mean*proj_vector_ij; - Jacobian_j[0][0] = Thermal_Diffusivity_Mean*proj_vector_ij; - } - - AD::SetPreaccOut(val_residual, nVar); - AD::EndPreacc(); - -} - -CAvgGradCorrected_Heat::CAvgGradCorrected_Heat(unsigned short val_nDim, unsigned short val_nVar, CConfig *config) : - CNumerics(val_nDim, val_nVar, config) { - - implicit = (config->GetKind_TimeIntScheme_Heat() == EULER_IMPLICIT); - - Edge_Vector = new su2double [nDim]; - Proj_Mean_GradHeatVar_Edge = new su2double [nVar]; - Proj_Mean_GradHeatVar_Kappa = new su2double [nVar]; - Proj_Mean_GradHeatVar_Corrected = new su2double [nVar]; - Mean_GradHeatVar = new su2double* [nVar]; - for (iVar = 0; iVar < nVar; iVar++) - Mean_GradHeatVar[iVar] = new su2double [nDim]; - -} - -CAvgGradCorrected_Heat::~CAvgGradCorrected_Heat(void) { - - delete [] Edge_Vector; - delete [] Proj_Mean_GradHeatVar_Edge; - delete [] Proj_Mean_GradHeatVar_Kappa; - delete [] Proj_Mean_GradHeatVar_Corrected; - for (iVar = 0; iVar < nVar; iVar++) - delete [] Mean_GradHeatVar[iVar]; - delete [] Mean_GradHeatVar; - -} - -void CAvgGradCorrected_Heat::ComputeResidual(su2double *val_residual, su2double **Jacobian_i, - su2double **Jacobian_j, CConfig *config) { - - AD::StartPreacc(); - AD::SetPreaccIn(Coord_i, nDim); AD::SetPreaccIn(Coord_j, nDim); - AD::SetPreaccIn(Normal, nDim); - AD::SetPreaccIn(Temp_i); AD::SetPreaccIn(Temp_j); - AD::SetPreaccIn(ConsVar_Grad_i[0],nDim); AD::SetPreaccIn(ConsVar_Grad_j[0],nDim); - AD::SetPreaccIn(Thermal_Diffusivity_i); AD::SetPreaccIn(Thermal_Diffusivity_j); - - Thermal_Diffusivity_Mean = 0.5*(Thermal_Diffusivity_i + Thermal_Diffusivity_j); - - /*--- Compute vector going from iPoint to jPoint ---*/ - - dist_ij_2 = 0; proj_vector_ij = 0; - for (iDim = 0; iDim < nDim; iDim++) { - Edge_Vector[iDim] = Coord_j[iDim]-Coord_i[iDim]; - dist_ij_2 += Edge_Vector[iDim]*Edge_Vector[iDim]; - proj_vector_ij += Edge_Vector[iDim]*Normal[iDim]; - } - if (dist_ij_2 == 0.0) proj_vector_ij = 0.0; - else proj_vector_ij = proj_vector_ij/dist_ij_2; - - /*--- Mean gradient approximation. Projection of the mean gradient - in the direction of the edge ---*/ - - for (iVar = 0; iVar < nVar; iVar++) { - Proj_Mean_GradHeatVar_Edge[iVar] = 0.0; - Proj_Mean_GradHeatVar_Kappa[iVar] = 0.0; - - for (iDim = 0; iDim < nDim; iDim++) { - Mean_GradHeatVar[iVar][iDim] = 0.5*(ConsVar_Grad_i[iVar][iDim] + ConsVar_Grad_j[iVar][iDim]); - Proj_Mean_GradHeatVar_Kappa[iVar] += Mean_GradHeatVar[iVar][iDim]*Normal[iDim]; - Proj_Mean_GradHeatVar_Edge[iVar] += Mean_GradHeatVar[iVar][iDim]*Edge_Vector[iDim]; - } - Proj_Mean_GradHeatVar_Corrected[iVar] = Proj_Mean_GradHeatVar_Kappa[iVar]; - Proj_Mean_GradHeatVar_Corrected[iVar] -= Proj_Mean_GradHeatVar_Edge[iVar]*proj_vector_ij - - (Temp_j-Temp_i)*proj_vector_ij; - } - - val_residual[0] = Thermal_Diffusivity_Mean*Proj_Mean_GradHeatVar_Corrected[0]; - - /*--- For Jacobians -> Use of TSL approx. to compute derivatives of the gradients ---*/ - - if (implicit) { - Jacobian_i[0][0] = -Thermal_Diffusivity_Mean*proj_vector_ij; - Jacobian_j[0][0] = Thermal_Diffusivity_Mean*proj_vector_ij; - } - - AD::SetPreaccOut(val_residual, nVar); - AD::EndPreacc(); -} diff --git a/SU2_CFD/src/solvers/CPBIncNSSolver.cpp b/SU2_CFD/src/solvers/CPBIncNSSolver.cpp index 01d7a216261..e1d193cb702 100644 --- a/SU2_CFD/src/solvers/CPBIncNSSolver.cpp +++ b/SU2_CFD/src/solvers/CPBIncNSSolver.cpp @@ -33,7 +33,7 @@ CPBIncNSSolver::CPBIncNSSolver(CGeometry *geometry, CConfig *config, unsigned short iMesh) : CPBIncEulerSolver(geometry, config, iMesh, true) { - + Viscosity_Inf = config->GetViscosity_FreeStreamND(); Tke_Inf = config->GetTke_FreeStreamND(); @@ -41,7 +41,10 @@ CPBIncNSSolver::CPBIncNSSolver(CGeometry *geometry, CConfig *config, unsigned sh void CPBIncNSSolver::Preprocessing(CGeometry *geometry, CSolver **solver_container, CConfig *config, unsigned short iMesh, unsigned short iRKStep, unsigned short RunTime_EqSystem, bool Output) { -unsigned long iPoint, ErrorCounter = 0; + auto* flowNodes = su2staticcast_p(solver_container[FLOW_SOL]->GetNodes()); + + + unsigned long iPoint, ErrorCounter = 0; su2double StrainMag = 0.0, Omega = 0.0, *Vorticity; unsigned long InnerIter = config->GetInnerIter(); @@ -82,13 +85,12 @@ unsigned long iPoint, ErrorCounter = 0; && !Output && !van_albada) { SetPrimitive_Limiter(geometry, config); } /*--- Evaluate the vorticity and strain rate magnitude ---*/ - - solver_container[FLOW_SOL]->GetNodes()->SetVorticity_StrainMag(); // TODO:PBFlow + ComputeVorticityAndStrainMag(*config, geometry, iMesh); StrainMag_Max = 0.0; Omega_Max = 0.0; for (iPoint = 0; iPoint < nPoint; iPoint++) { - StrainMag = solver_container[FLOW_SOL]->GetNodes()->GetStrainMag(iPoint); + StrainMag = flowNodes->GetStrainMag(iPoint); Vorticity = solver_container[FLOW_SOL]->GetNodes()->GetVorticity(iPoint); Omega = sqrt(Vorticity[0]*Vorticity[0]+ Vorticity[1]*Vorticity[1]+ Vorticity[2]*Vorticity[2]); @@ -189,8 +191,7 @@ unsigned long iPoint, ErrorCounter = 0; } /*--- Evaluate the vorticity and strain rate magnitude ---*/ - - solver_container[FLOW_SOL]->GetNodes()->SetVorticity_StrainMag(); + ComputeVorticityAndStrainMag(*config, geometry, iMesh); } diff --git a/SU2_CFD/src/variables/CIncNSVariable.cpp b/SU2_CFD/src/variables/CIncNSVariable.cpp index 7c8f820ef29..e32b8d6dd73 100644 --- a/SU2_CFD/src/variables/CIncNSVariable.cpp +++ b/SU2_CFD/src/variables/CIncNSVariable.cpp @@ -126,55 +126,3 @@ bool CIncNSVariable::SetPrimVar(unsigned long iPoint, su2double eddy_visc, su2do return physical; } - -bool CIncNSVariable::SetVorticity_StrainMag() { - - for (unsigned long iPoint = 0; iPoint < nPoint; ++iPoint) { - - /*--- Vorticity ---*/ - - Vorticity(iPoint,0) = 0.0; Vorticity(iPoint,1) = 0.0; - - Vorticity(iPoint,2) = Gradient_Primitive(iPoint,2,0)-Gradient_Primitive(iPoint,1,1); - - if (nDim == 3) { - Vorticity(iPoint,0) = Gradient_Primitive(iPoint,3,1)-Gradient_Primitive(iPoint,2,2); - Vorticity(iPoint,1) = -(Gradient_Primitive(iPoint,3,0)-Gradient_Primitive(iPoint,1,2)); - } - - /*--- Strain Magnitude ---*/ - - AD::StartPreacc(); - AD::SetPreaccIn(Gradient_Primitive[iPoint], nDim+1, nDim); - - su2double Div = 0.0; - for (unsigned long iDim = 0; iDim < nDim; iDim++) - Div += Gradient_Primitive(iPoint,iDim+1,iDim); - - StrainMag(iPoint) = 0.0; - - /*--- Add diagonal part ---*/ - - for (unsigned long iDim = 0; iDim < nDim; iDim++) { - StrainMag(iPoint) += pow(Gradient_Primitive(iPoint,iDim+1,iDim) - 1.0/3.0*Div, 2.0); - } - if (nDim == 2) { - StrainMag(iPoint) += pow(1.0/3.0*Div, 2.0); - } - - /*--- Add off diagonals ---*/ - - StrainMag(iPoint) += 2.0*pow(0.5*(Gradient_Primitive(iPoint,1,1) + Gradient_Primitive(iPoint,2,0)), 2); - - if (nDim == 3) { - StrainMag(iPoint) += 2.0*pow(0.5*(Gradient_Primitive(iPoint,1,2) + Gradient_Primitive(iPoint,3,0)), 2); - StrainMag(iPoint) += 2.0*pow(0.5*(Gradient_Primitive(iPoint,2,2) + Gradient_Primitive(iPoint,3,1)), 2); - } - - StrainMag(iPoint) = sqrt(2.0*StrainMag(iPoint)); - - AD::SetPreaccOut(StrainMag(iPoint)); - AD::EndPreacc(); - } - return false; -} \ No newline at end of file diff --git a/SU2_CFD/src/variables/CNSVariable.cpp b/SU2_CFD/src/variables/CNSVariable.cpp index dfce3ac29fb..b2fe70b35cc 100644 --- a/SU2_CFD/src/variables/CNSVariable.cpp +++ b/SU2_CFD/src/variables/CNSVariable.cpp @@ -45,59 +45,6 @@ CNSVariable::CNSVariable(su2double density, const su2double *velocity, su2double } -bool CNSVariable::SetVorticity_StrainMag() { - - SU2_OMP_FOR_STAT(256) - for (unsigned long iPoint = 0; iPoint < nPoint; ++iPoint) { - - /*--- Vorticity ---*/ - - Vorticity(iPoint,0) = 0.0; Vorticity(iPoint,1) = 0.0; - - Vorticity(iPoint,2) = Gradient_Primitive(iPoint,2,0)-Gradient_Primitive(iPoint,1,1); - - if (nDim == 3) { - Vorticity(iPoint,0) = Gradient_Primitive(iPoint,3,1)-Gradient_Primitive(iPoint,2,2); - Vorticity(iPoint,1) = -(Gradient_Primitive(iPoint,3,0)-Gradient_Primitive(iPoint,1,2)); - } - - /*--- Strain Magnitude ---*/ - - AD::StartPreacc(); - AD::SetPreaccIn(Gradient_Primitive[iPoint], nDim+1, nDim); - - su2double Div = 0.0; - for (unsigned long iDim = 0; iDim < nDim; iDim++) - Div += Gradient_Primitive(iPoint,iDim+1,iDim); - - StrainMag(iPoint) = 0.0; - - /*--- Add diagonal part ---*/ - - for (unsigned long iDim = 0; iDim < nDim; iDim++) { - StrainMag(iPoint) += pow(Gradient_Primitive(iPoint,iDim+1,iDim) - 1.0/3.0*Div, 2.0); - } - if (nDim == 2) { - StrainMag(iPoint) += pow(1.0/3.0*Div, 2.0); - } - - /*--- Add off diagonals ---*/ - - StrainMag(iPoint) += 2.0*pow(0.5*(Gradient_Primitive(iPoint,1,1) + Gradient_Primitive(iPoint,2,0)), 2); - - if (nDim == 3) { - StrainMag(iPoint) += 2.0*pow(0.5*(Gradient_Primitive(iPoint,1,2) + Gradient_Primitive(iPoint,3,0)), 2); - StrainMag(iPoint) += 2.0*pow(0.5*(Gradient_Primitive(iPoint,2,2) + Gradient_Primitive(iPoint,3,1)), 2); - } - - StrainMag(iPoint) = sqrt(2.0*StrainMag(iPoint)); - - AD::SetPreaccOut(StrainMag(iPoint)); - AD::EndPreacc(); - } - return false; -} - void CNSVariable::SetRoe_Dissipation_NTS(unsigned long iPoint, su2double val_delta, su2double val_const_DES){ diff --git a/SU2_CFD/src/variables/CPBIncNSVariable.cpp b/SU2_CFD/src/variables/CPBIncNSVariable.cpp index bd5587416d1..1480cae012b 100644 --- a/SU2_CFD/src/variables/CPBIncNSVariable.cpp +++ b/SU2_CFD/src/variables/CPBIncNSVariable.cpp @@ -40,7 +40,7 @@ CPBIncNSVariable::CPBIncNSVariable(su2double val_density, su2double val_pressure, su2double *val_velocity, unsigned long nPoint, unsigned short nDim, unsigned short nvar, CConfig *config) : CPBIncEulerVariable(val_density, val_pressure, val_velocity, nPoint, nDim, nvar, config) { - + Vorticity.resize(nPoint,3); StrainMag.resize(nPoint); DES_LengthScale.resize(nPoint) = su2double(0.0); @@ -48,71 +48,24 @@ CPBIncNSVariable::CPBIncNSVariable(su2double val_density, su2double val_pressure } bool CPBIncNSVariable::SetPrimVar(unsigned long iPoint, su2double Density_Inf, su2double Viscosity_Inf, su2double eddy_visc, su2double turb_ke, CConfig *config) { - + unsigned short iVar; bool check_dens = false, physical = true; - + /*--- Set eddy viscosity locally and in the fluid model. ---*/ SetEddyViscosity(iPoint, eddy_visc); - + /*--- Set the value of the density ---*/ - + check_dens = SetDensity(iPoint, Density_Inf); /*--- Set the value of the velocity and velocity^2 (requires density) ---*/ SetVelocity(iPoint); - + /*--- Set laminar viscosity ---*/ - + SetLaminarViscosity(iPoint, Viscosity_Inf); return physical; } - - -bool CPBIncNSVariable::SetVorticity_StrainMag(void) { - - for (unsigned long iPoint = 0; iPoint < nPoint; ++iPoint) { - - /*--- Vorticity ---*/ - - Vorticity(iPoint,0) = 0.0; Vorticity(iPoint,1) = 0.0; - - Vorticity(iPoint,2) = Gradient_Primitive(iPoint,2,0)-Gradient_Primitive(iPoint,1,1); - - if (nDim == 3) { - Vorticity(iPoint,0) = Gradient_Primitive(iPoint,3,1)-Gradient_Primitive(iPoint,2,2); - Vorticity(iPoint,1) = -(Gradient_Primitive(iPoint,3,0)-Gradient_Primitive(iPoint,1,2)); - } - - /*--- Strain Magnitude ---*/ - su2double Div = 0.0; - for (unsigned long iDim = 0; iDim < nDim; iDim++) - Div += Gradient_Primitive(iPoint,iDim+1,iDim); - - StrainMag(iPoint) = 0.0; - - /*--- Add diagonal part ---*/ - - for (unsigned long iDim = 0; iDim < nDim; iDim++) { - StrainMag(iPoint) += pow(Gradient_Primitive(iPoint,iDim+1,iDim) - 1.0/3.0*Div, 2.0); - } - if (nDim == 2) { - StrainMag(iPoint) += pow(1.0/3.0*Div, 2.0); - } - - /*--- Add off diagonals ---*/ - - StrainMag(iPoint) += 2.0*pow(0.5*(Gradient_Primitive(iPoint,1,1) + Gradient_Primitive(iPoint,2,0)), 2); - - if (nDim == 3) { - StrainMag(iPoint) += 2.0*pow(0.5*(Gradient_Primitive(iPoint,1,2) + Gradient_Primitive(iPoint,3,0)), 2); - StrainMag(iPoint) += 2.0*pow(0.5*(Gradient_Primitive(iPoint,2,2) + Gradient_Primitive(iPoint,3,1)), 2); - } - - StrainMag(iPoint) = sqrt(2.0*StrainMag(iPoint)); - - } - return false; -}