From cc64d25add4ea06a43685fb49bb6492b78902d9c Mon Sep 17 00:00:00 2001 From: koodlyakshay Date: Sun, 19 Jan 2020 23:04:57 +0100 Subject: [PATCH 01/29] Initial commit with roughness corrections for SA and SST models. --- Common/include/config_structure.hpp | 34 +++++- Common/include/dual_grid_structure.hpp | 12 ++ Common/include/dual_grid_structure.inl | 4 + Common/include/option_structure.hpp | 11 ++ Common/src/config_structure.cpp | 101 +++++++++++++++- Common/src/geometry/CPhysicalGeometry.cpp | 11 +- SU2_CFD/include/numerics_structure.hpp | 11 ++ SU2_CFD/include/numerics_structure.inl | 5 + SU2_CFD/src/numerics_direct_turbulent.cpp | 25 +++- SU2_CFD/src/solver_direct_turbulent.cpp | 141 ++++++++++++++++++---- 10 files changed, 325 insertions(+), 30 deletions(-) diff --git a/Common/include/config_structure.hpp b/Common/include/config_structure.hpp index 8bc4e3438695..9699c026a0fc 100644 --- a/Common/include/config_structure.hpp +++ b/Common/include/config_structure.hpp @@ -314,6 +314,7 @@ class CConfig { su2double *Outlet_Pressure; /*!< \brief Specified back pressures (static) for outlet boundaries. */ su2double *Isothermal_Temperature; /*!< \brief Specified isothermal wall temperatures (static). */ su2double *Heat_Flux; /*!< \brief Specified wall heat fluxes. */ + su2double *Roughness_Height; /*!< \brief Equivalent sand grain roughness for the marker. */ su2double *Displ_Value; /*!< \brief Specified displacement for displacement boundaries. */ su2double *Load_Value; /*!< \brief Specified force for load boundaries. */ su2double *Damper_Constant; /*!< \brief Specified constant for damper boundaries. */ @@ -565,6 +566,8 @@ class CConfig { *Kind_Inc_Outlet, *Kind_Data_Riemann, *Kind_Data_Giles; /*!< \brief Kind of inlet boundary treatment. */ + unsigned short *Kind_Wall; /*!< \brief Type of wall treatment. */ + unsigned short nWall_Types; /*!< \brief Number of wall treatment types listed. */ unsigned short nInc_Inlet; /*!< \brief Number of inlet boundary treatment types listed. */ unsigned short nInc_Outlet; /*!< \brief Number of inlet boundary treatment types listed. */ su2double Inc_Inlet_Damping; /*!< \brief Damping factor applied to the iterative updates to the velocity at a pressure inlet in incompressible flow. */ @@ -933,7 +936,8 @@ class CConfig { nMarkerPitching_Ampl, /*!< \brief Number of values provided for pitching amplitude of marker. */ nMarkerPitching_Phase, /*!< \brief Number of values provided for pitching phase offset of marker. */ nMarkerPlunging_Omega, /*!< \brief Number of values provided for angular frequency of marker. */ - nMarkerPlunging_Ampl; /*!< \brief Number of values provided for plunging amplitude of marker. */ + nMarkerPlunging_Ampl, /*!< \brief Number of values provided for plunging amplitude of marker. */ + nRoughWall; /*!< \brief Number of rough walls. */ su2double *Omega_HB; /*!< \brief Frequency for Harmonic Balance Operator (in rad/s). */ unsigned short nOmega_HB, /*!< \brief Number of frequencies in Harmonic Balance Operator. */ @@ -1034,7 +1038,8 @@ class CConfig { *default_ffd_axis, /*!< \brief Default FFD axis for the COption class. */ *default_inc_crit, /*!< \brief Default incremental criteria array for the COption class. */ *default_extrarelfac, /*!< \brief Default extra relaxation factor for Giles BC in the COption class. */ - *default_sineload_coeff; /*!< \brief Default values for a sine load. */ + *default_sineload_coeff, /*!< \brief Default values for a sine load. */ + *default_roughness; /*!< \brief Default values for wall roughness. */ unsigned short nSpanWiseSections; /*!< \brief number of span-wise sections */ unsigned short nSpanMaxAllZones; /*!< \brief number of maximum span-wise sections for all zones */ @@ -5082,6 +5087,17 @@ class CConfig { */ unsigned short GetKind_ActDisk(void); + /*! + * \brief Get the kind of wall. + * \return Kind of wall - smooth or rough. + */ + unsigned short GetKindWall(string val_marker); + + /*! + * \brief Set the kind of wall - rough or smooth. + */ + void SetKindWall(string val_marker, unsigned short val_kindwall); + /*! * \brief Get the number of sections. * \return Number of sections @@ -6801,6 +6817,20 @@ class CConfig { */ su2double* GetWallFunction_DoubleInfo(string val_marker); + /*! + * \brief Get the wall roughness height on a wall boundary (Heatflux or Isothermal). + * \param[in] val_index - Index corresponding to the boundary. + * \return The wall roughness height. + */ + su2double GetWall_RoughnessHeight(string val_marker); + + /*! + * \brief Set the wall roughness height on a wall boundary (Heatflux or Isothermal). + * \param[in] val_index - Index corresponding to the boundary. + * \return The wall roughness height. + */ + void SetWall_RoughnessHeight(string val_marker, su2double val_roughness); + /*! * \brief Get the target (pressure, massflow, etc) at an engine inflow boundary. * \param[in] val_index - Index corresponding to the engine inflow boundary. diff --git a/Common/include/dual_grid_structure.hpp b/Common/include/dual_grid_structure.hpp index 200b7bd3f693..2df7215969b0 100644 --- a/Common/include/dual_grid_structure.hpp +++ b/Common/include/dual_grid_structure.hpp @@ -162,6 +162,7 @@ class CPoint : public CDualGrid { su2double Wall_Distance; /*!< \brief Distance to the nearest wall. */ su2double SharpEdge_Distance; /*!< \brief Distance to a sharp edge. */ su2double Curvature; /*!< \brief Value of the surface curvature (SU2_GEO). */ + su2double RoughnessHeight; /*!< \brief Roughness of the nearest wall. */ unsigned long GlobalIndex; /*!< \brief Global index in the parallel simulation. */ unsigned short nNeighbor; /*!< \brief Number of neighbors. */ bool Flip_Orientation; /*!< \brief Flip the orientation of the normal. */ @@ -250,6 +251,17 @@ class CPoint : public CDualGrid { * \return Value of the distance to the nearest wall. */ su2double GetSharpEdge_Distance(void); + + /*! + * \brief Set the roughness height of the nearest wall. + */ + void SetRoughnessHeight(su2double val_roughness); + + /*! + * \brief Get the roughness height of the nearest wall. + * \return Value of the roughness at the nearest wall. + */ + su2double GetRoughnessHeight(); /*! * \brief Set the number of elements that compose the control volume. diff --git a/Common/include/dual_grid_structure.inl b/Common/include/dual_grid_structure.inl index 7e4b8b8609a3..c29da55c9346 100644 --- a/Common/include/dual_grid_structure.inl +++ b/Common/include/dual_grid_structure.inl @@ -232,6 +232,10 @@ inline void CPoint::SetCurvature(su2double val_curvature) { Curvature = val_curv inline void CPoint::SetSharpEdge_Distance(su2double val_distance) { SharpEdge_Distance = val_distance; } +inline void CPoint::SetRoughnessHeight(su2double val_roughness) {RoughnessHeight = val_roughness; } + +inline su2double CPoint::GetRoughnessHeight() { return RoughnessHeight; } + inline su2double CPoint::GetWall_Distance(void) { return Wall_Distance; } inline su2double CPoint::GetCurvature(void) { return Curvature; } diff --git a/Common/include/option_structure.hpp b/Common/include/option_structure.hpp index 759b64129317..969ce6a778b7 100644 --- a/Common/include/option_structure.hpp +++ b/Common/include/option_structure.hpp @@ -1289,6 +1289,17 @@ static const map ActDisk_Map = CCreateMap WallType_Map = CCreateMap +("SMOOTH", SMOOTH) +("ROUGH", ROUGH); + /*! * \brief types of geometric entities based on VTK nomenclature */ diff --git a/Common/src/config_structure.cpp b/Common/src/config_structure.cpp index 68007b40b7fa..9482fc0c1856 100644 --- a/Common/src/config_structure.cpp +++ b/Common/src/config_structure.cpp @@ -571,7 +571,7 @@ void CConfig::SetPointersNull(void) { Inlet_Velocity = NULL; Inflow_Mach = NULL; Inflow_Pressure = NULL; Exhaust_Pressure = NULL; Outlet_Pressure = NULL; Isothermal_Temperature= NULL; Heat_Flux = NULL; Displ_Value = NULL; Load_Value = NULL; - FlowLoad_Value = NULL; + FlowLoad_Value = NULL; Roughness_Height = NULL; ElasticityMod = NULL; PoissonRatio = NULL; MaterialDensity = NULL; @@ -825,6 +825,7 @@ void CConfig::SetConfig_Options() { default_sineload_coeff = new su2double[3]; default_nacelle_location = new su2double[5]; default_wrt_freq = new su2double[3]; + default_roughness = new su2double[1]; /*--- All temperature polynomial fits for the fluid models currently assume a quartic form (5 coefficients). For example, @@ -1281,6 +1282,11 @@ void CConfig::SetConfig_Options() { /*!\brief MARKER_HEATFLUX \n DESCRIPTION: Specified heat flux wall boundary marker(s) Format: ( Heat flux marker, wall heat flux (static), ... ) \ingroup Config*/ addStringDoubleListOption("MARKER_HEATFLUX", nMarker_HeatFlux, Marker_HeatFlux, Heat_Flux); + /*!\brief WALL_TYPE \n DESCRIPTION: List of wall types - SMOOTH or ROUGH (All walls are SMOOTH by default). List length must match number of wall markers or shouldn't appear at all. \ingroup Config*/ + addEnumListOption("WALL_TYPE", nWall_Types, Kind_Wall, WallType_Map); + /*!\brief Wall roughness height default value = 0 (smooth). \ingroup Config */ + default_roughness[0] = 0.0; //default_roughness[1] = 0.0; + addDoubleListOption("WALL_ROUGHNESS", nRoughWall, Roughness_Height); /*!\brief MARKER_ENGINE_INFLOW \n DESCRIPTION: Engine inflow boundary marker(s) Format: ( nacelle inflow marker, fan face Mach, ... ) \ingroup Config*/ addStringDoubleListOption("MARKER_ENGINE_INFLOW", nMarker_EngineInflow, Marker_EngineInflow, EngineInflow_Target); @@ -4354,6 +4360,55 @@ void CConfig::SetPostprocessing(unsigned short val_software, unsigned short val_ } } + /*--- Wall roughness handling begins here ---*/ + /*--- Check that if the wall roughness array are compatible and set deafult values if needed. ---*/ + + if ((nMarker_HeatFlux > 0) || (nMarker_Isothermal > 0)) { + /*--- If the config option is not declared, then assume smooth walls and set roughness to zero. ---*/ + if (nWall_Types != nMarker_HeatFlux + nMarker_Isothermal) { + /*--- If this option is not used at all, assume all walls are smooth and set roughness height to zero. ---*/ + if (nWall_Types == 0) { + if (Roughness_Height != NULL) Roughness_Height = NULL; + Roughness_Height = new su2double [nWall_Types]; + if (Kind_Wall != NULL) Kind_Wall = NULL; + Kind_Wall = new unsigned short [nWall_Types]; + for (unsigned short iMarker = 0; iMarker < nMarker_HeatFlux; iMarker++) { + Roughness_Height[iMarker] = 0.0; + Kind_Wall[iMarker] = 1; + } + for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) { + Roughness_Height[nMarker_HeatFlux + iMarker] = 0.0; + Kind_Wall[nMarker_HeatFlux + iMarker] = 1; + } + } + else + SU2_MPI::Error("Mismatch in number of wall markers. Specify type of wall for all markers (even if smooth).", CURRENT_FUNCTION); + } else { + if (nRoughWall != nWall_Types) + SU2_MPI::Error("Mismatch in number of roughness height definition and wall type definition.", CURRENT_FUNCTION); + else { + su2double *temp_rough; + unsigned short *temp_kindrough; + temp_rough = new su2double [nWall_Types]; + temp_kindrough = new unsigned short [nWall_Types]; + for (unsigned short iMarker = 0; iMarker < nWall_Types; iMarker++) { + temp_rough[iMarker] = Roughness_Height[iMarker]; + temp_kindrough[iMarker] = Kind_Wall[iMarker]; + } + for (iMarker = 0; iMarker < nMarker_HeatFlux; iMarker++) { + Roughness_Height[iMarker] = temp_rough[iMarker]; + Kind_Wall[iMarker] = temp_kindrough[iMarker]; + } + for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) { + Roughness_Height[nMarker_HeatFlux + iMarker] = temp_rough[nMarker_HeatFlux + iMarker]; + Kind_Wall[nMarker_HeatFlux + iMarker] = temp_kindrough[nMarker_HeatFlux + iMarker]; + } + } + } + + } + /*--- Wall roughness handling ends here. ---*/ + /*--- Handle default options for topology optimization ---*/ if (topology_optimization && top_optim_nKernel==0) { @@ -7432,6 +7487,7 @@ CConfig::~CConfig(void) { if (Load_Sine_Amplitude != NULL) delete[] Load_Sine_Amplitude; if (Load_Sine_Frequency != NULL) delete[] Load_Sine_Frequency; if (FlowLoad_Value != NULL) delete[] FlowLoad_Value; + if (Roughness_Height != NULL) delete[] Roughness_Height; /*--- related to periodic boundary conditions ---*/ @@ -8600,6 +8656,49 @@ su2double CConfig::GetWall_HeatFlux(string val_marker) { return Heat_Flux[iMarker_HeatFlux]; } +unsigned short CConfig::GetKindWall(string val_marker) { + unsigned short iMarker; + short flag = -1; + for (iMarker = 0; iMarker < nMarker_HeatFlux; iMarker++) + if (Marker_HeatFlux[iMarker] == val_marker) { + flag = iMarker; + break; + } + if (flag == -1) + for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) + if (Marker_Isothermal[iMarker] == val_marker) { + flag = iMarker; + break; + } + if (flag == -1) return 1; + else return Kind_Wall[flag]; +} + +void CConfig::SetWall_RoughnessHeight(string val_marker, su2double val_roughness) { + //Roughness_Height[val_marker] = val_roughness; +} + +su2double CConfig::GetWall_RoughnessHeight(string val_marker) { + unsigned short iMarker = 0; + short flag = -1; + + if (nMarker_HeatFlux > 0 || nMarker_Isothermal > 0) { + for (iMarker = 0; iMarker < nMarker_HeatFlux; iMarker++) + if (Marker_HeatFlux[iMarker] == val_marker) { + flag = iMarker; + break; + } + if (flag == -1) { + for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) + if (Marker_Isothermal[iMarker] == val_marker) { + flag = iMarker; + break; + } + } + } + return Roughness_Height[flag]; +} + unsigned short CConfig::GetWallFunction_Treatment(string val_marker) { unsigned short WallFunction = NO_WALL_FUNCTION; diff --git a/Common/src/geometry/CPhysicalGeometry.cpp b/Common/src/geometry/CPhysicalGeometry.cpp index 93f63d4c121c..2267ec6dde91 100644 --- a/Common/src/geometry/CPhysicalGeometry.cpp +++ b/Common/src/geometry/CPhysicalGeometry.cpp @@ -5186,11 +5186,20 @@ void CPhysicalGeometry::ComputeWall_Distance(CConfig *config) { unsigned short markerID; unsigned long elemID; int rankID; - su2double dist; + su2double dist, roughness; + string Marker_Tag; WallADT.DetermineNearestElement(node[iPoint]->GetCoord(), dist, markerID, elemID, rankID); + node[iPoint]->SetWall_Distance(dist); + + /*--- Use the markerID to find the corresponding wall roughness height. ---*/ + + Marker_Tag = config->GetMarker_All_TagBound(markerID); + roughness = config->GetWall_RoughnessHeight(Marker_Tag); + + node[iPoint]->SetRoughnessHeight(roughness); } } diff --git a/SU2_CFD/include/numerics_structure.hpp b/SU2_CFD/include/numerics_structure.hpp index 68e61cf9b009..46d469e9613c 100644 --- a/SU2_CFD/include/numerics_structure.hpp +++ b/SU2_CFD/include/numerics_structure.hpp @@ -240,6 +240,8 @@ class CNumerics { su2double StrainMag_i, StrainMag_j; /*!< \brief Strain rate magnitude. */ su2double Dissipation_i, Dissipation_j; /*!< \brief Dissipation. */ su2double Dissipation_ij; + su2double roughness_i, /*!< \brief Roughness of the wall nearest to point i. */ + roughness_j; /*!< \brief Roughness of the wall nearest to point j. */ su2double *l, *m; @@ -560,6 +562,13 @@ class CNumerics { */ void SetDistance(su2double val_dist_i, su2double val_dist_j); + /*! + * \brief Set the value of the roughness from the nearest wall. + * \param[in] val_dist_i - Value of of the roughness of the nearest wall from point i + * \param[in] val_dist_j - Value of of the roughness of the nearest wall from point j + */ + void SetRoughness(su2double val_roughness_i, su2double val_roughness_j); + /*! * \brief Set coordinates of the points. * \param[in] val_coord_i - Coordinates of the point i. @@ -4147,6 +4156,7 @@ class CSourcePieceWise_TurbSA : public CNumerics { su2double sigma; su2double cb2; su2double cw1; + su2double cr1; unsigned short iDim; su2double nu, Ji, fv1, fv2, ft2, Omega, S, Shat, inv_Shat, dist_i_2, Ji_2, Ji_3, inv_k2_d2; su2double r, g, g_6, glim, fw; @@ -4156,6 +4166,7 @@ class CSourcePieceWise_TurbSA : public CNumerics { bool incompressible; bool rotating_frame; bool transition; + bool roughwall; su2double gamma_BC; su2double intermittency; su2double Production, Destruction, CrossProduction; diff --git a/SU2_CFD/include/numerics_structure.inl b/SU2_CFD/include/numerics_structure.inl index 676a63759cbd..c6af3552eeca 100644 --- a/SU2_CFD/include/numerics_structure.inl +++ b/SU2_CFD/include/numerics_structure.inl @@ -223,6 +223,11 @@ inline void CNumerics::SetDistance(su2double val_dist_i, su2double val_dist_j) { dist_j = val_dist_j; } +inline void CNumerics::SetRoughness(su2double val_roughness_i, su2double val_roughness_j) { + roughness_i = val_roughness_i; + roughness_j = val_roughness_j; +} + inline void CNumerics::SetAdjointVar(su2double *val_psi_i, su2double *val_psi_j) { Psi_i = val_psi_i; Psi_j = val_psi_j; diff --git a/SU2_CFD/src/numerics_direct_turbulent.cpp b/SU2_CFD/src/numerics_direct_turbulent.cpp index 853b5291c83a..36121f5a872a 100644 --- a/SU2_CFD/src/numerics_direct_turbulent.cpp +++ b/SU2_CFD/src/numerics_direct_turbulent.cpp @@ -315,6 +315,7 @@ CSourcePieceWise_TurbSA::CSourcePieceWise_TurbSA(unsigned short val_nDim, unsign cb2 = 0.622; cb2_sigma = cb2/sigma; cw1 = cb1/k2+(1.0+cb2)/sigma; + cr1 = 0.5; } @@ -333,6 +334,9 @@ void CSourcePieceWise_TurbSA::ComputeResidual(su2double *val_residual, su2double // BC Transition Model variables su2double vmag, rey, re_theta, re_theta_t, re_v; su2double tu , nu_cr, nu_t, nu_BC, chi_1, chi_2, term1, term2, term_exponential; + + // Set the boolean here depending on whether the point is closest to a rough wall or not. + roughwall = (roughness_i > 0.0); if (incompressible) { Density_i = V_i[nDim+2]; @@ -373,13 +377,28 @@ void CSourcePieceWise_TurbSA::ComputeResidual(su2double *val_residual, su2double /*--- Production term ---*/ + dist_i_2 = dist_i*dist_i; - nu = Laminar_Viscosity_i/Density_i; - Ji = TurbVar_i[0]/nu; + nu = Laminar_Viscosity_i/Density_i; + + /*--- Old values without roughness ---*/ + //Ji = TurbVar_i[0]/nu; + //Ji_2 = Ji*Ji; + //Ji_3 = Ji_2*Ji; + //fv1 = Ji_3/(Ji_3+cv1_3); + //fv2 = 1.0 - Ji/(1.0+Ji*fv1); + + /*--- Modified values for roughness ---*/ + /*--- Ref: Aupoix, B. and Spalart, P. R., "Extensions of the Spalart-Allmaras Turbulence Model to Account for Wall Roughness," + * International Journal of Heat and Fluid Flow, Vol. 24, 2003, pp. 454-462. ---*/ + /* --- See https://turbmodels.larc.nasa.gov/spalart.html#sarough for detailed explanation. ---*/ + + Ji = TurbVar_i[0]/nu + cr1*(roughness_i/dist_i); //roughness_i = 0 for smooth walls and Ji remains the same, changes only if roughness is specified. Ji_2 = Ji*Ji; Ji_3 = Ji_2*Ji; fv1 = Ji_3/(Ji_3+cv1_3); - fv2 = 1.0 - Ji/(1.0+Ji*fv1); + fv2 = 1.0 - TurbVar_i[0]/(nu+TurbVar_i[0]*fv1); //Using a modified relation so as to not change the Shat that depends on fv2. + ft2 = ct3*exp(-ct4*Ji_2); S = Omega; inv_k2_d2 = 1.0/(k2*dist_i_2); diff --git a/SU2_CFD/src/solver_direct_turbulent.cpp b/SU2_CFD/src/solver_direct_turbulent.cpp index 92e797b764a8..8a3a0a12c055 100644 --- a/SU2_CFD/src/solver_direct_turbulent.cpp +++ b/SU2_CFD/src/solver_direct_turbulent.cpp @@ -230,6 +230,10 @@ void CTurbSolver::Viscous_Residual(CGeometry *geometry, CSolver **solver_contain if ((config->GetKind_Turb_Model() == SST) || (config->GetKind_Turb_Model() == SST_SUST)) numerics->SetF1blending(nodes->GetF1blending(iPoint), nodes->GetF1blending(jPoint)); + /*--- Roughness heights. ---*/ + if (config->GetKind_Turb_Model() == SA) + numerics->SetRoughness(geometry->node[iPoint]->GetRoughnessHeight(),geometry->node[jPoint]->GetRoughnessHeight()); + /*--- Compute residual, and Jacobians ---*/ numerics->ComputeResidual(Residual, Jacobian_i, Jacobian_j, config); @@ -1279,6 +1283,7 @@ void CTurbSASolver::Source_Residual(CGeometry *geometry, CSolver **solver_contai bool harmonic_balance = (config->GetTime_Marching() == HARMONIC_BALANCE); bool transition = (config->GetKind_Trans_Model() == LM); bool transition_BC = (config->GetKind_Trans_Model() == BC); + su2double modifiedWallDistance = 0.0; for (iPoint = 0; iPoint < nPointDomain; iPoint++) { @@ -1314,10 +1319,23 @@ void CTurbSASolver::Source_Residual(CGeometry *geometry, CSolver **solver_contai /*--- Get Hybrid RANS/LES Type and set the appropriate wall distance ---*/ if (config->GetKind_HybridRANSLES() == NO_HYBRIDRANSLES) { + + /*--- For the SA model, wall roughness is accounted by modifying the computed wall distance + * d_new = d + 0.03 k_s + * where k_s is the equivalent sand grain roughness height that is specified in cfg file. + * For smooth walls, wall roughness is zero and computed wall distance remains the same. */ + + modifiedWallDistance = geometry->node[iPoint]->GetWall_Distance(); + + modifiedWallDistance = modifiedWallDistance + 0.03*geometry->node[iPoint]->GetRoughnessHeight(); /*--- Set distance to the surface ---*/ - numerics->SetDistance(geometry->node[iPoint]->GetWall_Distance(), 0.0); + numerics->SetDistance(modifiedWallDistance, 0.0); + + /*--- Set the roughness of the closest wall. ---*/ + + numerics->SetRoughness(geometry->node[iPoint]->GetRoughnessHeight(), 0.0 ); } else { @@ -1381,7 +1399,16 @@ void CTurbSASolver::Source_Template(CGeometry *geometry, CSolver **solver_contai void CTurbSASolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_container, CNumerics *conv_numerics, CNumerics *visc_numerics, CConfig *config, unsigned short val_marker) { unsigned long iPoint, iVertex; - unsigned short iVar; + unsigned short iVar, iDim; + bool rough_wall = false; + su2double RoughWallBC, Roughness_Height; + su2double *Res_Wall = new su2double [nVar]; + su2double *Normal, Area; + + string Marker_Tag = config->GetMarker_All_TagBound(val_marker); + + Roughness_Height = config->GetWall_RoughnessHeight(Marker_Tag); + if (Roughness_Height > 0.0 ) rough_wall = true; /*--- The dirichlet condition is used only without wall function, otherwise the convergence is compromised as we are providing nu tilde values for the @@ -1395,8 +1422,9 @@ void CTurbSASolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_conta /*--- Check if the node belongs to the domain (i.e, not a halo node) ---*/ if (geometry->node[iPoint]->GetDomain()) { - - /*--- Get the velocity vector ---*/ + + if (!rough_wall) { + /*--- Get the solution vector ---*/ for (iVar = 0; iVar < nVar; iVar++) Solution[iVar] = 0.0; @@ -1407,6 +1435,28 @@ void CTurbSASolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_conta /*--- Includes 1 in the diagonal ---*/ Jacobian.DeleteValsRowi(iPoint); + } + else { + /*--- Compute dual-grid area and boundary normal ---*/ + + Normal = geometry->vertex[val_marker][iVertex]->GetNormal(); + + Area = 0.0; + for (iDim = 0; iDim < nDim; iDim++) + Area += Normal[iDim]*Normal[iDim]; + Area = sqrt (Area); + + /*--- For rough walls, the boundary condition is given by + * (\frac{\partial \nu}{\partial n})_wall = \frac{\nu}{0.03*k_s} + * where \nu is the solution variable, $n$ is the wall normal direction + * and k_s is the equivalent sand grain roughness specified. ---*/ + RoughWallBC = nodes->GetSolution(iPoint,0)/(0.03*Roughness_Height); + Res_Wall[0] = RoughWallBC*Area; + LinSysRes.SubtractBlock(iPoint, Res_Wall); + + Jacobian_i[0][0] = Area/(0.03*Roughness_Height); + Jacobian.SubtractBlock(iPoint,iPoint,Jacobian_i); + } } } } @@ -3671,31 +3721,77 @@ void CTurbSSTSolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_cont unsigned long iPoint, jPoint, iVertex, total_index; unsigned short iDim, iVar; su2double distance, density = 0.0, laminar_viscosity = 0.0, beta_1; + bool rough_wall = false; + su2double RoughWallBC, Roughness_Height, S_R, FrictionVel, kPlus, WallShearStress; + string Marker_Tag = config->GetMarker_All_TagBound(val_marker); + + + if (config->GetKindWall(Marker_Tag) == ROUGH ) rough_wall = true; for (iVertex = 0; iVertex < geometry->nVertex[val_marker]; iVertex++) { iPoint = geometry->vertex[val_marker][iVertex]->GetNode(); /*--- Check if the node belongs to the domain (i.e, not a halo node) ---*/ if (geometry->node[iPoint]->GetDomain()) { - - /*--- distance to closest neighbor ---*/ - jPoint = geometry->vertex[val_marker][iVertex]->GetNormal_Neighbor(); - distance = 0.0; - for (iDim = 0; iDim < nDim; iDim++) { - distance += (geometry->node[iPoint]->GetCoord(iDim) - geometry->node[jPoint]->GetCoord(iDim))* - (geometry->node[iPoint]->GetCoord(iDim) - geometry->node[jPoint]->GetCoord(iDim)); - } - distance = sqrt(distance); - - /*--- Set wall values ---*/ - - density = solver_container[FLOW_SOL]->GetNodes()->GetDensity(jPoint); - laminar_viscosity = solver_container[FLOW_SOL]->GetNodes()->GetLaminarViscosity(jPoint); - - beta_1 = constants[4]; - + + /*--- Set wall values ---*/ + + density = solver_container[FLOW_SOL]->GetNodes()->GetDensity(iPoint); + laminar_viscosity = solver_container[FLOW_SOL]->GetNodes()->GetLaminarViscosity(iPoint); + + if (rough_wall) { + + WallShearStress = solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 0)*solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 0); + WallShearStress += solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 1)*solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 1); + if (nDim == 3) WallShearStress += solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 2)*solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 2); + + WallShearStress = sqrt(WallShearStress); + /*--- Compute non-dimensional velocity ---*/ + FrictionVel = sqrt(fabs(WallShearStress)/density); + + /*--- Compute roughness in wall units. ---*/ + Roughness_Height = config->GetWall_RoughnessHeight(Marker_Tag); + kPlus = FrictionVel*Roughness_Height*density/laminar_viscosity; + + /*--- Reference 1 from Sandia and Aupoix ---*/ + if (kPlus <= 25) + S_R = (50/(kPlus+EPS))*(50/(kPlus+EPS)); + else + S_R = 100/(kPlus+EPS); + + /*--- Reference 2 from D.C. Wilcox Turbulence Modeling for CFD ---*/ + /*if (kPlus <= 5) + S_R = (200/(kPlus+EPS))*(200/(kPlus+EPS)); + else + S_R = 100/(kPlus+EPS) + ((200/(kPlus+EPS))*(200/(kPlus+EPS)) - 100/(kPlus+EPS))*exp(5-kPlus);*/ + + + /*--- Modify the omega to account for a rough wall. ---*/ + Solution[1] = FrictionVel*FrictionVel*S_R/(laminar_viscosity/density); + + + + + } else { + /*--- distance to closest neighbor ---*/ + jPoint = geometry->vertex[val_marker][iVertex]->GetNormal_Neighbor(); + distance = 0.0; + for (iDim = 0; iDim < nDim; iDim++) { + distance += (geometry->node[iPoint]->GetCoord(iDim) - geometry->node[jPoint]->GetCoord(iDim))* + (geometry->node[iPoint]->GetCoord(iDim) - geometry->node[jPoint]->GetCoord(iDim)); + } + distance = sqrt(distance); + + /*--- Set wall values ---*/ + + beta_1 = constants[4]; + + + Solution[1] = 60.0*laminar_viscosity/(density*beta_1*distance*distance); + } + + /*--- Only the \omega solution is modified at the wall. ---*/ Solution[0] = 0.0; - Solution[1] = 60.0*laminar_viscosity/(density*beta_1*distance*distance); /*--- Set the solution values and zero the residual ---*/ nodes->SetSolution_Old(iPoint,Solution); @@ -3707,7 +3803,6 @@ void CTurbSSTSolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_cont total_index = iPoint*nVar+iVar; Jacobian.DeleteValsRowi(total_index); } - } } From 75e68fcdce91a405db0c0374cb8979e51bb2836f Mon Sep 17 00:00:00 2001 From: koodlyakshay Date: Fri, 31 Jan 2020 10:00:13 +0100 Subject: [PATCH 02/29] Small changes. --- SU2_CFD/src/numerics_direct_turbulent.cpp | 2 +- SU2_CFD/src/solvers/CTurbSSTSolver.cpp | 44 +++++++++++------------ 2 files changed, 21 insertions(+), 25 deletions(-) diff --git a/SU2_CFD/src/numerics_direct_turbulent.cpp b/SU2_CFD/src/numerics_direct_turbulent.cpp index 36121f5a872a..9f0a705e9e42 100644 --- a/SU2_CFD/src/numerics_direct_turbulent.cpp +++ b/SU2_CFD/src/numerics_direct_turbulent.cpp @@ -425,7 +425,7 @@ void CSourcePieceWise_TurbSA::ComputeResidual(su2double *val_residual, su2double re_theta_t = (803.73 * pow((tu + 0.6067),-1.027)); //MENTER correlation //re_theta_t = 163.0 + exp(6.91-tu); //ABU-GHANNAM & SHAW correlation - term1 = sqrt(max(re_theta-re_theta_t,0.)/(chi_1*re_theta_t)); + term1 = sqrt(max(re_theta-re_theta_t , 0.0)/(chi_1*re_theta_t)); term2 = sqrt(max(nu_BC-nu_cr,0.)/(nu_cr)); term_exponential = (term1 + term2); gamma_BC = 1.0 - exp(-term_exponential); diff --git a/SU2_CFD/src/solvers/CTurbSSTSolver.cpp b/SU2_CFD/src/solvers/CTurbSSTSolver.cpp index 87a4b3767d38..09158a004cd7 100644 --- a/SU2_CFD/src/solvers/CTurbSSTSolver.cpp +++ b/SU2_CFD/src/solvers/CTurbSSTSolver.cpp @@ -468,37 +468,33 @@ void CTurbSSTSolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_cont laminar_viscosity = solver_container[FLOW_SOL]->GetNodes()->GetLaminarViscosity(iPoint); if (rough_wall) { - - WallShearStress = solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 0)*solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 0); - WallShearStress += solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 1)*solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 1); - if (nDim == 3) WallShearStress += solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 2)*solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 2); - - WallShearStress = sqrt(WallShearStress); - /*--- Compute non-dimensional velocity ---*/ - FrictionVel = sqrt(fabs(WallShearStress)/density); - /*--- Compute roughness in wall units. ---*/ - Roughness_Height = config->GetWall_RoughnessHeight(Marker_Tag); - kPlus = FrictionVel*Roughness_Height*density/laminar_viscosity; + WallShearStress = solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 0)*solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 0); + WallShearStress += solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 1)*solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 1); + if (nDim == 3) WallShearStress += solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 2)*solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 2); - /*--- Reference 1 from Sandia and Aupoix ---*/ - if (kPlus <= 25) - S_R = (50/(kPlus+EPS))*(50/(kPlus+EPS)); - else - S_R = 100/(kPlus+EPS); + WallShearStress = sqrt(WallShearStress); + /*--- Compute non-dimensional velocity ---*/ + FrictionVel = sqrt(fabs(WallShearStress)/density); + + /*--- Compute roughness in wall units. ---*/ + Roughness_Height = config->GetWall_RoughnessHeight(Marker_Tag); + kPlus = FrictionVel*Roughness_Height*density/laminar_viscosity; - /*--- Reference 2 from D.C. Wilcox Turbulence Modeling for CFD ---*/ - /*if (kPlus <= 5) - S_R = (200/(kPlus+EPS))*(200/(kPlus+EPS)); + /*--- Reference 1 original Wilcox (1998) ---*/ + /*if (kPlus <= 25) + S_R = (50/(kPlus+EPS))*(50/(kPlus+EPS)); else - S_R = 100/(kPlus+EPS) + ((200/(kPlus+EPS))*(200/(kPlus+EPS)) - 100/(kPlus+EPS))*exp(5-kPlus);*/ + S_R = 100/(kPlus+EPS);*/ + /*--- Reference 2 from D.C. Wilcox Turbulence Modeling for CFD (2006) ---*/ + if (kPlus <= 5) + S_R = (200/(kPlus+EPS))*(200/(kPlus+EPS)); + else + S_R = 100/(kPlus+EPS) + ((200/(kPlus+EPS))*(200/(kPlus+EPS)) - 100/(kPlus+EPS))*exp(5-kPlus); /*--- Modify the omega to account for a rough wall. ---*/ - Solution[1] = FrictionVel*FrictionVel*S_R/(laminar_viscosity/density); - - - + Solution[1] = FrictionVel*FrictionVel*S_R/(laminar_viscosity/density); } else { /*--- distance to closest neighbor ---*/ From 015b46bd3845d03f0741b352cf9f8f9be2308ebb Mon Sep 17 00:00:00 2001 From: koodlyakshay Date: Sun, 2 Feb 2020 22:20:32 +0100 Subject: [PATCH 03/29] Fix boundary condition implementation. --- SU2_CFD/src/solvers/CTurbSASolver.cpp | 28 +++++++++++++++++++-------- 1 file changed, 20 insertions(+), 8 deletions(-) diff --git a/SU2_CFD/src/solvers/CTurbSASolver.cpp b/SU2_CFD/src/solvers/CTurbSASolver.cpp index 51fc62e8161b..cd549cf3193d 100644 --- a/SU2_CFD/src/solvers/CTurbSASolver.cpp +++ b/SU2_CFD/src/solvers/CTurbSASolver.cpp @@ -502,7 +502,7 @@ void CTurbSASolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_conta unsigned long iPoint, iVertex; unsigned short iVar, iDim; bool rough_wall = false; - su2double RoughWallBC, Roughness_Height; + su2double RoughWallBC, Roughness_Height, laminar_viscosity, density, sigma = 2.0/3.0, nu_total, coeff; su2double *Res_Wall = new su2double [nVar]; su2double *Normal, Area; @@ -538,24 +538,36 @@ void CTurbSASolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_conta Jacobian.DeleteValsRowi(iPoint); } else { + /*--- For rough walls, the boundary condition is given by + * (\frac{\partial \nu}{\partial n})_wall = \frac{\nu}{0.03*k_s} + * where \nu is the solution variable, $n$ is the wall normal direction + * and k_s is the equivalent sand grain roughness specified. ---*/ + /*--- Compute dual-grid area and boundary normal ---*/ Normal = geometry->vertex[val_marker][iVertex]->GetNormal(); - Area = 0.0; + Area = 0.0; for (iDim = 0; iDim < nDim; iDim++) Area += Normal[iDim]*Normal[iDim]; Area = sqrt (Area); + + /*--- Get laminar_viscosity and density ---*/ - /*--- For rough walls, the boundary condition is given by - * (\frac{\partial \nu}{\partial n})_wall = \frac{\nu}{0.03*k_s} - * where \nu is the solution variable, $n$ is the wall normal direction - * and k_s is the equivalent sand grain roughness specified. ---*/ + laminar_viscosity = solver_container[FLOW_SOL]->GetNodes()->GetLaminarViscosity(iPoint); + density = solver_container[FLOW_SOL]->GetNodes()->GetDensity(iPoint); + + nu_total = (laminar_viscosity/density + nodes->GetSolution(iPoint,0)); + + coeff = (nu_total/sigma); + RoughWallBC = nodes->GetSolution(iPoint,0)/(0.03*Roughness_Height); - Res_Wall[0] = RoughWallBC*Area; + + Res_Wall[0] = coeff*RoughWallBC*Area; LinSysRes.SubtractBlock(iPoint, Res_Wall); - Jacobian_i[0][0] = Area/(0.03*Roughness_Height); + Jacobian_i[0][0] = (laminar_viscosity*Area)/(0.03*Roughness_Height*sigma); + Jacobian_i[0][0] += 2.0*RoughWallBC*Area/sigma; Jacobian.SubtractBlock(iPoint,iPoint,Jacobian_i); } } From f1f35e259852b62a4f1cad0c27cf439ec04d89d5 Mon Sep 17 00:00:00 2001 From: koodlyakshay Date: Thu, 6 Feb 2020 21:50:21 +0100 Subject: [PATCH 04/29] Temporary fix for wall distance modification when run in parallel. --- Common/include/adt_structure.hpp | 5 +++- Common/src/adt_structure.cpp | 31 ++++++++++++++-------- Common/src/fem_geometry_structure.cpp | 5 +++- Common/src/fem_wall_distance.cpp | 17 ++++++++---- Common/src/geometry/CPhysicalGeometry.cpp | 22 ++++++++++----- Common/src/geometry_structure_fem_part.cpp | 5 +++- SU2_CFD/src/solvers/CTurbSASolver.cpp | 5 ++-- 7 files changed, 63 insertions(+), 27 deletions(-) diff --git a/Common/include/adt_structure.hpp b/Common/include/adt_structure.hpp index 9b89441436cb..49eac288fdf0 100644 --- a/Common/include/adt_structure.hpp +++ b/Common/include/adt_structure.hpp @@ -350,6 +350,7 @@ class CADTElemClass : public CADTBaseClass { vector BBoxTargets; /*!< \brief Vector, used to store possible bounding box candidates during the nearest element search. */ + vector markerRoughness; public: /*! * \brief Constructor of the class. @@ -370,6 +371,7 @@ class CADTElemClass : public CADTBaseClass { vector &val_VTKElem, vector &val_markerID, vector &val_elemID, + vector &val_roughheight, const bool globalTree); /*! @@ -410,7 +412,8 @@ class CADTElemClass : public CADTBaseClass { su2double &dist, unsigned short &markerID, unsigned long &elemID, - int &rankID); + int &rankID, + su2double &localRoughness); private: /*! diff --git a/Common/src/adt_structure.cpp b/Common/src/adt_structure.cpp index a75c4b235a5b..d1e5a7b91551 100644 --- a/Common/src/adt_structure.cpp +++ b/Common/src/adt_structure.cpp @@ -484,6 +484,7 @@ CADTElemClass::CADTElemClass(unsigned short val_nDim, vector &val_VTKElem, vector &val_markerID, vector &val_elemID, + vector &val_roughheight, const bool globalTree) { /* Copy the dimension of the problem into nDim. */ @@ -552,6 +553,7 @@ CADTElemClass::CADTElemClass(unsigned short val_nDim, elemVTK_Type.resize(sizeGlobal); localMarkers.resize(sizeGlobal); + markerRoughness.resize(sizeGlobal); localElemIDs.resize(sizeGlobal); SU2_MPI::Allgatherv(val_VTKElem.data(), sizeLocal, MPI_UNSIGNED_SHORT, elemVTK_Type.data(), @@ -559,6 +561,9 @@ CADTElemClass::CADTElemClass(unsigned short val_nDim, SU2_MPI::Allgatherv(val_markerID.data(), sizeLocal, MPI_UNSIGNED_SHORT, localMarkers.data(), recvCounts.data(), displs.data(), MPI_UNSIGNED_SHORT, MPI_COMM_WORLD); + + SU2_MPI::Allgatherv(val_roughheight.data(), sizeLocal, MPI_DOUBLE, markerRoughness.data(), + recvCounts.data(), displs.data(), MPI_DOUBLE, MPI_COMM_WORLD); SU2_MPI::Allgatherv(val_elemID.data(), sizeLocal, MPI_UNSIGNED_LONG, localElemIDs.data(), recvCounts.data(), displs.data(), MPI_UNSIGNED_LONG, MPI_COMM_WORLD); @@ -596,11 +601,12 @@ CADTElemClass::CADTElemClass(unsigned short val_nDim, int rank; SU2_MPI::Comm_rank(MPI_COMM_WORLD, &rank); - coorPoints = val_coor; - elemConns = val_connElem; - elemVTK_Type = val_VTKElem; - localMarkers = val_markerID; - localElemIDs = val_elemID; + coorPoints = val_coor; + elemConns = val_connElem; + elemVTK_Type = val_VTKElem; + localMarkers = val_markerID; + localElemIDs = val_elemID; + markerRoughness = val_roughheight; ranksOfElems.assign(elemVTK_Type.size(), rank); } @@ -608,11 +614,12 @@ CADTElemClass::CADTElemClass(unsigned short val_nDim, /*--- Sequential mode. Copy the data from the arguments into the member variables and set the ranks to MASTER_NODE. ---*/ - coorPoints = val_coor; - elemConns = val_connElem; - elemVTK_Type = val_VTKElem; - localMarkers = val_markerID; - localElemIDs = val_elemID; + coorPoints = val_coor; + elemConns = val_connElem; + elemVTK_Type = val_VTKElem; + localMarkers = val_markerID; + localElemIDs = val_elemID; + markerRoughness = val_roughheight; ranksOfElems.assign(elemVTK_Type.size(), MASTER_NODE); @@ -789,7 +796,8 @@ void CADTElemClass::DetermineNearestElement(const su2double *coor, su2double &dist, unsigned short &markerID, unsigned long &elemID, - int &rankID) { + int &rankID, + su2double &localRoughness) { AD_BEGIN_PASSIVE @@ -965,6 +973,7 @@ void CADTElemClass::DetermineNearestElement(const su2double *coor, markerID = localMarkers[ii]; elemID = localElemIDs[ii]; rankID = ranksOfElems[ii]; + localRoughness = markerRoughness[ii]; } } diff --git a/Common/src/fem_geometry_structure.cpp b/Common/src/fem_geometry_structure.cpp index fa33a7d08da7..c3bcd0180e68 100644 --- a/Common/src/fem_geometry_structure.cpp +++ b/Common/src/fem_geometry_structure.cpp @@ -6220,6 +6220,7 @@ void CMeshFEM_DG::WallFunctionPreprocessing(CConfig *config) { vector subElementIDInParent; vector VTK_TypeElem; vector elemConn; + vector dummyRough; /* Loop over the locally stored volume elements (including halo elements) to create the connectivity of the subelements. */ @@ -6256,6 +6257,7 @@ void CMeshFEM_DG::WallFunctionPreprocessing(CConfig *config) { for(unsigned short j=0; j().swap(parentElement); vector().swap(elemConn); vector().swap(volCoor); + vector().swap(dummyRough); /*--------------------------------------------------------------------------*/ /*--- Step 3. Search for donor elements at the exchange locations in ---*/ diff --git a/Common/src/fem_wall_distance.cpp b/Common/src/fem_wall_distance.cpp index 88b9b8f2d5de..a630fa791c42 100644 --- a/Common/src/fem_wall_distance.cpp +++ b/Common/src/fem_wall_distance.cpp @@ -48,6 +48,7 @@ void CMeshFEM_DG::ComputeWall_Distance(CConfig *config) { vector elemIDs; vector VTK_TypeElem; vector markerIDs; + vector dummyRoughs; /* Loop over the boundary markers. */ for(unsigned short iMarker=0; iMarker().swap(elemIDs); vector().swap(surfaceConn); vector().swap(surfaceCoor); + vector().swap(dummyRoughs); /*--------------------------------------------------------------------------*/ /*--- Step 3: Determine the wall distance of the integration points of ---*/ @@ -162,7 +165,8 @@ void CMeshFEM_DG::ComputeWall_Distance(CConfig *config) { unsigned long elemID; int rankID; su2double dist; - WallADT.DetermineNearestElement(coor, dist, markerID, elemID, rankID); + su2double dumRough; + WallADT.DetermineNearestElement(coor, dist, markerID, elemID, rankID, dumRough); volElem[l].wallDistance[i] = dist; } @@ -202,7 +206,8 @@ void CMeshFEM_DG::ComputeWall_Distance(CConfig *config) { unsigned long elemID; int rankID; su2double dist; - WallADT.DetermineNearestElement(coor, dist, markerID, elemID, rankID); + su2double dumRough; + WallADT.DetermineNearestElement(coor, dist, markerID, elemID, rankID, dumRough); volElem[l].wallDistanceSolDOFs[i] = dist; } @@ -242,7 +247,8 @@ void CMeshFEM_DG::ComputeWall_Distance(CConfig *config) { unsigned long elemID; int rankID; su2double dist; - WallADT.DetermineNearestElement(coor, dist, markerID, elemID, rankID); + su2double dumRough; + WallADT.DetermineNearestElement(coor, dist, markerID, elemID, rankID, dumRough); matchingFaces[l].wallDistance[i] = dist; } @@ -296,7 +302,8 @@ void CMeshFEM_DG::ComputeWall_Distance(CConfig *config) { unsigned long elemID; int rankID; su2double dist; - WallADT.DetermineNearestElement(coor, dist, markerID, elemID, rankID); + su2double dumRough; + WallADT.DetermineNearestElement(coor, dist, markerID, elemID, rankID, dumRough); surfElem[l].wallDistance[i] = dist; } diff --git a/Common/src/geometry/CPhysicalGeometry.cpp b/Common/src/geometry/CPhysicalGeometry.cpp index e595edac45d9..6391f8151586 100644 --- a/Common/src/geometry/CPhysicalGeometry.cpp +++ b/Common/src/geometry/CPhysicalGeometry.cpp @@ -5098,6 +5098,9 @@ void CPhysicalGeometry::ComputeWall_Distance(CConfig *config) { vector elemIDs; vector VTK_TypeElem; vector markerIDs; + vector roughheights; + string Marker_Tag; + su2double roughness; /* Loop over the boundary markers. */ @@ -5106,6 +5109,9 @@ void CPhysicalGeometry::ComputeWall_Distance(CConfig *config) { /* Check for a viscous wall. */ if( config->GetViscous_Wall(iMarker)) { + + Marker_Tag = config->GetMarker_All_TagBound(iMarker); + roughness = config->GetWall_RoughnessHeight(Marker_Tag); /* Loop over the surface elements of this marker. */ for(unsigned long iElem=0; iElem < nElem_Bound[iMarker]; iElem++) { @@ -5126,7 +5132,8 @@ void CPhysicalGeometry::ComputeWall_Distance(CConfig *config) { markerIDs.push_back(iMarker); VTK_TypeElem.push_back(VTK_Type); elemIDs.push_back(iElem); - + roughheights.push_back(roughness); + for (unsigned short iNode = 0; iNode < nDOFsPerElem; iNode++) surfaceConn.push_back(bound[iMarker][iElem]->GetNode(iNode)); } @@ -5163,7 +5170,7 @@ void CPhysicalGeometry::ComputeWall_Distance(CConfig *config) { /* Build the ADT. */ CADTElemClass WallADT(nDim, surfaceCoor, surfaceConn, VTK_TypeElem, - markerIDs, elemIDs, true); + markerIDs, elemIDs, roughheights, true); /* Release the memory of the vectors used to build the ADT. To make sure that all the memory is deleted, the swap function is used. */ @@ -5172,6 +5179,7 @@ void CPhysicalGeometry::ComputeWall_Distance(CConfig *config) { vector().swap(elemIDs); vector().swap(surfaceConn); vector().swap(surfaceCoor); + vector().swap(roughheights); /*--------------------------------------------------------------------------*/ /*--- Step 3: Loop over all interior mesh nodes and compute minimum ---*/ @@ -5196,11 +5204,11 @@ void CPhysicalGeometry::ComputeWall_Distance(CConfig *config) { unsigned short markerID; unsigned long elemID; int rankID; - su2double dist, roughness; - string Marker_Tag; + su2double dist; + su2double localRoughness; WallADT.DetermineNearestElement(node[iPoint]->GetCoord(), dist, markerID, - elemID, rankID); + elemID, rankID, localRoughness); node[iPoint]->SetWall_Distance(dist); @@ -5209,7 +5217,9 @@ void CPhysicalGeometry::ComputeWall_Distance(CConfig *config) { Marker_Tag = config->GetMarker_All_TagBound(markerID); roughness = config->GetWall_RoughnessHeight(Marker_Tag); - node[iPoint]->SetRoughnessHeight(roughness); + node[iPoint]->SetRoughnessHeight(localRoughness); + + //if (rank ==2) cout<GetCoord(0)<<"\t"<GetCoord(1)<<"\t"< subElementIDInParent; vector VTK_TypeElem; vector elemConn; + vector dummyRough; /* Loop over the local volume elements to create the connectivity of the linear sub-elements. */ @@ -3478,6 +3479,7 @@ void CPhysicalGeometry::DetermineDonorElementsWallFunctions(CConfig *config) { parentElement.push_back(elem[l]->GetGlobalElemID()); subElementIDInParent.push_back(jj); VTK_TypeElem.push_back(VTK_Type[i]); + dummyRough.push_back(0.0); for(unsigned short k=0; kGetNode(connSubElems[i][kk]); @@ -3502,7 +3504,7 @@ void CPhysicalGeometry::DetermineDonorElementsWallFunctions(CConfig *config) { /* Build the local ADT. */ CADTElemClass localVolumeADT(nDim, volCoor, elemConn, VTK_TypeElem, - subElementIDInParent, parentElement, false); + subElementIDInParent, parentElement, dummyRough, false); /* Release the memory of the vectors used to build the ADT. To make sure that all the memory is deleted, the swap function is used. */ @@ -3511,6 +3513,7 @@ void CPhysicalGeometry::DetermineDonorElementsWallFunctions(CConfig *config) { vector().swap(parentElement); vector().swap(elemConn); vector().swap(volCoor); + vector().swap(dummyRough); /*--------------------------------------------------------------------------*/ /*--- Step 3. Search for donor elements at the exchange locations in ---*/ diff --git a/SU2_CFD/src/solvers/CTurbSASolver.cpp b/SU2_CFD/src/solvers/CTurbSASolver.cpp index cd549cf3193d..468c092ffd1b 100644 --- a/SU2_CFD/src/solvers/CTurbSASolver.cpp +++ b/SU2_CFD/src/solvers/CTurbSASolver.cpp @@ -510,6 +510,7 @@ void CTurbSASolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_conta Roughness_Height = config->GetWall_RoughnessHeight(Marker_Tag); if (Roughness_Height > 0.0 ) rough_wall = true; + unsigned short rank = SU2_MPI::GetRank(); /*--- The dirichlet condition is used only without wall function, otherwise the convergence is compromised as we are providing nu tilde values for the @@ -532,7 +533,7 @@ void CTurbSASolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_conta nodes->SetSolution_Old(iPoint,Solution); LinSysRes.SetBlock_Zero(iPoint); - + /*--- Includes 1 in the diagonal ---*/ Jacobian.DeleteValsRowi(iPoint); @@ -550,7 +551,7 @@ void CTurbSASolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_conta Area = 0.0; for (iDim = 0; iDim < nDim; iDim++) Area += Normal[iDim]*Normal[iDim]; - Area = sqrt (Area); + Area = sqrt(Area); /*--- Get laminar_viscosity and density ---*/ From a211446f563b8d2979e9d0bc03f1ced6a29a4a42 Mon Sep 17 00:00:00 2001 From: koodlyakshay Date: Sat, 8 Feb 2020 23:45:46 +0100 Subject: [PATCH 05/29] Fix postprocessing routines. --- Common/src/CConfig.cpp | 3 +-- Common/src/geometry/CPhysicalGeometry.cpp | 3 +-- SU2_CFD/src/numerics_direct_turbulent.cpp | 4 +++- SU2_CFD/src/solvers/CIncNSSolver.cpp | 5 ++++- SU2_CFD/src/solvers/CTurbSASolver.cpp | 8 ++++++-- 5 files changed, 15 insertions(+), 8 deletions(-) diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index b220dc86571b..af0ed50e0649 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -4642,9 +4642,8 @@ void CConfig::SetPostprocessing(unsigned short val_software, unsigned short val_ /*--- Check that if the wall roughness array are compatible and set deafult values if needed. ---*/ if ((nMarker_HeatFlux > 0) || (nMarker_Isothermal > 0)) { - /*--- If the config option is not declared, then assume smooth walls and set roughness to zero. ---*/ if (nWall_Types != nMarker_HeatFlux + nMarker_Isothermal) { - /*--- If this option is not used at all, assume all walls are smooth and set roughness height to zero. ---*/ + /*--- If the config option is not declared, then assume smooth walls and set roughness to zero. ---*/ if (nWall_Types == 0) { if (Roughness_Height != NULL) Roughness_Height = NULL; Roughness_Height = new su2double [nWall_Types]; diff --git a/Common/src/geometry/CPhysicalGeometry.cpp b/Common/src/geometry/CPhysicalGeometry.cpp index 6391f8151586..5a76c2e82bd9 100644 --- a/Common/src/geometry/CPhysicalGeometry.cpp +++ b/Common/src/geometry/CPhysicalGeometry.cpp @@ -5218,8 +5218,7 @@ void CPhysicalGeometry::ComputeWall_Distance(CConfig *config) { roughness = config->GetWall_RoughnessHeight(Marker_Tag); node[iPoint]->SetRoughnessHeight(localRoughness); - - //if (rank ==2) cout<GetCoord(0)<<"\t"<GetCoord(1)<<"\t"<GetGradient_Primitive(iPoint,nDim+1, iDim); } + TurbViscosity = nodes->GetEddyViscosity(iPoint); Viscosity = nodes->GetLaminarViscosity(iPoint); + Viscosity += TurbViscosity; + Density = nodes->GetDensity(iPoint); Area = 0.0; for (iDim = 0; iDim < nDim; iDim++) Area += Normal[iDim]*Normal[iDim]; Area = sqrt(Area); diff --git a/SU2_CFD/src/solvers/CTurbSASolver.cpp b/SU2_CFD/src/solvers/CTurbSASolver.cpp index 468c092ffd1b..8cb372f34fbf 100644 --- a/SU2_CFD/src/solvers/CTurbSASolver.cpp +++ b/SU2_CFD/src/solvers/CTurbSASolver.cpp @@ -347,7 +347,7 @@ void CTurbSASolver::Preprocessing(CGeometry *geometry, CSolver **solver_containe void CTurbSASolver::Postprocessing(CGeometry *geometry, CSolver **solver_container, CConfig *config, unsigned short iMesh) { su2double rho = 0.0, mu = 0.0, nu, *nu_hat, muT, Ji, Ji_3, fv1; - su2double cv1_3 = 7.1*7.1*7.1; + su2double cv1_3 = 7.1*7.1*7.1, cR1 = 0.5, roughness, dist; unsigned long iPoint; bool neg_spalart_allmaras = (config->GetKind_Turb_Model() == SA_NEG); @@ -361,8 +361,12 @@ void CTurbSASolver::Postprocessing(CGeometry *geometry, CSolver **solver_contain nu = mu/rho; nu_hat = nodes->GetSolution(iPoint); + + roughness = geometry->node[iPoint]->GetRoughnessHeight(); + dist = geometry->node[iPoint]->GetWall_Distance(); + dist += 0.03*roughness; - Ji = nu_hat[0]/nu; + Ji = nu_hat[0]/nu + cR1*roughness/dist; Ji_3 = Ji*Ji*Ji; fv1 = Ji_3/(Ji_3+cv1_3); From 305667bef75b09779db69195e48f4471d0a3159a Mon Sep 17 00:00:00 2001 From: koodlyakshay Date: Sun, 9 Feb 2020 00:30:59 +0100 Subject: [PATCH 06/29] Fix for zero wall distance and zero roughness by adding EPS. --- Common/src/CConfig.cpp | 3 ++- SU2_CFD/src/numerics_direct_turbulent.cpp | 2 +- SU2_CFD/src/solvers/CTurbSASolver.cpp | 2 +- 3 files changed, 4 insertions(+), 3 deletions(-) diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index af0ed50e0649..b220dc86571b 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -4642,8 +4642,9 @@ void CConfig::SetPostprocessing(unsigned short val_software, unsigned short val_ /*--- Check that if the wall roughness array are compatible and set deafult values if needed. ---*/ if ((nMarker_HeatFlux > 0) || (nMarker_Isothermal > 0)) { + /*--- If the config option is not declared, then assume smooth walls and set roughness to zero. ---*/ if (nWall_Types != nMarker_HeatFlux + nMarker_Isothermal) { - /*--- If the config option is not declared, then assume smooth walls and set roughness to zero. ---*/ + /*--- If this option is not used at all, assume all walls are smooth and set roughness height to zero. ---*/ if (nWall_Types == 0) { if (Roughness_Height != NULL) Roughness_Height = NULL; Roughness_Height = new su2double [nWall_Types]; diff --git a/SU2_CFD/src/numerics_direct_turbulent.cpp b/SU2_CFD/src/numerics_direct_turbulent.cpp index 608b594a4d20..25b8766ec2db 100644 --- a/SU2_CFD/src/numerics_direct_turbulent.cpp +++ b/SU2_CFD/src/numerics_direct_turbulent.cpp @@ -393,7 +393,7 @@ void CSourcePieceWise_TurbSA::ComputeResidual(su2double *val_residual, su2double * International Journal of Heat and Fluid Flow, Vol. 24, 2003, pp. 454-462. ---*/ /* --- See https://turbmodels.larc.nasa.gov/spalart.html#sarough for detailed explanation. ---*/ - Ji = TurbVar_i[0]/nu + cr1*(roughness_i/dist_i); //roughness_i = 0 for smooth walls and Ji remains the same, changes only if roughness is specified. + Ji = TurbVar_i[0]/nu + cr1*(roughness_i/(dist_i+EPS)); //roughness_i = 0 for smooth walls and Ji remains the same, changes only if roughness is specified. Ji_2 = Ji*Ji; Ji_3 = Ji_2*Ji; fv1 = Ji_3/(Ji_3+cv1_3); diff --git a/SU2_CFD/src/solvers/CTurbSASolver.cpp b/SU2_CFD/src/solvers/CTurbSASolver.cpp index 8cb372f34fbf..424c7888b859 100644 --- a/SU2_CFD/src/solvers/CTurbSASolver.cpp +++ b/SU2_CFD/src/solvers/CTurbSASolver.cpp @@ -366,7 +366,7 @@ void CTurbSASolver::Postprocessing(CGeometry *geometry, CSolver **solver_contain dist = geometry->node[iPoint]->GetWall_Distance(); dist += 0.03*roughness; - Ji = nu_hat[0]/nu + cR1*roughness/dist; + Ji = nu_hat[0]/nu + cR1*roughness/(dist+EPS); Ji_3 = Ji*Ji*Ji; fv1 = Ji_3/(Ji_3+cv1_3); From 2fb7cc1954f3bc59f95bd61818fd57cc452fa759 Mon Sep 17 00:00:00 2001 From: koodlyakshay Date: Tue, 11 Feb 2020 13:02:21 +0100 Subject: [PATCH 07/29] Add some comments. --- SU2_CFD/src/solvers/CIncNSSolver.cpp | 2 ++ 1 file changed, 2 insertions(+) diff --git a/SU2_CFD/src/solvers/CIncNSSolver.cpp b/SU2_CFD/src/solvers/CIncNSSolver.cpp index ed944f70fc9d..38b0fe83a6fd 100644 --- a/SU2_CFD/src/solvers/CIncNSSolver.cpp +++ b/SU2_CFD/src/solvers/CIncNSSolver.cpp @@ -1265,6 +1265,8 @@ void CIncNSSolver::Friction_Forces(CGeometry *geometry, CConfig *config) { Grad_Temp[iDim] = nodes->GetGradient_Primitive(iPoint,nDim+1, iDim); } + /*--- TurbViscosity is usually zero on the walls, except when the SA roughness model is used. + * If the SA roughness model is not active, this term will be zero and will not affect the solution---*/ TurbViscosity = nodes->GetEddyViscosity(iPoint); Viscosity = nodes->GetLaminarViscosity(iPoint); Viscosity += TurbViscosity; From 456db04ebf3f9358092eff8e6e86dbf1f82f1dc8 Mon Sep 17 00:00:00 2001 From: koodlyakshay Date: Mon, 17 Feb 2020 22:45:14 +0100 Subject: [PATCH 08/29] Fix some issues in roughness values and change back to jPoint in SST BCHeatflux routine. --- SU2_CFD/src/solvers/CIncNSSolver.cpp | 3 ++- SU2_CFD/src/solvers/CNSSolver.cpp | 7 ++++- SU2_CFD/src/solvers/CTurbSASolver.cpp | 7 +++-- SU2_CFD/src/solvers/CTurbSSTSolver.cpp | 36 +++++++++++++------------- 4 files changed, 31 insertions(+), 22 deletions(-) diff --git a/SU2_CFD/src/solvers/CIncNSSolver.cpp b/SU2_CFD/src/solvers/CIncNSSolver.cpp index d1f1a2223928..0ca2ee470fe2 100644 --- a/SU2_CFD/src/solvers/CIncNSSolver.cpp +++ b/SU2_CFD/src/solvers/CIncNSSolver.cpp @@ -1267,7 +1267,8 @@ void CIncNSSolver::Friction_Forces(CGeometry *geometry, CConfig *config) { /*--- TurbViscosity is usually zero on the walls, except when the SA roughness model is used. * If the SA roughness model is not active, this term will be zero and will not affect the solution---*/ - TurbViscosity = nodes->GetEddyViscosity(iPoint); + TurbViscosity = 0.0; + if ((config->GetKindWall(Marker_Tag) == ROUGH ) && (config->GetKind_Turb_Model() == SA)) TurbViscosity = nodes->GetEddyViscosity(iPoint); Viscosity = nodes->GetLaminarViscosity(iPoint); Viscosity += TurbViscosity; diff --git a/SU2_CFD/src/solvers/CNSSolver.cpp b/SU2_CFD/src/solvers/CNSSolver.cpp index 96065b8dea4d..b421e13b0882 100644 --- a/SU2_CFD/src/solvers/CNSSolver.cpp +++ b/SU2_CFD/src/solvers/CNSSolver.cpp @@ -1408,7 +1408,7 @@ void CNSSolver::Friction_Forces(CGeometry *geometry, CConfig *config) { unsigned long iVertex, iPoint, iPointNormal; unsigned short Boundary, Monitoring, iMarker, iMarker_Monitoring, iDim, jDim; - su2double Viscosity = 0.0, div_vel, *Normal, MomentDist[3] = {0.0, 0.0, 0.0}, WallDist[3] = {0.0, 0.0, 0.0}, + su2double Viscosity = 0.0, TurbViscosity = 0.0, div_vel, *Normal, MomentDist[3] = {0.0, 0.0, 0.0}, WallDist[3] = {0.0, 0.0, 0.0}, *Coord, *Coord_Normal, Area, WallShearStress, TauNormal, factor, RefTemp, RefVel2, RefDensity, GradTemperature, Density = 0.0, WallDistMod, FrictionVel, Mach2Vel, Mach_Motion, UnitNormal[3] = {0.0, 0.0, 0.0}, TauElem[3] = {0.0, 0.0, 0.0}, TauTangent[3] = {0.0, 0.0, 0.0}, @@ -1529,7 +1529,12 @@ void CNSSolver::Friction_Forces(CGeometry *geometry, CConfig *config) { Grad_Temp[iDim] = nodes->GetGradient_Primitive(iPoint,0, iDim); } + /*--- TurbViscosity is usually zero on the walls, except when the SA roughness model is used. + * If the SA roughness model is not active, this term will be zero and will not affect the solution---*/ + TurbViscosity = 0.0; + if ((config->GetKindWall(Marker_Tag) == ROUGH ) && (config->GetKind_Turb_Model() == SA)) TurbViscosity = nodes->GetEddyViscosity(iPoint); Viscosity = nodes->GetLaminarViscosity(iPoint); + Viscosity += TurbViscosity; Density = nodes->GetDensity(iPoint); Area = 0.0; for (iDim = 0; iDim < nDim; iDim++) Area += Normal[iDim]*Normal[iDim]; Area = sqrt(Area); diff --git a/SU2_CFD/src/solvers/CTurbSASolver.cpp b/SU2_CFD/src/solvers/CTurbSASolver.cpp index 0a39eff37a10..502206fa9025 100644 --- a/SU2_CFD/src/solvers/CTurbSASolver.cpp +++ b/SU2_CFD/src/solvers/CTurbSASolver.cpp @@ -366,7 +366,10 @@ void CTurbSASolver::Postprocessing(CGeometry *geometry, CSolver **solver_contain dist = geometry->node[iPoint]->GetWall_Distance(); dist += 0.03*roughness; - Ji = nu_hat[0]/nu + cR1*roughness/(dist+EPS); + Ji = nu_hat[0]/nu ; + if (roughness > 1.0e-10) + Ji+= cR1*roughness/(dist+EPS); + Ji_3 = Ji*Ji*Ji; fv1 = Ji_3/(Ji_3+cv1_3); @@ -541,7 +544,7 @@ void CTurbSASolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_conta /*--- Includes 1 in the diagonal ---*/ Jacobian.DeleteValsRowi(iPoint); - } + } else { /*--- For rough walls, the boundary condition is given by * (\frac{\partial \nu}{\partial n})_wall = \frac{\nu}{0.03*k_s} diff --git a/SU2_CFD/src/solvers/CTurbSSTSolver.cpp b/SU2_CFD/src/solvers/CTurbSSTSolver.cpp index 54d921b6e06a..0d6d254b25e1 100644 --- a/SU2_CFD/src/solvers/CTurbSSTSolver.cpp +++ b/SU2_CFD/src/solvers/CTurbSSTSolver.cpp @@ -452,8 +452,7 @@ void CTurbSSTSolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_cont bool rough_wall = false; su2double RoughWallBC, Roughness_Height, S_R, FrictionVel, kPlus, WallShearStress; string Marker_Tag = config->GetMarker_All_TagBound(val_marker); - - + if (config->GetKindWall(Marker_Tag) == ROUGH ) rough_wall = true; for (iVertex = 0; iVertex < geometry->nVertex[val_marker]; iVertex++) { @@ -461,13 +460,13 @@ void CTurbSSTSolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_cont /*--- Check if the node belongs to the domain (i.e, not a halo node) ---*/ if (geometry->node[iPoint]->GetDomain()) { - + + if (rough_wall) { + /*--- Set wall values ---*/ - - density = solver_container[FLOW_SOL]->GetNodes()->GetDensity(iPoint); - laminar_viscosity = solver_container[FLOW_SOL]->GetNodes()->GetLaminarViscosity(iPoint); - - if (rough_wall) { + + density = solver_container[FLOW_SOL]->GetNodes()->GetDensity(iPoint); + laminar_viscosity = solver_container[FLOW_SOL]->GetNodes()->GetLaminarViscosity(iPoint); WallShearStress = solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 0)*solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 0); WallShearStress += solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 1)*solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 1); @@ -494,9 +493,8 @@ void CTurbSSTSolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_cont S_R = 100/(kPlus+EPS) + ((200/(kPlus+EPS))*(200/(kPlus+EPS)) - 100/(kPlus+EPS))*exp(5-kPlus); /*--- Modify the omega to account for a rough wall. ---*/ - Solution[1] = FrictionVel*FrictionVel*S_R/(laminar_viscosity/density); - - } else { + Solution[1] = FrictionVel*FrictionVel*S_R/(laminar_viscosity/density); + } else { /*--- distance to closest neighbor ---*/ jPoint = geometry->vertex[val_marker][iVertex]->GetNormal_Neighbor(); distance = 0.0; @@ -505,23 +503,25 @@ void CTurbSSTSolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_cont (geometry->node[iPoint]->GetCoord(iDim) - geometry->node[jPoint]->GetCoord(iDim)); } distance = sqrt(distance); - + /*--- Set wall values ---*/ - + + density = solver_container[FLOW_SOL]->GetNodes()->GetDensity(jPoint); + laminar_viscosity = solver_container[FLOW_SOL]->GetNodes()->GetLaminarViscosity(jPoint); + beta_1 = constants[4]; - - + Solution[1] = 60.0*laminar_viscosity/(density*beta_1*distance*distance); } - + /*--- Only the \omega solution is modified at the wall. ---*/ Solution[0] = 0.0; - + /*--- Set the solution values and zero the residual ---*/ nodes->SetSolution_Old(iPoint,Solution); nodes->SetSolution(iPoint,Solution); LinSysRes.SetBlock_Zero(iPoint); - + /*--- Change rows of the Jacobian (includes 1 in the diagonal) ---*/ for (iVar = 0; iVar < nVar; iVar++) { total_index = iPoint*nVar+iVar; From 5a9ad66b59da95b6de495c9dbe8400aa5506ba6d Mon Sep 17 00:00:00 2001 From: koodlyakshay Date: Tue, 18 Feb 2020 10:28:37 +0100 Subject: [PATCH 09/29] Fix nWallTypes and initialize Roughness in constructor. --- Common/src/CConfig.cpp | 5 +++-- Common/src/geometry/dual_grid/CPoint.cpp | 4 ++++ SU2_CFD/src/solvers/CTurbSASolver.cpp | 5 ++++- 3 files changed, 11 insertions(+), 3 deletions(-) diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index 6367c97ed936..cc38eeabaa42 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -4646,17 +4646,18 @@ void CConfig::SetPostprocessing(unsigned short val_software, unsigned short val_ if (nWall_Types != nMarker_HeatFlux + nMarker_Isothermal) { /*--- If this option is not used at all, assume all walls are smooth and set roughness height to zero. ---*/ if (nWall_Types == 0) { + nWall_Types = nMarker_HeatFlux + nMarker_Isothermal; if (Roughness_Height != NULL) Roughness_Height = NULL; Roughness_Height = new su2double [nWall_Types]; if (Kind_Wall != NULL) Kind_Wall = NULL; Kind_Wall = new unsigned short [nWall_Types]; for (unsigned short iMarker = 0; iMarker < nMarker_HeatFlux; iMarker++) { Roughness_Height[iMarker] = 0.0; - Kind_Wall[iMarker] = 1; + Kind_Wall[iMarker] = SMOOTH; } for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) { Roughness_Height[nMarker_HeatFlux + iMarker] = 0.0; - Kind_Wall[nMarker_HeatFlux + iMarker] = 1; + Kind_Wall[nMarker_HeatFlux + iMarker] = SMOOTH; } } else diff --git a/Common/src/geometry/dual_grid/CPoint.cpp b/Common/src/geometry/dual_grid/CPoint.cpp index 2b5c12c94da0..1ebd0714fdf2 100644 --- a/Common/src/geometry/dual_grid/CPoint.cpp +++ b/Common/src/geometry/dual_grid/CPoint.cpp @@ -140,6 +140,10 @@ CPoint::CPoint(unsigned short val_nDim, unsigned long val_globalindex, CConfig * /*--- Init walldistance ---*/ Wall_Distance = 0.0; + + /*--- Init walldistance ---*/ + + RoughnessHeight = 0.0; } CPoint::CPoint(su2double val_coord_0, su2double val_coord_1, unsigned long val_globalindex, CConfig *config) : CDualGrid(2) { diff --git a/SU2_CFD/src/solvers/CTurbSASolver.cpp b/SU2_CFD/src/solvers/CTurbSASolver.cpp index 502206fa9025..1d79d6e95aa2 100644 --- a/SU2_CFD/src/solvers/CTurbSASolver.cpp +++ b/SU2_CFD/src/solvers/CTurbSASolver.cpp @@ -349,6 +349,7 @@ void CTurbSASolver::Postprocessing(CGeometry *geometry, CSolver **solver_contain su2double rho = 0.0, mu = 0.0, nu, *nu_hat, muT, Ji, Ji_3, fv1; su2double cv1_3 = 7.1*7.1*7.1, cR1 = 0.5, roughness, dist; unsigned long iPoint; + unsigned short rank = SU2_MPI::GetRank(); bool neg_spalart_allmaras = (config->GetKind_Turb_Model() == SA_NEG); @@ -367,8 +368,10 @@ void CTurbSASolver::Postprocessing(CGeometry *geometry, CSolver **solver_contain dist += 0.03*roughness; Ji = nu_hat[0]/nu ; - if (roughness > 1.0e-10) + + if (roughness > 1.0e-10) { Ji+= cR1*roughness/(dist+EPS); + } Ji_3 = Ji*Ji*Ji; fv1 = Ji_3/(Ji_3+cv1_3); From 529772a17632608a2a118fb5ef6a4d844178a665 Mon Sep 17 00:00:00 2001 From: koodlyakshay Date: Wed, 26 Feb 2020 10:42:25 +0100 Subject: [PATCH 10/29] Fix tabs and spaces. --- Common/include/CConfig.hpp | 16 +-- Common/include/geometry/dual_grid/CPoint.hpp | 2 +- Common/src/CConfig.cpp | 112 +++++++-------- Common/src/geometry/CPhysicalGeometry.cpp | 142 +++++++++---------- Common/src/geometry/dual_grid/CPoint.cpp | 2 +- SU2_CFD/src/solvers/CTurbSASolver.cpp | 78 +++++----- SU2_CFD/src/solvers/CTurbSSTSolver.cpp | 20 +-- SU2_CFD/src/solvers/CTurbSolver.cpp | 2 +- 8 files changed, 186 insertions(+), 188 deletions(-) diff --git a/Common/include/CConfig.hpp b/Common/include/CConfig.hpp index 2d86d163bc97..2799f020b7c8 100644 --- a/Common/include/CConfig.hpp +++ b/Common/include/CConfig.hpp @@ -928,9 +928,9 @@ class CConfig { nMarkerTranslation, /*!< \brief Number of values provided for translational velocity of marker. */ nMarkerRotation_Rate, /*!< \brief Number of values provided for angular velocity of marker. */ nMarkerPitching_Omega, /*!< \brief Number of values provided for angular frequency of marker. */ - nMarkerPitching_Ampl, /*!< \brief Number of values provided for pitching amplitude of marker. */ - nMarkerPitching_Phase, /*!< \brief Number of values provided for pitching phase offset of marker. */ - nMarkerPlunging_Omega, /*!< \brief Number of values provided for angular frequency of marker. */ + nMarkerPitching_Ampl, /*!< \brief Number of values provided for pitching amplitude of marker. */ + nMarkerPitching_Phase, /*!< \brief Number of values provided for pitching phase offset of marker. */ + nMarkerPlunging_Omega, /*!< \brief Number of values provided for angular frequency of marker. */ nMarkerPlunging_Ampl, /*!< \brief Number of values provided for plunging amplitude of marker. */ nRoughWall; /*!< \brief Number of rough walls. */ su2double *Omega_HB; /*!< \brief Frequency for Harmonic Balance Operator (in rad/s). */ @@ -4940,12 +4940,12 @@ class CConfig { * \return Kind of wall - smooth or rough. */ unsigned short GetKindWall(string val_marker); - + /*! * \brief Set the kind of wall - rough or smooth. */ void SetKindWall(string val_marker, unsigned short val_kindwall); - + /*! * \brief Get the number of sections. * \return Number of sections @@ -6654,21 +6654,21 @@ class CConfig { * \return Pointer to the double info for the given marker. */ su2double* GetWallFunction_DoubleInfo(string val_marker); - + /*! * \brief Get the wall roughness height on a wall boundary (Heatflux or Isothermal). * \param[in] val_index - Index corresponding to the boundary. * \return The wall roughness height. */ su2double GetWall_RoughnessHeight(string val_marker); - + /*! * \brief Set the wall roughness height on a wall boundary (Heatflux or Isothermal). * \param[in] val_index - Index corresponding to the boundary. * \return The wall roughness height. */ void SetWall_RoughnessHeight(string val_marker, su2double val_roughness); - + /*! * \brief Get the target (pressure, massflow, etc) at an engine inflow boundary. * \param[in] val_index - Index corresponding to the engine inflow boundary. diff --git a/Common/include/geometry/dual_grid/CPoint.hpp b/Common/include/geometry/dual_grid/CPoint.hpp index a8d19dbdfa23..17dd311df2fe 100644 --- a/Common/include/geometry/dual_grid/CPoint.hpp +++ b/Common/include/geometry/dual_grid/CPoint.hpp @@ -168,7 +168,7 @@ class CPoint final : public CDualGrid { * \return Value of the roughness at the nearest wall. */ inline su2double GetRoughnessHeight() const { return RoughnessHeight; } - + /*! * \brief Set the number of elements that compose the control volume. * \param[in] val_nElem - Number of elements that make the control volume around a node. diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index 884c495b11aa..c307e721a66a 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -2531,7 +2531,7 @@ void CConfig::SetConfig_Options() { /*!\par INLETINTERPOLATION \n * DESCRIPTION: Type of spanwise interpolation to use for the inlet face. \n OPTIONS: see \link Inlet_SpanwiseInterpolation_Map \endlink * Sets Kind_InletInterpolation \ingroup Config - */ + */ addEnumOption("INLET_INTERPOLATION_FUNCTION",Kind_InletInterpolationFunction, Inlet_SpanwiseInterpolation_Map, NO_INTERPOLATION); /*!\par INLETINTERPOLATION \n @@ -4657,7 +4657,7 @@ void CConfig::SetPostprocessing(unsigned short val_software, unsigned short val_ SU2_MPI::Error("Must list two markers for the pressure drop objective function.\n Expected format: MARKER_ANALYZE= (outlet_name, inlet_name).", CURRENT_FUNCTION); } } - /*--- Wall roughness handling begins here ---*/ + /*--- Check that if the wall roughness array are compatible and set deafult values if needed. ---*/ if ((nMarker_HeatFlux > 0) || (nMarker_Isothermal > 0)) { @@ -4665,47 +4665,45 @@ void CConfig::SetPostprocessing(unsigned short val_software, unsigned short val_ if (nWall_Types != nMarker_HeatFlux + nMarker_Isothermal) { /*--- If this option is not used at all, assume all walls are smooth and set roughness height to zero. ---*/ if (nWall_Types == 0) { - nWall_Types = nMarker_HeatFlux + nMarker_Isothermal; - if (Roughness_Height != NULL) Roughness_Height = NULL; - Roughness_Height = new su2double [nWall_Types]; - if (Kind_Wall != NULL) Kind_Wall = NULL; - Kind_Wall = new unsigned short [nWall_Types]; - for (unsigned short iMarker = 0; iMarker < nMarker_HeatFlux; iMarker++) { - Roughness_Height[iMarker] = 0.0; - Kind_Wall[iMarker] = SMOOTH; - } - for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) { - Roughness_Height[nMarker_HeatFlux + iMarker] = 0.0; - Kind_Wall[nMarker_HeatFlux + iMarker] = SMOOTH; - } - } - else - SU2_MPI::Error("Mismatch in number of wall markers. Specify type of wall for all markers (even if smooth).", CURRENT_FUNCTION); + nWall_Types = nMarker_HeatFlux + nMarker_Isothermal; + if (Roughness_Height != NULL) Roughness_Height = NULL; + Roughness_Height = new su2double [nWall_Types]; + if (Kind_Wall != NULL) Kind_Wall = NULL; + Kind_Wall = new unsigned short [nWall_Types]; + for (unsigned short iMarker = 0; iMarker < nMarker_HeatFlux; iMarker++) { + Roughness_Height[iMarker] = 0.0; + Kind_Wall[iMarker] = SMOOTH; + } + for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) { + Roughness_Height[nMarker_HeatFlux + iMarker] = 0.0; + Kind_Wall[nMarker_HeatFlux + iMarker] = SMOOTH; + } + } + else + SU2_MPI::Error("Mismatch in number of wall markers. Specify type of wall for all markers (even if smooth).", CURRENT_FUNCTION); } else { - if (nRoughWall != nWall_Types) - SU2_MPI::Error("Mismatch in number of roughness height definition and wall type definition.", CURRENT_FUNCTION); - else { - su2double *temp_rough; - unsigned short *temp_kindrough; - temp_rough = new su2double [nWall_Types]; - temp_kindrough = new unsigned short [nWall_Types]; - for (unsigned short iMarker = 0; iMarker < nWall_Types; iMarker++) { - temp_rough[iMarker] = Roughness_Height[iMarker]; - temp_kindrough[iMarker] = Kind_Wall[iMarker]; - } - for (iMarker = 0; iMarker < nMarker_HeatFlux; iMarker++) { - Roughness_Height[iMarker] = temp_rough[iMarker]; - Kind_Wall[iMarker] = temp_kindrough[iMarker]; - } - for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) { - Roughness_Height[nMarker_HeatFlux + iMarker] = temp_rough[nMarker_HeatFlux + iMarker]; - Kind_Wall[nMarker_HeatFlux + iMarker] = temp_kindrough[nMarker_HeatFlux + iMarker]; - } - } - } - - } - /*--- Wall roughness handling ends here. ---*/ + if (nRoughWall != nWall_Types) + SU2_MPI::Error("Mismatch in number of roughness height definition and wall type definition.", CURRENT_FUNCTION); + else { + su2double *temp_rough; + unsigned short *temp_kindrough; + temp_rough = new su2double [nWall_Types]; + temp_kindrough = new unsigned short [nWall_Types]; + for (unsigned short iMarker = 0; iMarker < nWall_Types; iMarker++) { + temp_rough[iMarker] = Roughness_Height[iMarker]; + temp_kindrough[iMarker] = Kind_Wall[iMarker]; + } + for (iMarker = 0; iMarker < nMarker_HeatFlux; iMarker++) { + Roughness_Height[iMarker] = temp_rough[iMarker]; + Kind_Wall[iMarker] = temp_kindrough[iMarker]; + } + for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) { + Roughness_Height[nMarker_HeatFlux + iMarker] = temp_rough[nMarker_HeatFlux + iMarker]; + Kind_Wall[nMarker_HeatFlux + iMarker] = temp_kindrough[nMarker_HeatFlux + iMarker]; + } + } + } + } /*--- Handle default options for topology optimization ---*/ @@ -8910,21 +8908,21 @@ unsigned short CConfig::GetKindWall(string val_marker) { short flag = -1; for (iMarker = 0; iMarker < nMarker_HeatFlux; iMarker++) if (Marker_HeatFlux[iMarker] == val_marker) { - flag = iMarker; - break; - } + flag = iMarker; + break; + } if (flag == -1) for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) if (Marker_Isothermal[iMarker] == val_marker) { - flag = iMarker; - break; - } + flag = iMarker; + break; + } if (flag == -1) return 1; else return Kind_Wall[flag]; } void CConfig::SetWall_RoughnessHeight(string val_marker, su2double val_roughness) { - //Roughness_Height[val_marker] = val_roughness; + //Roughness_Height[val_marker] = val_roughness; } su2double CConfig::GetWall_RoughnessHeight(string val_marker) { @@ -8934,16 +8932,16 @@ su2double CConfig::GetWall_RoughnessHeight(string val_marker) { if (nMarker_HeatFlux > 0 || nMarker_Isothermal > 0) { for (iMarker = 0; iMarker < nMarker_HeatFlux; iMarker++) if (Marker_HeatFlux[iMarker] == val_marker) { - flag = iMarker; - break; - } + flag = iMarker; + break; + } if (flag == -1) { - for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) + for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) if (Marker_Isothermal[iMarker] == val_marker) { - flag = iMarker; - break; - } - } + flag = iMarker; + break; + } + } } return Roughness_Height[flag]; } diff --git a/Common/src/geometry/CPhysicalGeometry.cpp b/Common/src/geometry/CPhysicalGeometry.cpp index 55db0f542699..14fc2317744e 100644 --- a/Common/src/geometry/CPhysicalGeometry.cpp +++ b/Common/src/geometry/CPhysicalGeometry.cpp @@ -6,7 +6,7 @@ * * SU2 Project Website: https://su2code.github.io * - * The SU2 Project is maintained by the SU2 Foundation + * The SU2 Project is maintained by the SU2 Foundation * (http://su2foundation.org) * * Copyright 2012-2019, SU2 Contributors (cf. AUTHORS.md) @@ -257,7 +257,7 @@ CPhysicalGeometry::CPhysicalGeometry(CConfig *config, unsigned short val_iZone, if (bound[iMarker][iElem_Bound]->GetVTK_Type() == VERTEX) { boundary_file << bound[iMarker][iElem_Bound]->GetRotation_Type() << "\t"; } - boundary_file << iElem_Bound << endl; + boundary_file << iElem_Bound << endl; } } @@ -270,7 +270,7 @@ CPhysicalGeometry::CPhysicalGeometry(CConfig *config, unsigned short val_iZone, if (bound[iMarker][iElem_Bound]->GetVTK_Type() == VERTEX) { boundary_file << bound[iMarker][iElem_Bound]->GetRotation_Type() << "\t"; } - boundary_file << iElem_Bound << endl; + boundary_file << iElem_Bound << endl; } } @@ -3071,7 +3071,7 @@ void CPhysicalGeometry::LoadSurfaceElements(CConfig *config, CGeometry *geometry /*--- Initialize pointers for turbomachinery computations ---*/ nSpanWiseSections = new unsigned short[2]; - nSpanSectionsByMarker = new unsigned short[nMarker]; + nSpanSectionsByMarker = new unsigned short[nMarker]; SpanWiseValue = new su2double*[2]; for (unsigned short iMarker = 0; iMarker < 2; iMarker++){ nSpanWiseSections[iMarker] = 0; @@ -3092,7 +3092,7 @@ void CPhysicalGeometry::LoadSurfaceElements(CConfig *config, CGeometry *geometry MinRelAngularCoord = new su2double*[nMarker]; for (unsigned short iMarker = 0; iMarker < nMarker; iMarker++){ - nSpanSectionsByMarker[iMarker] = 0; + nSpanSectionsByMarker[iMarker] = 0; nVertexSpan[iMarker] = NULL; nTotVertexSpan[iMarker] = NULL; turbovertex[iMarker] = NULL; @@ -3118,8 +3118,8 @@ void CPhysicalGeometry::LoadSurfaceElements(CConfig *config, CGeometry *geometry TurboRadiusOut = new su2double*[config->GetnMarker_TurboPerformance()]; for (unsigned short iMarker = 0; iMarker < config->GetnMarker_TurboPerformance(); iMarker++){ - TangGridVelIn[iMarker] = NULL; - SpanAreaIn[iMarker] = NULL; + TangGridVelIn[iMarker] = NULL; + SpanAreaIn[iMarker] = NULL; TurboRadiusIn[iMarker] = NULL; TangGridVelOut[iMarker] = NULL; SpanAreaOut[iMarker] = NULL; @@ -3596,7 +3596,7 @@ void CPhysicalGeometry::SetSendReceive(CConfig *config) { } /*--- First compute the Send/Receive boundaries ---*/ - Counter_Send = 0; Counter_Receive = 0; + Counter_Send = 0; Counter_Receive = 0; for (iDomain = 0; iDomain < nDomain; iDomain++) if (SendDomainLocal[iDomain].size() != 0) Counter_Send++; @@ -4667,8 +4667,8 @@ void CPhysicalGeometry::Check_IntElem_Orientation(CConfig *config) { test = a[0]*b[1]-b[0]*a[1]; if (test < 0.0) { - elem[iElem]->Change_Orientation(); - triangle_flip++; + elem[iElem]->Change_Orientation(); + triangle_flip++; } } @@ -4726,8 +4726,8 @@ void CPhysicalGeometry::Check_IntElem_Orientation(CConfig *config) { test = n[0]*c[0]+n[1]*c[1]+n[2]*c[2]; if (test < 0.0) { - elem[iElem]->Change_Orientation(); - tet_flip++; + elem[iElem]->Change_Orientation(); + tet_flip++; } } @@ -4844,8 +4844,8 @@ void CPhysicalGeometry::Check_IntElem_Orientation(CConfig *config) { if ((test_1 < 0.0) || (test_2 < 0.0) || (test_3 < 0.0) || (test_4 < 0.0)) { - elem[iElem]->Change_Orientation(); - hexa_flip++; + elem[iElem]->Change_Orientation(); + hexa_flip++; } } @@ -4884,7 +4884,7 @@ void CPhysicalGeometry::Check_IntElem_Orientation(CConfig *config) { if ((test_1 < 0.0) || (test_2 < 0.0)) { elem[iElem]->Change_Orientation(); - pyram_flip++; + pyram_flip++; } } @@ -5109,7 +5109,7 @@ void CPhysicalGeometry::ComputeWall_Distance(CConfig *config) { /* Check for a viscous wall. */ if( config->GetViscous_Wall(iMarker)) { - + Marker_Tag = config->GetMarker_All_TagBound(iMarker); roughness = config->GetWall_RoughnessHeight(Marker_Tag); @@ -5133,7 +5133,7 @@ void CPhysicalGeometry::ComputeWall_Distance(CConfig *config) { VTK_TypeElem.push_back(VTK_Type); elemIDs.push_back(iElem); roughheights.push_back(roughness); - + for (unsigned short iNode = 0; iNode < nDOFsPerElem; iNode++) surfaceConn.push_back(bound[iMarker][iElem]->GetNode(iNode)); } @@ -5209,13 +5209,13 @@ void CPhysicalGeometry::ComputeWall_Distance(CConfig *config) { WallADT.DetermineNearestElement(node[iPoint]->GetCoord(), dist, markerID, elemID, rankID, localRoughness); - + node[iPoint]->SetWall_Distance(dist); - + /*--- Use the markerID to find the corresponding wall roughness height. ---*/ - - Marker_Tag = config->GetMarker_All_TagBound(markerID); - roughness = config->GetWall_RoughnessHeight(Marker_Tag); + + Marker_Tag = config->GetMarker_All_TagBound(markerID); + roughness = config->GetWall_RoughnessHeight(Marker_Tag); node[iPoint]->SetRoughnessHeight(localRoughness); @@ -5376,7 +5376,7 @@ void CPhysicalGeometry::SetPositive_ZArea(CConfig *config) { if (config->GetSystemMeasurements() == SI) cout <<" m"; else cout <<" ft"; if (nDim == 3) { - cout << ", z-direction = "<< TotalMaxCoordZ; + cout << ", z-direction = "<< TotalMaxCoordZ; if (config->GetSystemMeasurements() == SI) cout <<" m." << endl; else cout <<" ft."<< endl; } else cout << "." << endl; @@ -5388,7 +5388,7 @@ void CPhysicalGeometry::SetPositive_ZArea(CConfig *config) { if (config->GetSystemMeasurements() == SI) cout <<" m"; else cout <<" ft"; if (nDim == 3) { - cout << ", z-direction = "<< TotalMinCoordZ; + cout << ", z-direction = "<< TotalMinCoordZ; if (config->GetSystemMeasurements() == SI) cout <<" m." << endl; else cout <<" ft."<< endl; } else cout << "." << endl; @@ -5796,7 +5796,7 @@ void CPhysicalGeometry::ComputeNSpan(CConfig *config, unsigned short val_iZone, if (nDim == 2){ nSpanWiseSections[marker_flag-1] = 1; //TODO (turbo) make it more genral - if(marker_flag == OUTFLOW) config->SetnSpanWiseSections(1); + if(marker_flag == OUTFLOW) config->SetnSpanWiseSections(1); /*---Initilize the vector of span-wise values that will be ordered ---*/ SpanWiseValue[marker_flag -1] = new su2double[1]; @@ -5862,7 +5862,7 @@ void CPhysicalGeometry::ComputeNSpan(CConfig *config, unsigned short val_iZone, for (jMarker = 0; jMarker < nMarker; jMarker++){ if (config->GetMarker_All_KindBC(jMarker) == PERIODIC_BOUNDARY) { PeriodicBoundary = config->GetMarker_All_PerBound(jMarker); - jVertex = node[iPoint]->GetVertex(jMarker); + jVertex = node[iPoint]->GetVertex(jMarker); if ((jVertex != -1) && (PeriodicBoundary == (val_iZone + 1))){ coord = node[iPoint]->GetCoord(); switch (config->GetKind_TurboMachinery(val_iZone)){ @@ -5942,9 +5942,9 @@ void CPhysicalGeometry::ComputeNSpan(CConfig *config, unsigned short val_iZone, // check if the value are gathered correctly // // for (iSpan = 0; iSpan < nSpan; iSpan++){ - // if(rank == MASTER_NODE){ - // cout << setprecision(16)<< iSpan +1 << " with a value of " <GetMarker_All_Turbomachinery(iMarker) == iMarkerTP){ if (config->GetMarker_All_TurbomachineryFlag(iMarker) == marker_flag){ turboNormal2 = 0.0; - Normal2 = 0.0; + Normal2 = 0.0; for (iDim = 0; iDim < nDim; iDim++){ turboNormal2 += AverageTurboNormal[iMarker][nSpanWiseSections[marker_flag-1]][iDim]*AverageTurboNormal[iMarker][nSpanWiseSections[marker_flag-1]][iDim]; @@ -8003,11 +8003,11 @@ void CPhysicalGeometry::MatchPeriodic(CConfig *config, su2double rotMatrix[3][3] = {{1.0,0.0,0.0},{0.0,1.0,0.0},{0.0,0.0,1.0}}; su2double Theta, Phi, Psi, cosTheta, sinTheta, cosPhi, sinPhi, cosPsi, sinPsi; su2double rotCoord[3] = {0.0, 0.0, 0.0}; - + bool pointOnAxis = false; - + bool chkSamePoint = false; - + su2double distToAxis = 0.0; /*--- Tolerance for distance-based match to report warning. ---*/ @@ -8217,7 +8217,7 @@ void CPhysicalGeometry::MatchPeriodic(CConfig *config, rotMatrix[2][1]*dy + rotMatrix[2][2]*dz + translation[2]); - /*--- Check if the point lies on the axis of rotation. If it does, + /*--- Check if the point lies on the axis of rotation. If it does, the rotated coordinate and the original coordinate are the same. ---*/ pointOnAxis = false; @@ -8272,7 +8272,7 @@ void CPhysicalGeometry::MatchPeriodic(CConfig *config, and also perform checks just to be sure that this is an independent periodic point (even if on the same rank), unless it lies on the axis of rotation. ---*/ - + chkSamePoint = false; chkSamePoint = (((dist < mindist) && (iProcessor != rank)) || ((dist < mindist) && (iProcessor == rank) && @@ -9606,7 +9606,7 @@ void CPhysicalGeometry::SetBoundSensitivity(CConfig *config) { /*--- Write file name with extension if unsteady or steady ---*/ if (config->GetTime_Marching() == HARMONIC_BALANCE) - SPRINTF (buffer, "_%d.csv", SU2_TYPE::Int(iTimeIter)); + SPRINTF (buffer, "_%d.csv", SU2_TYPE::Int(iTimeIter)); if ((config->GetTime_Marching() && config->GetTime_Domain()) || (config->GetTime_Marching() == HARMONIC_BALANCE)) { @@ -10430,7 +10430,7 @@ su2double CPhysicalGeometry::Compute_MaxThickness(su2double *Plane_P0, su2double unsigned short index = 2; /*--- Removing the trailing edge from list of points that we are going to use in the interpolation, - to be sure that a blunt trailing edge do not affect the interpolation ---*/ + to be sure that a blunt trailing edge do not affect the interpolation ---*/ if ((Normal[index] >= 0.0) && (fabs(Xcoord_Airfoil_[iVertex]) > MaxDistance*0.01)) { Xcoord.push_back(Xcoord_Airfoil_[iVertex]); @@ -10656,8 +10656,8 @@ su2double CPhysicalGeometry::Compute_WaterLineWidth(su2double *Plane_P0, su2doub for (iVertex = 0; iVertex < Xcoord_Airfoil.size(); iVertex++) { Distance = fabs(Zcoord_Airfoil[iVertex] - WaterLine); if (Distance < MinDistance) { - MinDistance = Distance; - WaterLineWidth = fabs(Xcoord_Airfoil[iVertex] - Xcoord_Airfoil[Trailing_Point]); + MinDistance = Distance; + WaterLineWidth = fabs(Xcoord_Airfoil[iVertex] - Xcoord_Airfoil[Trailing_Point]); } } @@ -11116,21 +11116,21 @@ void CPhysicalGeometry::Compute_Wing(CConfig *config, bool original_surface, /*--- Plot the geometrical quatities ---*/ for (iPlane = 0; iPlane < nPlane; iPlane++) { - if (Xcoord_Airfoil[iPlane].size() > 1) { - if (config->GetTabular_FileFormat() == TAB_CSV) { - Wing_File << Ycoord_Airfoil[iPlane][0]/SemiSpan <<", "<< Area[iPlane] <<", "<< MaxThickness[iPlane] <<", "<< Chord[iPlane] <<", "<< LERadius[iPlane] <<", "<< ToC[iPlane] - <<", "<< Twist[iPlane] <<", "<< Curvature[iPlane] <<", "<< Dihedral[iPlane] - <<", "<< LeadingEdge[iPlane][0]/SemiSpan <<", "<< LeadingEdge[iPlane][2]/SemiSpan - <<", "<< TrailingEdge[iPlane][0]/SemiSpan <<", "<< TrailingEdge[iPlane][2]/SemiSpan << endl; - } - else { - Wing_File << Ycoord_Airfoil[iPlane][0]/SemiSpan <<" "<< Area[iPlane] <<" "<< MaxThickness[iPlane] <<" "<< Chord[iPlane] <<" "<< LERadius[iPlane] <<" "<< ToC[iPlane] - <<" "<< Twist[iPlane] <<" "<< Curvature[iPlane] <<" "<< Dihedral[iPlane] - <<" "<< LeadingEdge[iPlane][0]/SemiSpan <<" "<< LeadingEdge[iPlane][2]/SemiSpan - <<" "<< TrailingEdge[iPlane][0]/SemiSpan <<" "<< TrailingEdge[iPlane][2]/SemiSpan << endl; - - } - } + if (Xcoord_Airfoil[iPlane].size() > 1) { + if (config->GetTabular_FileFormat() == TAB_CSV) { + Wing_File << Ycoord_Airfoil[iPlane][0]/SemiSpan <<", "<< Area[iPlane] <<", "<< MaxThickness[iPlane] <<", "<< Chord[iPlane] <<", "<< LERadius[iPlane] <<", "<< ToC[iPlane] + <<", "<< Twist[iPlane] <<", "<< Curvature[iPlane] <<", "<< Dihedral[iPlane] + <<", "<< LeadingEdge[iPlane][0]/SemiSpan <<", "<< LeadingEdge[iPlane][2]/SemiSpan + <<", "<< TrailingEdge[iPlane][0]/SemiSpan <<", "<< TrailingEdge[iPlane][2]/SemiSpan << endl; + } + else { + Wing_File << Ycoord_Airfoil[iPlane][0]/SemiSpan <<" "<< Area[iPlane] <<" "<< MaxThickness[iPlane] <<" "<< Chord[iPlane] <<" "<< LERadius[iPlane] <<" "<< ToC[iPlane] + <<" "<< Twist[iPlane] <<" "<< Curvature[iPlane] <<" "<< Dihedral[iPlane] + <<" "<< LeadingEdge[iPlane][0]/SemiSpan <<" "<< LeadingEdge[iPlane][2]/SemiSpan + <<" "<< TrailingEdge[iPlane][0]/SemiSpan <<" "<< TrailingEdge[iPlane][2]/SemiSpan << endl; + + } + } } Wing_File.close(); @@ -11402,18 +11402,18 @@ void CPhysicalGeometry::Compute_Fuselage(CConfig *config, bool original_surface, /*--- Plot the geometrical quatities ---*/ for (iPlane = 0; iPlane < nPlane; iPlane++) { - if (Xcoord_Airfoil[iPlane].size() > 1) { - if (config->GetTabular_FileFormat() == TAB_CSV) { - Fuselage_File << -Ycoord_Airfoil[iPlane][0] <<", "<< Area[iPlane] <<", "<< Length[iPlane] <<", "<< Width[iPlane] <<", "<< WaterLineWidth[iPlane] <<", "<< Height[iPlane] <<", "<< Curvature[iPlane] - <<", "<< -LeadingEdge[iPlane][1] <<", "<< LeadingEdge[iPlane][0] <<", "<< LeadingEdge[iPlane][2] - <<", "<< -TrailingEdge[iPlane][1] <<", "<< TrailingEdge[iPlane][0] <<", "<< TrailingEdge[iPlane][2] << endl; - } - else { - Fuselage_File << -Ycoord_Airfoil[iPlane][0] <<" "<< Area[iPlane] <<" "<< Length[iPlane] <<" "<< Width[iPlane] <<" "<< WaterLineWidth[iPlane] <<" "<< Height[iPlane] <<" "<< Curvature[iPlane] - <<" "<< -LeadingEdge[iPlane][1] <<" "<< LeadingEdge[iPlane][0] <<" "<< LeadingEdge[iPlane][2] - <<" "<< -TrailingEdge[iPlane][1] <<" "<< TrailingEdge[iPlane][0] <<" "<< TrailingEdge[iPlane][2] << endl; - } - } + if (Xcoord_Airfoil[iPlane].size() > 1) { + if (config->GetTabular_FileFormat() == TAB_CSV) { + Fuselage_File << -Ycoord_Airfoil[iPlane][0] <<", "<< Area[iPlane] <<", "<< Length[iPlane] <<", "<< Width[iPlane] <<", "<< WaterLineWidth[iPlane] <<", "<< Height[iPlane] <<", "<< Curvature[iPlane] + <<", "<< -LeadingEdge[iPlane][1] <<", "<< LeadingEdge[iPlane][0] <<", "<< LeadingEdge[iPlane][2] + <<", "<< -TrailingEdge[iPlane][1] <<", "<< TrailingEdge[iPlane][0] <<", "<< TrailingEdge[iPlane][2] << endl; + } + else { + Fuselage_File << -Ycoord_Airfoil[iPlane][0] <<" "<< Area[iPlane] <<" "<< Length[iPlane] <<" "<< Width[iPlane] <<" "<< WaterLineWidth[iPlane] <<" "<< Height[iPlane] <<" "<< Curvature[iPlane] + <<" "<< -LeadingEdge[iPlane][1] <<" "<< LeadingEdge[iPlane][0] <<" "<< LeadingEdge[iPlane][2] + <<" "<< -TrailingEdge[iPlane][1] <<" "<< TrailingEdge[iPlane][0] <<" "<< TrailingEdge[iPlane][2] << endl; + } + } } Fuselage_File.close(); @@ -11451,7 +11451,7 @@ void CPhysicalGeometry::Compute_Fuselage(CConfig *config, bool original_surface, Fuselage_Volume = 0.0; for (iPlane = 0; iPlane < nPlane-2; iPlane+=2) { if (Xcoord_Airfoil[iPlane].size() > 1) { - Fuselage_Volume += (1.0/3.0)*dPlane*(Area[iPlane] + 4.0*Area[iPlane+1] + Area[iPlane+2]); + Fuselage_Volume += (1.0/3.0)*dPlane*(Area[iPlane] + 4.0*Area[iPlane+1] + Area[iPlane+2]); } } @@ -11461,7 +11461,7 @@ void CPhysicalGeometry::Compute_Fuselage(CConfig *config, bool original_surface, if (Xcoord_Airfoil[0].size() > 1) Fuselage_WettedArea += (1.0/2.0)*dPlane*Length[0]; for (iPlane = 1; iPlane < nPlane-1; iPlane++) { if (Xcoord_Airfoil[iPlane].size() > 1) { - Fuselage_WettedArea += dPlane*Length[iPlane]; + Fuselage_WettedArea += dPlane*Length[iPlane]; } } if (Xcoord_Airfoil[nPlane-1].size() > 1) Fuselage_WettedArea += (1.0/2.0)*dPlane*Length[nPlane-1]; diff --git a/Common/src/geometry/dual_grid/CPoint.cpp b/Common/src/geometry/dual_grid/CPoint.cpp index 1ebd0714fdf2..56544e928e00 100644 --- a/Common/src/geometry/dual_grid/CPoint.cpp +++ b/Common/src/geometry/dual_grid/CPoint.cpp @@ -140,7 +140,7 @@ CPoint::CPoint(unsigned short val_nDim, unsigned long val_globalindex, CConfig * /*--- Init walldistance ---*/ Wall_Distance = 0.0; - + /*--- Init walldistance ---*/ RoughnessHeight = 0.0; diff --git a/SU2_CFD/src/solvers/CTurbSASolver.cpp b/SU2_CFD/src/solvers/CTurbSASolver.cpp index 1d79d6e95aa2..8a9fc5780bb6 100644 --- a/SU2_CFD/src/solvers/CTurbSASolver.cpp +++ b/SU2_CFD/src/solvers/CTurbSASolver.cpp @@ -362,17 +362,17 @@ void CTurbSASolver::Postprocessing(CGeometry *geometry, CSolver **solver_contain nu = mu/rho; nu_hat = nodes->GetSolution(iPoint); - + roughness = geometry->node[iPoint]->GetRoughnessHeight(); dist = geometry->node[iPoint]->GetWall_Distance(); dist += 0.03*roughness; Ji = nu_hat[0]/nu ; - + if (roughness > 1.0e-10) { - Ji+= cR1*roughness/(dist+EPS); + Ji+= cR1*roughness/(dist+EPS); } - + Ji_3 = Ji*Ji*Ji; fv1 = Ji_3/(Ji_3+cv1_3); @@ -430,20 +430,20 @@ void CTurbSASolver::Source_Residual(CGeometry *geometry, CSolver **solver_contai /*--- Get Hybrid RANS/LES Type and set the appropriate wall distance ---*/ if (config->GetKind_HybridRANSLES() == NO_HYBRIDRANSLES) { - - /*--- For the SA model, wall roughness is accounted by modifying the computed wall distance + + /*--- For the SA model, wall roughness is accounted by modifying the computed wall distance * d_new = d + 0.03 k_s * where k_s is the equivalent sand grain roughness height that is specified in cfg file. * For smooth walls, wall roughness is zero and computed wall distance remains the same. */ - + modifiedWallDistance = geometry->node[iPoint]->GetWall_Distance(); - + modifiedWallDistance += 0.03*geometry->node[iPoint]->GetRoughnessHeight(); - + /*--- Set distance to the surface ---*/ - + numerics->SetDistance(modifiedWallDistance, 0.0); - + /*--- Set the roughness of the closest wall. ---*/ numerics->SetRoughness(geometry->node[iPoint]->GetRoughnessHeight(), 0.0 ); @@ -543,44 +543,44 @@ void CTurbSASolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_conta nodes->SetSolution_Old(iPoint,Solution); LinSysRes.SetBlock_Zero(iPoint); - + /*--- Includes 1 in the diagonal ---*/ Jacobian.DeleteValsRowi(iPoint); } - else { - /*--- For rough walls, the boundary condition is given by - * (\frac{\partial \nu}{\partial n})_wall = \frac{\nu}{0.03*k_s} - * where \nu is the solution variable, $n$ is the wall normal direction - * and k_s is the equivalent sand grain roughness specified. ---*/ - - /*--- Compute dual-grid area and boundary normal ---*/ - - Normal = geometry->vertex[val_marker][iVertex]->GetNormal(); - + else { + /*--- For rough walls, the boundary condition is given by + * (\frac{\partial \nu}{\partial n})_wall = \frac{\nu}{0.03*k_s} + * where \nu is the solution variable, $n$ is the wall normal direction + * and k_s is the equivalent sand grain roughness specified. ---*/ + + /*--- Compute dual-grid area and boundary normal ---*/ + + Normal = geometry->vertex[val_marker][iVertex]->GetNormal(); + Area = 0.0; - for (iDim = 0; iDim < nDim; iDim++) - Area += Normal[iDim]*Normal[iDim]; - Area = sqrt(Area); + for (iDim = 0; iDim < nDim; iDim++) + Area += Normal[iDim]*Normal[iDim]; + Area = sqrt(Area); + + /*--- Get laminar_viscosity and density ---*/ - /*--- Get laminar_viscosity and density ---*/ - - laminar_viscosity = solver_container[FLOW_SOL]->GetNodes()->GetLaminarViscosity(iPoint); + laminar_viscosity = solver_container[FLOW_SOL]->GetNodes()->GetLaminarViscosity(iPoint); density = solver_container[FLOW_SOL]->GetNodes()->GetDensity(iPoint); - + nu_total = (laminar_viscosity/density + nodes->GetSolution(iPoint,0)); - + coeff = (nu_total/sigma); - RoughWallBC = nodes->GetSolution(iPoint,0)/(0.03*Roughness_Height); + RoughWallBC = nodes->GetSolution(iPoint,0)/(0.03*Roughness_Height); - Res_Wall[0] = coeff*RoughWallBC*Area; - LinSysRes.SubtractBlock(iPoint, Res_Wall); - - Jacobian_i[0][0] = (laminar_viscosity*Area)/(0.03*Roughness_Height*sigma); - Jacobian_i[0][0] += 2.0*RoughWallBC*Area/sigma; - Jacobian.SubtractBlock(iPoint,iPoint,Jacobian_i); - } + Res_Wall[0] = coeff*RoughWallBC*Area; + LinSysRes.SubtractBlock(iPoint, Res_Wall); + + Jacobian_i[0][0] = (laminar_viscosity*Area)/(0.03*Roughness_Height*sigma); + Jacobian_i[0][0] += 2.0*RoughWallBC*Area/sigma; + Jacobian.SubtractBlock(iPoint,iPoint,Jacobian_i); + } } } } @@ -2422,5 +2422,5 @@ void CTurbSASolver::SetUniformInlet(CConfig* config, unsigned short iMarker) { for(unsigned long iVertex=0; iVertex < nVertex[iMarker]; iVertex++){ Inlet_TurbVars[iMarker][iVertex][0] = nu_tilde_Inf; } - + } diff --git a/SU2_CFD/src/solvers/CTurbSSTSolver.cpp b/SU2_CFD/src/solvers/CTurbSSTSolver.cpp index 0d6d254b25e1..5ca3a0584fe1 100644 --- a/SU2_CFD/src/solvers/CTurbSSTSolver.cpp +++ b/SU2_CFD/src/solvers/CTurbSSTSolver.cpp @@ -445,25 +445,25 @@ void CTurbSSTSolver::Source_Template(CGeometry *geometry, CSolver **solver_conta } void CTurbSSTSolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_container, CNumerics *conv_numerics, CNumerics *visc_numerics, CConfig *config, unsigned short val_marker) { - + unsigned long iPoint, jPoint, iVertex, total_index; unsigned short iDim, iVar; su2double distance, density = 0.0, laminar_viscosity = 0.0, beta_1; bool rough_wall = false; su2double RoughWallBC, Roughness_Height, S_R, FrictionVel, kPlus, WallShearStress; string Marker_Tag = config->GetMarker_All_TagBound(val_marker); - + if (config->GetKindWall(Marker_Tag) == ROUGH ) rough_wall = true; - + for (iVertex = 0; iVertex < geometry->nVertex[val_marker]; iVertex++) { iPoint = geometry->vertex[val_marker][iVertex]->GetNode(); - + /*--- Check if the node belongs to the domain (i.e, not a halo node) ---*/ if (geometry->node[iPoint]->GetDomain()) { if (rough_wall) { - /*--- Set wall values ---*/ + /*--- Set wall values ---*/ density = solver_container[FLOW_SOL]->GetNodes()->GetDensity(iPoint); laminar_viscosity = solver_container[FLOW_SOL]->GetNodes()->GetLaminarViscosity(iPoint); @@ -481,13 +481,13 @@ void CTurbSSTSolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_cont kPlus = FrictionVel*Roughness_Height*density/laminar_viscosity; /*--- Reference 1 original Wilcox (1998) ---*/ - /*if (kPlus <= 25) + /*if (kPlus <= 25) S_R = (50/(kPlus+EPS))*(50/(kPlus+EPS)); else S_R = 100/(kPlus+EPS);*/ /*--- Reference 2 from D.C. Wilcox Turbulence Modeling for CFD (2006) ---*/ - if (kPlus <= 5) + if (kPlus <= 5) S_R = (200/(kPlus+EPS))*(200/(kPlus+EPS)); else S_R = 100/(kPlus+EPS) + ((200/(kPlus+EPS))*(200/(kPlus+EPS)) - 100/(kPlus+EPS))*exp(5-kPlus); @@ -505,14 +505,14 @@ void CTurbSSTSolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_cont distance = sqrt(distance); /*--- Set wall values ---*/ - + density = solver_container[FLOW_SOL]->GetNodes()->GetDensity(jPoint); laminar_viscosity = solver_container[FLOW_SOL]->GetNodes()->GetLaminarViscosity(jPoint); beta_1 = constants[4]; Solution[1] = 60.0*laminar_viscosity/(density*beta_1*distance*distance); - } + } /*--- Only the \omega solution is modified at the wall. ---*/ Solution[0] = 0.0; @@ -528,7 +528,7 @@ void CTurbSSTSolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_cont Jacobian.DeleteValsRowi(total_index); } } - } + } } void CTurbSSTSolver::BC_Isothermal_Wall(CGeometry *geometry, CSolver **solver_container, CNumerics *conv_numerics, CNumerics *visc_numerics, CConfig *config, diff --git a/SU2_CFD/src/solvers/CTurbSolver.cpp b/SU2_CFD/src/solvers/CTurbSolver.cpp index c3071afc2729..ac23c00388d2 100644 --- a/SU2_CFD/src/solvers/CTurbSolver.cpp +++ b/SU2_CFD/src/solvers/CTurbSolver.cpp @@ -232,7 +232,7 @@ void CTurbSolver::Viscous_Residual(CGeometry *geometry, CSolver **solver_contain /*--- Roughness heights. ---*/ if (config->GetKind_Turb_Model() == SA) numerics->SetRoughness(geometry->node[iPoint]->GetRoughnessHeight(),geometry->node[jPoint]->GetRoughnessHeight()); - + /*--- Compute residual, and Jacobians ---*/ numerics->ComputeResidual(Residual, Jacobian_i, Jacobian_j, config); From f0b88e3f63674a06d5054f9a5b4fb37fc4ba9c5b Mon Sep 17 00:00:00 2001 From: koodlyakshay Date: Thu, 27 Feb 2020 14:56:44 +0100 Subject: [PATCH 11/29] Fix comments mentioned in review and update config tempelate --- Common/include/CConfig.hpp | 4 +-- Common/src/CConfig.cpp | 61 +++++++++++++++++++------------------- config_template.cfg | 15 ++++++++++ 3 files changed, 48 insertions(+), 32 deletions(-) diff --git a/Common/include/CConfig.hpp b/Common/include/CConfig.hpp index 2799f020b7c8..8cce3476fac7 100644 --- a/Common/include/CConfig.hpp +++ b/Common/include/CConfig.hpp @@ -4939,7 +4939,7 @@ class CConfig { * \brief Get the kind of wall. * \return Kind of wall - smooth or rough. */ - unsigned short GetKindWall(string val_marker); + unsigned short GetKindWall(string val_marker) const; /*! * \brief Set the kind of wall - rough or smooth. @@ -6660,7 +6660,7 @@ class CConfig { * \param[in] val_index - Index corresponding to the boundary. * \return The wall roughness height. */ - su2double GetWall_RoughnessHeight(string val_marker); + su2double GetWall_RoughnessHeight(string val_marker) const; /*! * \brief Set the wall roughness height on a wall boundary (Heatflux or Isothermal). diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index c307e721a66a..b7952d603ba2 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -792,6 +792,7 @@ void CConfig::SetPointersNull(void) { Kind_WallFunctions = NULL; IntInfo_WallFunctions = NULL; DoubleInfo_WallFunctions = NULL; + Kind_Wall = NULL; Config_Filenames = NULL; @@ -4666,9 +4667,7 @@ void CConfig::SetPostprocessing(unsigned short val_software, unsigned short val_ /*--- If this option is not used at all, assume all walls are smooth and set roughness height to zero. ---*/ if (nWall_Types == 0) { nWall_Types = nMarker_HeatFlux + nMarker_Isothermal; - if (Roughness_Height != NULL) Roughness_Height = NULL; Roughness_Height = new su2double [nWall_Types]; - if (Kind_Wall != NULL) Kind_Wall = NULL; Kind_Wall = new unsigned short [nWall_Types]; for (unsigned short iMarker = 0; iMarker < nMarker_HeatFlux; iMarker++) { Roughness_Height[iMarker] = 0.0; @@ -4685,13 +4684,11 @@ void CConfig::SetPostprocessing(unsigned short val_software, unsigned short val_ if (nRoughWall != nWall_Types) SU2_MPI::Error("Mismatch in number of roughness height definition and wall type definition.", CURRENT_FUNCTION); else { - su2double *temp_rough; - unsigned short *temp_kindrough; - temp_rough = new su2double [nWall_Types]; - temp_kindrough = new unsigned short [nWall_Types]; + vector temp_rough; + vector temp_kindrough; for (unsigned short iMarker = 0; iMarker < nWall_Types; iMarker++) { - temp_rough[iMarker] = Roughness_Height[iMarker]; - temp_kindrough[iMarker] = Kind_Wall[iMarker]; + temp_rough.push_back(Roughness_Height[iMarker]); + temp_kindrough.push_back(Kind_Wall[iMarker]); } for (iMarker = 0; iMarker < nMarker_HeatFlux; iMarker++) { Roughness_Height[iMarker] = temp_rough[iMarker]; @@ -4701,6 +4698,9 @@ void CConfig::SetPostprocessing(unsigned short val_software, unsigned short val_ Roughness_Height[nMarker_HeatFlux + iMarker] = temp_rough[nMarker_HeatFlux + iMarker]; Kind_Wall[nMarker_HeatFlux + iMarker] = temp_kindrough[nMarker_HeatFlux + iMarker]; } + /*--- Release memory allocated earlier. ---*/ + vector ().swap(temp_rough); + vector ().swap(temp_kindrough); } } } @@ -7566,6 +7566,7 @@ CConfig::~CConfig(void) { if (Kind_Inc_Outlet != NULL) delete[] Kind_Inc_Outlet; if (Kind_WallFunctions != NULL) delete[] Kind_WallFunctions; + if (Kind_Wall != NULL) delete[] Kind_Wall; if (Config_Filenames != NULL) delete[] Config_Filenames; @@ -8903,20 +8904,20 @@ su2double CConfig::GetWall_HeatFlux(string val_marker) { return Heat_Flux[iMarker_HeatFlux]; } -unsigned short CConfig::GetKindWall(string val_marker) { +unsigned short CConfig::GetKindWall(string val_marker) const { unsigned short iMarker; short flag = -1; for (iMarker = 0; iMarker < nMarker_HeatFlux; iMarker++) if (Marker_HeatFlux[iMarker] == val_marker) { - flag = iMarker; - break; - } + flag = iMarker; + break; + } if (flag == -1) - for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) - if (Marker_Isothermal[iMarker] == val_marker) { - flag = iMarker; - break; - } + for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) + if (Marker_Isothermal[iMarker] == val_marker) { + flag = iMarker; + break; + } if (flag == -1) return 1; else return Kind_Wall[flag]; } @@ -8925,23 +8926,23 @@ void CConfig::SetWall_RoughnessHeight(string val_marker, su2double val_roughness //Roughness_Height[val_marker] = val_roughness; } -su2double CConfig::GetWall_RoughnessHeight(string val_marker) { +su2double CConfig::GetWall_RoughnessHeight(string val_marker) const { unsigned short iMarker = 0; short flag = -1; if (nMarker_HeatFlux > 0 || nMarker_Isothermal > 0) { - for (iMarker = 0; iMarker < nMarker_HeatFlux; iMarker++) - if (Marker_HeatFlux[iMarker] == val_marker) { - flag = iMarker; - break; - } - if (flag == -1) { - for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) - if (Marker_Isothermal[iMarker] == val_marker) { - flag = iMarker; - break; - } - } + for (iMarker = 0; iMarker < nMarker_HeatFlux; iMarker++) + if (Marker_HeatFlux[iMarker] == val_marker) { + flag = iMarker; + break; + } + if (flag == -1) { + for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) + if (Marker_Isothermal[iMarker] == val_marker) { + flag = iMarker; + break; + } + } } return Roughness_Height[flag]; } diff --git a/config_template.cfg b/config_template.cfg index 697ad7fbda5d..72d689a7c8e0 100644 --- a/config_template.cfg +++ b/config_template.cfg @@ -770,7 +770,17 @@ GILES_EXTRA_RELAXFACTOR= ( 0.05, 0.05) % % YES Non reflectivity activated, NO the Giles BC behaves as a normal 1D characteristic-based BC SPATIAL_FOURIER= NO +% ------------------------ WALL ROUGHNESS DEFINITION --------------------------% +% Type of wall (SMOOTH or ROUGH). By default all walls are assumed to be smooth. +% If multiple walls are present, list type for each of them (ex- SMOOTH, ROUGH,..) +% in the order they are listed in MARKER_HEATFLUX definition. If all walls are smooth +% this option can be removed. +%WALL_TYPE= SMOOTH +%The equivalent sand grain roughness height (k_s) on each of the wall. This must be in m. +% This is a list of doubles each element corresponding to the MARKER defined in WALL_TYPE. +% For SMOOTH walls, set it to 0.0. +%WALL_ROUGHNESS = 0.0 % ------------------------ SURFACES IDENTIFICATION ----------------------------% % % Marker(s) of the surface in the surface flow solution file @@ -801,6 +811,11 @@ MARKER_ANALYZE_AVERAGE = MASSFLUX % % Numerical method for spatial gradients (GREEN_GAUSS, WEIGHTED_LEAST_SQUARES) NUM_METHOD_GRAD= GREEN_GAUSS + +% Numerical method for spatial gradients to be used for MUSCL reconstruction +% Options are (GREEN_GAUSS, WEIGHTED_LEAST_SQUARES, LEAST_SQUARES). Default value is +% NONE and the method specified in NUM_METHOD_GRAD is used. +NUM_METHOD_GRAD_RECON = LEAST_SQUARES % % CFL number (initial value for the adaptive CFL number) CFL_NUMBER= 15.0 From 807191391ceb61783fc79dd6ab3a8c93a65be35c Mon Sep 17 00:00:00 2001 From: koodlyakshay Date: Fri, 28 Feb 2020 16:35:47 +0100 Subject: [PATCH 12/29] Fix comments, indentations and MARKER_ISOTHERMAL. --- Common/src/CConfig.cpp | 3 - SU2_CFD/src/solvers/CTurbSASolver.cpp | 118 +++++++++++++++++-------- SU2_CFD/src/solvers/CTurbSSTSolver.cpp | 68 +++++++++++--- config_template.cfg | 6 +- 4 files changed, 138 insertions(+), 57 deletions(-) diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index d2b7addad738..e0901b82c5fe 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -4703,9 +4703,6 @@ void CConfig::SetPostprocessing(unsigned short val_software, unsigned short val_ Roughness_Height[nMarker_HeatFlux + iMarker] = temp_rough[nMarker_HeatFlux + iMarker]; Kind_Wall[nMarker_HeatFlux + iMarker] = temp_kindrough[nMarker_HeatFlux + iMarker]; } - /*--- Release memory allocated earlier. ---*/ - vector ().swap(temp_rough); - vector ().swap(temp_kindrough); } } } diff --git a/SU2_CFD/src/solvers/CTurbSASolver.cpp b/SU2_CFD/src/solvers/CTurbSASolver.cpp index 617171b0bc13..27b732d21a7f 100644 --- a/SU2_CFD/src/solvers/CTurbSASolver.cpp +++ b/SU2_CFD/src/solvers/CTurbSASolver.cpp @@ -324,8 +324,6 @@ void CTurbSASolver::Preprocessing(CGeometry *geometry, CSolver **solver_containe void CTurbSASolver::Postprocessing(CGeometry *geometry, CSolver **solver_container, CConfig *config, unsigned short iMesh) { -// su2double rho = 0.0, mu = 0.0, nu, *nu_hat, muT, Ji, Ji_3, fv1; -// su2double cv1_3 = 7.1*7.1*7.1, cR1 = 0.5, roughness, dist; const su2double cv1_3 = 7.1*7.1*7.1, cR1 = 0.5; const bool neg_spalart_allmaras = (config->GetKind_Turb_Model() == SA_NEG); @@ -491,8 +489,7 @@ void CTurbSASolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_conta string Marker_Tag = config->GetMarker_All_TagBound(val_marker); Roughness_Height = config->GetWall_RoughnessHeight(Marker_Tag); - if (Roughness_Height > 0.0 ) rough_wall = true; - unsigned short rank = SU2_MPI::GetRank(); + if (config->GetKindWall(Marker_Tag) == ROUGH ) rough_wall = true; /*--- The dirichlet condition is used only without wall function, otherwise the convergence is compromised as we are providing nu tilde values for the @@ -519,40 +516,40 @@ void CTurbSASolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_conta /*--- Includes 1 in the diagonal ---*/ Jacobian.DeleteValsRowi(iPoint); - } - else { - /*--- For rough walls, the boundary condition is given by - * (\frac{\partial \nu}{\partial n})_wall = \frac{\nu}{0.03*k_s} - * where \nu is the solution variable, $n$ is the wall normal direction - * and k_s is the equivalent sand grain roughness specified. ---*/ + } + else { + /*--- For rough walls, the boundary condition is given by + * (\frac{\partial \nu}{\partial n})_wall = \frac{\nu}{0.03*k_s} + * where \nu is the solution variable, $n$ is the wall normal direction + * and k_s is the equivalent sand grain roughness specified. ---*/ - /*--- Compute dual-grid area and boundary normal ---*/ + /*--- Compute dual-grid area and boundary normal ---*/ - Normal = geometry->vertex[val_marker][iVertex]->GetNormal(); + Normal = geometry->vertex[val_marker][iVertex]->GetNormal(); - Area = 0.0; - for (iDim = 0; iDim < nDim; iDim++) - Area += Normal[iDim]*Normal[iDim]; - Area = sqrt(Area); + Area = 0.0; + for (iDim = 0; iDim < nDim; iDim++) + Area += Normal[iDim]*Normal[iDim]; + Area = sqrt(Area); - /*--- Get laminar_viscosity and density ---*/ + /*--- Get laminar_viscosity and density ---*/ - laminar_viscosity = solver_container[FLOW_SOL]->GetNodes()->GetLaminarViscosity(iPoint); - density = solver_container[FLOW_SOL]->GetNodes()->GetDensity(iPoint); + laminar_viscosity = solver_container[FLOW_SOL]->GetNodes()->GetLaminarViscosity(iPoint); + density = solver_container[FLOW_SOL]->GetNodes()->GetDensity(iPoint); - nu_total = (laminar_viscosity/density + nodes->GetSolution(iPoint,0)); + nu_total = (laminar_viscosity/density + nodes->GetSolution(iPoint,0)); - coeff = (nu_total/sigma); + coeff = (nu_total/sigma); - RoughWallBC = nodes->GetSolution(iPoint,0)/(0.03*Roughness_Height); + RoughWallBC = nodes->GetSolution(iPoint,0)/(0.03*Roughness_Height); - Res_Wall[0] = coeff*RoughWallBC*Area; - LinSysRes.SubtractBlock(iPoint, Res_Wall); + Res_Wall[0] = coeff*RoughWallBC*Area; + LinSysRes.SubtractBlock(iPoint, Res_Wall); - Jacobian_i[0][0] = (laminar_viscosity*Area)/(0.03*Roughness_Height*sigma); - Jacobian_i[0][0] += 2.0*RoughWallBC*Area/sigma; - Jacobian.SubtractBlock(iPoint,iPoint,Jacobian_i); - } + Jacobian_i[0][0] = (laminar_viscosity*Area)/(0.03*Roughness_Height*sigma); + Jacobian_i[0][0] += 2.0*RoughWallBC*Area/sigma; + Jacobian.SubtractBlock(iPoint,iPoint,Jacobian_i); + } } } } @@ -569,7 +566,16 @@ void CTurbSASolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_conta void CTurbSASolver::BC_Isothermal_Wall(CGeometry *geometry, CSolver **solver_container, CNumerics *conv_numerics, CNumerics *visc_numerics, CConfig *config, unsigned short val_marker) { unsigned long iPoint, iVertex; - unsigned short iVar; + unsigned short iVar, iDim; + bool rough_wall = false; + su2double RoughWallBC, Roughness_Height, laminar_viscosity, density, sigma = 2.0/3.0, nu_total, coeff; + su2double *Res_Wall = new su2double [nVar]; + su2double *Normal, Area; + + string Marker_Tag = config->GetMarker_All_TagBound(val_marker); + + Roughness_Height = config->GetWall_RoughnessHeight(Marker_Tag); + if (config->GetKindWall(Marker_Tag) == ROUGH ) rough_wall = true; for (iVertex = 0; iVertex < geometry->nVertex[val_marker]; iVertex++) { iPoint = geometry->vertex[val_marker][iVertex]->GetNode(); @@ -578,18 +584,54 @@ void CTurbSASolver::BC_Isothermal_Wall(CGeometry *geometry, CSolver **solver_con if (geometry->node[iPoint]->GetDomain()) { - /*--- Get the velocity vector ---*/ - for (iVar = 0; iVar < nVar; iVar++) - Solution[iVar] = 0.0; + if (!rough_wall) { + + /*--- Get the solution vector ---*/ + for (iVar = 0; iVar < nVar; iVar++) + Solution[iVar] = 0.0; + + nodes->SetSolution_Old(iPoint,Solution); + LinSysRes.SetBlock_Zero(iPoint); + + /*--- Includes 1 in the diagonal ---*/ + + Jacobian.DeleteValsRowi(iPoint); + } + else { + /*--- For rough walls, the boundary condition is given by + * (\frac{\partial \nu}{\partial n})_wall = \frac{\nu}{0.03*k_s} + * where \nu is the solution variable, $n$ is the wall normal direction + * and k_s is the equivalent sand grain roughness specified. ---*/ - nodes->SetSolution_Old(iPoint,Solution); - LinSysRes.SetBlock_Zero(iPoint); + /*--- Compute dual-grid area and boundary normal ---*/ - /*--- Includes 1 in the diagonal ---*/ + Normal = geometry->vertex[val_marker][iVertex]->GetNormal(); - Jacobian.DeleteValsRowi(iPoint); - } - } + Area = 0.0; + for (iDim = 0; iDim < nDim; iDim++) + Area += Normal[iDim]*Normal[iDim]; + Area = sqrt(Area); + + /*--- Get laminar_viscosity and density ---*/ + + laminar_viscosity = solver_container[FLOW_SOL]->GetNodes()->GetLaminarViscosity(iPoint); + density = solver_container[FLOW_SOL]->GetNodes()->GetDensity(iPoint); + + nu_total = (laminar_viscosity/density + nodes->GetSolution(iPoint,0)); + + coeff = (nu_total/sigma); + + RoughWallBC = nodes->GetSolution(iPoint,0)/(0.03*Roughness_Height); + + Res_Wall[0] = coeff*RoughWallBC*Area; + LinSysRes.SubtractBlock(iPoint, Res_Wall); + + Jacobian_i[0][0] = (laminar_viscosity*Area)/(0.03*Roughness_Height*sigma); + Jacobian_i[0][0] += 2.0*RoughWallBC*Area/sigma; + Jacobian.SubtractBlock(iPoint,iPoint,Jacobian_i); + } + } + } } diff --git a/SU2_CFD/src/solvers/CTurbSSTSolver.cpp b/SU2_CFD/src/solvers/CTurbSSTSolver.cpp index 8ee1691ce0ea..b54986378c31 100644 --- a/SU2_CFD/src/solvers/CTurbSSTSolver.cpp +++ b/SU2_CFD/src/solvers/CTurbSSTSolver.cpp @@ -508,6 +508,11 @@ void CTurbSSTSolver::BC_Isothermal_Wall(CGeometry *geometry, CSolver **solver_co unsigned long iPoint, jPoint, iVertex, total_index; unsigned short iDim, iVar; su2double distance, density = 0.0, laminar_viscosity = 0.0, beta_1; + bool rough_wall = false; + su2double RoughWallBC, Roughness_Height, S_R, FrictionVel, kPlus, WallShearStress; + string Marker_Tag = config->GetMarker_All_TagBound(val_marker); + + if (config->GetKindWall(Marker_Tag) == ROUGH ) rough_wall = true; for (iVertex = 0; iVertex < geometry->nVertex[val_marker]; iVertex++) { iPoint = geometry->vertex[val_marker][iVertex]->GetNode(); @@ -515,24 +520,61 @@ void CTurbSSTSolver::BC_Isothermal_Wall(CGeometry *geometry, CSolver **solver_co /*--- Check if the node belongs to the domain (i.e, not a halo node) ---*/ if (geometry->node[iPoint]->GetDomain()) { - /*--- distance to closest neighbor ---*/ - jPoint = geometry->vertex[val_marker][iVertex]->GetNormal_Neighbor(); - distance = 0.0; - for (iDim = 0; iDim < nDim; iDim++) { - distance += (geometry->node[iPoint]->GetCoord(iDim) - geometry->node[jPoint]->GetCoord(iDim))* - (geometry->node[iPoint]->GetCoord(iDim) - geometry->node[jPoint]->GetCoord(iDim)); - } - distance = sqrt(distance); + if (rough_wall) { - /*--- Set wall values ---*/ + /*--- Set wall values ---*/ - density = solver_container[FLOW_SOL]->GetNodes()->GetDensity(jPoint); - laminar_viscosity = solver_container[FLOW_SOL]->GetNodes()->GetLaminarViscosity(jPoint); + density = solver_container[FLOW_SOL]->GetNodes()->GetDensity(iPoint); + laminar_viscosity = solver_container[FLOW_SOL]->GetNodes()->GetLaminarViscosity(iPoint); - beta_1 = constants[4]; + WallShearStress = solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 0)*solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 0); + WallShearStress += solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 1)*solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 1); + if (nDim == 3) WallShearStress += solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 2)*solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 2); + + WallShearStress = sqrt(WallShearStress); + /*--- Compute non-dimensional velocity ---*/ + FrictionVel = sqrt(fabs(WallShearStress)/density); + /*--- Compute roughness in wall units. ---*/ + Roughness_Height = config->GetWall_RoughnessHeight(Marker_Tag); + kPlus = FrictionVel*Roughness_Height*density/laminar_viscosity; + + /*--- Reference 1 original Wilcox (1998) ---*/ + /*if (kPlus <= 25) + S_R = (50/(kPlus+EPS))*(50/(kPlus+EPS)); + else + S_R = 100/(kPlus+EPS);*/ + + /*--- Reference 2 from D.C. Wilcox Turbulence Modeling for CFD (2006) ---*/ + if (kPlus <= 5) + S_R = (200/(kPlus+EPS))*(200/(kPlus+EPS)); + else + S_R = 100/(kPlus+EPS) + ((200/(kPlus+EPS))*(200/(kPlus+EPS)) - 100/(kPlus+EPS))*exp(5-kPlus); + + /*--- Modify the omega to account for a rough wall. ---*/ + Solution[1] = FrictionVel*FrictionVel*S_R/(laminar_viscosity/density); + } else { + /*--- distance to closest neighbor ---*/ + jPoint = geometry->vertex[val_marker][iVertex]->GetNormal_Neighbor(); + distance = 0.0; + for (iDim = 0; iDim < nDim; iDim++) { + distance += pow(geometry->node[iPoint]->GetCoord(iDim)- + geometry->node[jPoint]->GetCoord(iDim), 2); + } + distance = sqrt(distance); + + /*--- Set wall values ---*/ + + density = solver_container[FLOW_SOL]->GetNodes()->GetDensity(jPoint); + laminar_viscosity = solver_container[FLOW_SOL]->GetNodes()->GetLaminarViscosity(jPoint); + + beta_1 = constants[4]; + + Solution[1] = 60.0*laminar_viscosity/(density*beta_1*distance*distance); + } + + /*--- Only the \omega solution is modified at the wall. ---*/ Solution[0] = 0.0; - Solution[1] = 60.0*laminar_viscosity/(density*beta_1*distance*distance); /*--- Set the solution values and zero the residual ---*/ nodes->SetSolution_Old(iPoint,Solution); diff --git a/config_template.cfg b/config_template.cfg index 72d689a7c8e0..f018e8c7dc64 100644 --- a/config_template.cfg +++ b/config_template.cfg @@ -773,11 +773,11 @@ SPATIAL_FOURIER= NO % ------------------------ WALL ROUGHNESS DEFINITION --------------------------% % Type of wall (SMOOTH or ROUGH). By default all walls are assumed to be smooth. % If multiple walls are present, list type for each of them (ex- SMOOTH, ROUGH,..) -% in the order they are listed in MARKER_HEATFLUX definition. If all walls are smooth -% this option can be removed. +% in the order they are listed in MARKER_HEATFLUX followed by the MARKER_ISOTHERMAL +% definition. If all walls are smooth this option can be removed. %WALL_TYPE= SMOOTH -%The equivalent sand grain roughness height (k_s) on each of the wall. This must be in m. +% The equivalent sand grain roughness height (k_s) on each of the wall. This must be in m. % This is a list of doubles each element corresponding to the MARKER defined in WALL_TYPE. % For SMOOTH walls, set it to 0.0. %WALL_ROUGHNESS = 0.0 From 444698fb2923d32a5e1b19b4ea8d848d0e3fb7bd Mon Sep 17 00:00:00 2001 From: koodlyakshay Date: Mon, 2 Mar 2020 17:09:09 +0100 Subject: [PATCH 13/29] Merge with develop (III). --- Common/src/CConfig.cpp | 3 --- 1 file changed, 3 deletions(-) diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index c1fb4fd16053..759ce9eb98f9 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -7760,11 +7760,8 @@ CConfig::~CConfig(void) { if (Load_Sine_Amplitude != NULL) delete[] Load_Sine_Amplitude; if (Load_Sine_Frequency != NULL) delete[] Load_Sine_Frequency; if (FlowLoad_Value != NULL) delete[] FlowLoad_Value; -<<<<<<< HEAD if (Roughness_Height != NULL) delete[] Roughness_Height; -======= if (Wall_Emissivity != NULL) delete[] Wall_Emissivity; ->>>>>>> develop /*--- related to periodic boundary conditions ---*/ From 782ca72f51384b60ef5cb2c7ec51a4aa25ffcdf1 Mon Sep 17 00:00:00 2001 From: Akshay Date: Thu, 5 Mar 2020 10:44:59 +0100 Subject: [PATCH 14/29] Fix indentation and spaces. --- .../geometry/dual_grid/CTurboVertex.hpp | 2 +- Common/src/geometry/CPhysicalGeometry.cpp | 4 +- SU2_CFD/include/numerics/CNumerics.hpp | 12 ++-- .../include/variables/CTurbSSTVariable.hpp | 2 +- SU2_CFD/include/variables/CTurbVariable.hpp | 6 +- .../src/numerics/turbulent/turb_sources.cpp | 16 ++--- SU2_CFD/src/solvers/CIncNSSolver.cpp | 4 +- SU2_CFD/src/solvers/CNSSolver.cpp | 2 +- SU2_CFD/src/solvers/CTurbSASolver.cpp | 71 +++++++++---------- SU2_CFD/src/solvers/CTurbSSTSolver.cpp | 16 ++--- SU2_CFD/src/solvers/CTurbSolver.cpp | 3 +- SU2_CFD/src/variables/CTurbSSTVariable.cpp | 2 +- SU2_CFD/src/variables/CTurbVariable.cpp | 6 +- 13 files changed, 72 insertions(+), 74 deletions(-) diff --git a/Common/include/geometry/dual_grid/CTurboVertex.hpp b/Common/include/geometry/dual_grid/CTurboVertex.hpp index fa50790debb0..1e43fa2b2f62 100644 --- a/Common/include/geometry/dual_grid/CTurboVertex.hpp +++ b/Common/include/geometry/dual_grid/CTurboVertex.hpp @@ -117,7 +117,7 @@ class CTurboVertex final : public CVertex { * \brief get global index for ordered span-wise turbovertex. */ inline int GetGlobalVertexIndex(void) const {return GlobalIndex;} - + /*! * \brief set angular coord. */ diff --git a/Common/src/geometry/CPhysicalGeometry.cpp b/Common/src/geometry/CPhysicalGeometry.cpp index f5e8443733dc..2b1f153676ec 100644 --- a/Common/src/geometry/CPhysicalGeometry.cpp +++ b/Common/src/geometry/CPhysicalGeometry.cpp @@ -5218,8 +5218,8 @@ void CPhysicalGeometry::ComputeWall_Distance(CConfig *config) { /*--- Use the markerID to find the corresponding wall roughness height. ---*/ - Marker_Tag = config->GetMarker_All_TagBound(markerID); - roughness = config->GetWall_RoughnessHeight(Marker_Tag); + Marker_Tag = config->GetMarker_All_TagBound(markerID); + roughness = config->GetWall_RoughnessHeight(Marker_Tag); node[iPoint]->SetRoughnessHeight(localRoughness); diff --git a/SU2_CFD/include/numerics/CNumerics.hpp b/SU2_CFD/include/numerics/CNumerics.hpp index 792c13c1ccf7..705608461294 100644 --- a/SU2_CFD/include/numerics/CNumerics.hpp +++ b/SU2_CFD/include/numerics/CNumerics.hpp @@ -560,14 +560,14 @@ class CNumerics { /*! * \brief Set the value of the roughness from the nearest wall. - * \param[in] val_dist_i - Value of of the roughness of the nearest wall from point i - * \param[in] val_dist_j - Value of of the roughness of the nearest wall from point j + * \param[in] val_dist_i - Value of of the roughness of the nearest wall from point i + * \param[in] val_dist_j - Value of of the roughness of the nearest wall from point j */ void SetRoughness(su2double val_roughness_i, su2double val_roughness_j) { - roughness_i = val_roughness_i; - roughness_j = val_roughness_j; -} - + roughness_i = val_roughness_i; + roughness_j = val_roughness_j; + } + /*! * \brief Set coordinates of the points. * \param[in] val_coord_i - Coordinates of the point i. diff --git a/SU2_CFD/include/variables/CTurbSSTVariable.hpp b/SU2_CFD/include/variables/CTurbSSTVariable.hpp index f762664658d1..899bfb50cea4 100644 --- a/SU2_CFD/include/variables/CTurbSSTVariable.hpp +++ b/SU2_CFD/include/variables/CTurbSSTVariable.hpp @@ -6,7 +6,7 @@ * * SU2 Project Website: https://su2code.github.io * - * The SU2 Project is maintained by the SU2 Foundation + * The SU2 Project is maintained by the SU2 Foundation * (http://su2foundation.org) * * Copyright 2012-2020, SU2 Contributors (cf. AUTHORS.md) diff --git a/SU2_CFD/include/variables/CTurbVariable.hpp b/SU2_CFD/include/variables/CTurbVariable.hpp index 2d07630016ca..fcf462b5a09f 100644 --- a/SU2_CFD/include/variables/CTurbVariable.hpp +++ b/SU2_CFD/include/variables/CTurbVariable.hpp @@ -6,7 +6,7 @@ * * SU2 Project Website: https://su2code.github.io * - * The SU2 Project is maintained by the SU2 Foundation + * The SU2 Project is maintained by the SU2 Foundation * (http://su2foundation.org) * * Copyright 2012-2020, SU2 Contributors (cf. AUTHORS.md) @@ -82,7 +82,7 @@ class CTurbVariable : public CVariable { inline su2double GetGradient_Reconstruction(unsigned long iPoint, unsigned long iVar, unsigned long iDim) const final { return Gradient_Reconstruction(iPoint,iVar,iDim); } - + /*! * \brief Get the value of the reconstruction variables gradient at a node. * \param[in] iPoint - Index of the current node. @@ -93,7 +93,7 @@ class CTurbVariable : public CVariable { inline void SetGradient_Reconstruction(unsigned long iPoint, unsigned long iVar, unsigned long iDim, su2double value) final { Gradient_Reconstruction(iPoint,iVar,iDim) = value; } - + /*! * \brief Get the array of the reconstruction variables gradient at a node. * \param[in] iPoint - Index of the current node. diff --git a/SU2_CFD/src/numerics/turbulent/turb_sources.cpp b/SU2_CFD/src/numerics/turbulent/turb_sources.cpp index fc1abbaf350b..526473ab806e 100644 --- a/SU2_CFD/src/numerics/turbulent/turb_sources.cpp +++ b/SU2_CFD/src/numerics/turbulent/turb_sources.cpp @@ -79,7 +79,7 @@ CNumerics::ResidualType<> CSourcePieceWise_TurbSA::ComputeResidual(const CConfig // BC Transition Model variables su2double vmag, rey, re_theta, re_theta_t, re_v; su2double tu , nu_cr, nu_t, nu_BC, chi_1, chi_2, term1, term2, term_exponential; - + // Set the boolean here depending on whether the point is closest to a rough wall or not. roughwall = (roughness_i > 0.0); @@ -124,31 +124,31 @@ CNumerics::ResidualType<> CSourcePieceWise_TurbSA::ComputeResidual(const CConfig dist_i_2 = dist_i*dist_i; nu = Laminar_Viscosity_i/Density_i; - + /*--- Old values without roughness ---*/ /*Ji = TurbVar_i[0]/nu; Ji_2 = Ji*Ji; Ji_3 = Ji_2*Ji; fv1 = Ji_3/(Ji_3+cv1_3); fv2 = 1.0 - Ji/(1.0+Ji*fv1);*/ - + /*--- Modified values for roughness ---*/ - /*--- Ref: Aupoix, B. and Spalart, P. R., "Extensions of the Spalart-Allmaras Turbulence Model to Account for Wall Roughness," + /*--- Ref: Aupoix, B. and Spalart, P. R., "Extensions of the Spalart-Allmaras Turbulence Model to Account for Wall Roughness," * International Journal of Heat and Fluid Flow, Vol. 24, 2003, pp. 454-462. ---*/ /* --- See https://turbmodels.larc.nasa.gov/spalart.html#sarough for detailed explanation. ---*/ - + Ji = TurbVar_i[0]/nu + cr1*(roughness_i/(dist_i+EPS)); //roughness_i = 0 for smooth walls and Ji remains the same, changes only if roughness is specified. Ji_2 = Ji*Ji; Ji_3 = Ji_2*Ji; fv1 = Ji_3/(Ji_3+cv1_3); - + /*--- Using a modified relation so as to not change the Shat that depends on fv2. ---*/ fv2 = 1.0 - TurbVar_i[0]/(nu+TurbVar_i[0]*fv1); // From NASA turb modeling resource and 2003 paper - + ft2 = ct3*exp(-ct4*Ji_2); S = Omega; inv_k2_d2 = 1.0/(k2*dist_i_2); - + Shat = S + TurbVar_i[0]*fv2*inv_k2_d2; Shat = max(Shat, 1.0e-10); inv_Shat = 1.0/Shat; diff --git a/SU2_CFD/src/solvers/CIncNSSolver.cpp b/SU2_CFD/src/solvers/CIncNSSolver.cpp index 24f9d7cca356..00efce358857 100644 --- a/SU2_CFD/src/solvers/CIncNSSolver.cpp +++ b/SU2_CFD/src/solvers/CIncNSSolver.cpp @@ -1261,13 +1261,13 @@ void CIncNSSolver::Friction_Forces(CGeometry *geometry, CConfig *config) { Grad_Temp[iDim] = nodes->GetGradient_Primitive(iPoint,nDim+1, iDim); } - /*--- TurbViscosity is usually zero on the walls, except when the SA roughness model is used. + /*--- TurbViscosity is usually zero on the walls, except when the SA roughness model is used. * If the SA roughness model is not active, this term will be zero and will not affect the solution---*/ TurbViscosity = 0.0; if ((config->GetKindWall(Marker_Tag) == ROUGH ) && (config->GetKind_Turb_Model() == SA)) TurbViscosity = nodes->GetEddyViscosity(iPoint); Viscosity = nodes->GetLaminarViscosity(iPoint); Viscosity += TurbViscosity; - + Density = nodes->GetDensity(iPoint); Area = 0.0; for (iDim = 0; iDim < nDim; iDim++) Area += Normal[iDim]*Normal[iDim]; Area = sqrt(Area); diff --git a/SU2_CFD/src/solvers/CNSSolver.cpp b/SU2_CFD/src/solvers/CNSSolver.cpp index 09b253ebe354..e92199c5211d 100644 --- a/SU2_CFD/src/solvers/CNSSolver.cpp +++ b/SU2_CFD/src/solvers/CNSSolver.cpp @@ -562,7 +562,7 @@ void CNSSolver::Friction_Forces(CGeometry *geometry, CConfig *config) { Grad_Temp[iDim] = nodes->GetGradient_Primitive(iPoint,0, iDim); } - /*--- TurbViscosity is usually zero on the walls, except when the SA roughness model is used. + /*--- TurbViscosity is usually zero on the walls, except when the SA roughness model is used. * If the SA roughness model is not active, this term will be zero and will not affect the solution---*/ TurbViscosity = 0.0; if ((config->GetKindWall(Marker_Tag) == ROUGH ) && (config->GetKind_Turb_Model() == SA)) TurbViscosity = nodes->GetEddyViscosity(iPoint); diff --git a/SU2_CFD/src/solvers/CTurbSASolver.cpp b/SU2_CFD/src/solvers/CTurbSASolver.cpp index 0d2b7db29e50..01cb98c2d47c 100644 --- a/SU2_CFD/src/solvers/CTurbSASolver.cpp +++ b/SU2_CFD/src/solvers/CTurbSASolver.cpp @@ -340,13 +340,13 @@ void CTurbSASolver::Postprocessing(CGeometry *geometry, CSolver **solver_contain su2double nu_hat = nodes->GetSolution(iPoint,0); su2double roughness = geometry->node[iPoint]->GetRoughnessHeight(); su2double dist = geometry->node[iPoint]->GetWall_Distance(); - + dist += 0.03*roughness; su2double Ji = nu_hat/nu ; - if (roughness > 1.0e-10) + if (roughness > 1.0e-10) Ji+= cR1*roughness/(dist+EPS); - + su2double Ji_3 = Ji*Ji*Ji; su2double fv1 = Ji_3/(Ji_3+cv1_3); @@ -412,8 +412,8 @@ void CTurbSASolver::Source_Residual(CGeometry *geometry, CSolver **solver_contai if (config->GetKind_HybridRANSLES() == NO_HYBRIDRANSLES) { - /*--- For the SA model, wall roughness is accounted by modifying the computed wall distance - * d_new = d + 0.03 k_s + /*--- For the SA model, wall roughness is accounted by modifying the computed wall distance + * d_new = d + 0.03 k_s * where k_s is the equivalent sand grain roughness height that is specified in cfg file. * For smooth walls, wall roughness is zero and computed wall distance remains the same. */ @@ -507,8 +507,7 @@ void CTurbSASolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_conta if (!rough_wall) { /*--- Get the solution vector ---*/ - for (iVar = 0; iVar < nVar; iVar++) - Solution[iVar] = 0.0; + Solution[0] = 0.0; nodes->SetSolution_Old(iPoint,Solution); LinSysRes.SetBlock_Zero(iPoint); @@ -548,7 +547,7 @@ void CTurbSASolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_conta Jacobian_i[0][0] = (laminar_viscosity*Area)/(0.03*Roughness_Height*sigma); Jacobian_i[0][0] += 2.0*RoughWallBC*Area/sigma; - Jacobian.SubtractBlock(iPoint,iPoint,Jacobian_i); + Jacobian.SubtractBlock2Diag(iPoint,Jacobian_i); } } } @@ -587,48 +586,46 @@ void CTurbSASolver::BC_Isothermal_Wall(CGeometry *geometry, CSolver **solver_con if (!rough_wall) { /*--- Get the solution vector ---*/ - for (iVar = 0; iVar < nVar; iVar++) - Solution[iVar] = 0.0; + Solution[0] = 0.0; - nodes->SetSolution_Old(iPoint,Solution); - LinSysRes.SetBlock_Zero(iPoint); + nodes->SetSolution_Old(iPoint,Solution); + LinSysRes.SetBlock_Zero(iPoint); - /*--- Includes 1 in the diagonal ---*/ + /*--- Includes 1 in the diagonal ---*/ + Jacobian.DeleteValsRowi(iPoint); - Jacobian.DeleteValsRowi(iPoint); - } - else { - /*--- For rough walls, the boundary condition is given by - * (\frac{\partial \nu}{\partial n})_wall = \frac{\nu}{0.03*k_s} - * where \nu is the solution variable, $n$ is the wall normal direction - * and k_s is the equivalent sand grain roughness specified. ---*/ + } else { + /*--- For rough walls, the boundary condition is given by + * (\frac{\partial \nu}{\partial n})_wall = \frac{\nu}{0.03*k_s} + * where \nu is the solution variable, $n$ is the wall normal direction + * and k_s is the equivalent sand grain roughness specified. ---*/ - /*--- Compute dual-grid area and boundary normal ---*/ + /*--- Compute dual-grid area and boundary normal ---*/ - Normal = geometry->vertex[val_marker][iVertex]->GetNormal(); + Normal = geometry->vertex[val_marker][iVertex]->GetNormal(); - Area = 0.0; - for (iDim = 0; iDim < nDim; iDim++) - Area += Normal[iDim]*Normal[iDim]; - Area = sqrt(Area); + Area = 0.0; + for (iDim = 0; iDim < nDim; iDim++) + Area += Normal[iDim]*Normal[iDim]; + Area = sqrt(Area); - /*--- Get laminar_viscosity and density ---*/ + /*--- Get laminar_viscosity and density ---*/ - laminar_viscosity = solver_container[FLOW_SOL]->GetNodes()->GetLaminarViscosity(iPoint); - density = solver_container[FLOW_SOL]->GetNodes()->GetDensity(iPoint); + laminar_viscosity = solver_container[FLOW_SOL]->GetNodes()->GetLaminarViscosity(iPoint); + density = solver_container[FLOW_SOL]->GetNodes()->GetDensity(iPoint); - nu_total = (laminar_viscosity/density + nodes->GetSolution(iPoint,0)); + nu_total = (laminar_viscosity/density + nodes->GetSolution(iPoint,0)); - coeff = (nu_total/sigma); + coeff = (nu_total/sigma); - RoughWallBC = nodes->GetSolution(iPoint,0)/(0.03*Roughness_Height); + RoughWallBC = nodes->GetSolution(iPoint,0)/(0.03*Roughness_Height); - Res_Wall[0] = coeff*RoughWallBC*Area; - LinSysRes.SubtractBlock(iPoint, Res_Wall); + Res_Wall[0] = coeff*RoughWallBC*Area; + LinSysRes.SubtractBlock(iPoint, Res_Wall); - Jacobian_i[0][0] = (laminar_viscosity*Area)/(0.03*Roughness_Height*sigma); - Jacobian_i[0][0] += 2.0*RoughWallBC*Area/sigma; - Jacobian.SubtractBlock(iPoint,iPoint,Jacobian_i); + Jacobian_i[0][0] = (laminar_viscosity*Area)/(0.03*Roughness_Height*sigma); + Jacobian_i[0][0] += 2.0*RoughWallBC*Area/sigma; + Jacobian.SubtractBlock2Diag(iPoint,Jacobian_i); } } } diff --git a/SU2_CFD/src/solvers/CTurbSSTSolver.cpp b/SU2_CFD/src/solvers/CTurbSSTSolver.cpp index bdedd5c7b9e8..135bbe763753 100644 --- a/SU2_CFD/src/solvers/CTurbSSTSolver.cpp +++ b/SU2_CFD/src/solvers/CTurbSSTSolver.cpp @@ -433,15 +433,15 @@ void CTurbSSTSolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_cont if (geometry->node[iPoint]->GetDomain()) { if (rough_wall) { - + /*--- Set wall values ---*/ density = solver_container[FLOW_SOL]->GetNodes()->GetDensity(iPoint); laminar_viscosity = solver_container[FLOW_SOL]->GetNodes()->GetLaminarViscosity(iPoint); - WallShearStress = solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 0)*solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 0); - WallShearStress += solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 1)*solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 1); - if (nDim == 3) WallShearStress += solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 2)*solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 2); + WallShearStress = 0.0; + for (iDim = 0; iDim < nDim; iDim++) + WallShearStress += pow(solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, iDim),2.0); WallShearStress = sqrt(WallShearStress); /*--- Compute non-dimensional velocity ---*/ @@ -511,7 +511,7 @@ void CTurbSSTSolver::BC_Isothermal_Wall(CGeometry *geometry, CSolver **solver_co bool rough_wall = false; su2double RoughWallBC, Roughness_Height, S_R, FrictionVel, kPlus, WallShearStress; string Marker_Tag = config->GetMarker_All_TagBound(val_marker); - + if (config->GetKindWall(Marker_Tag) == ROUGH ) rough_wall = true; for (iVertex = 0; iVertex < geometry->nVertex[val_marker]; iVertex++) { @@ -527,9 +527,9 @@ void CTurbSSTSolver::BC_Isothermal_Wall(CGeometry *geometry, CSolver **solver_co density = solver_container[FLOW_SOL]->GetNodes()->GetDensity(iPoint); laminar_viscosity = solver_container[FLOW_SOL]->GetNodes()->GetLaminarViscosity(iPoint); - WallShearStress = solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 0)*solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 0); - WallShearStress += solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 1)*solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 1); - if (nDim == 3) WallShearStress += solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 2)*solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, 2); + WallShearStress = 0.0; + for (iDim = 0; iDim < nDim; iDim++) + WallShearStress += pow(solver_container[FLOW_SOL]->GetCSkinFriction(val_marker, iVertex, iDim),2.0); WallShearStress = sqrt(WallShearStress); /*--- Compute non-dimensional velocity ---*/ diff --git a/SU2_CFD/src/solvers/CTurbSolver.cpp b/SU2_CFD/src/solvers/CTurbSolver.cpp index 1763adc70722..a12d7c7c99ef 100644 --- a/SU2_CFD/src/solvers/CTurbSolver.cpp +++ b/SU2_CFD/src/solvers/CTurbSolver.cpp @@ -271,7 +271,8 @@ void CTurbSolver::Viscous_Residual(CGeometry *geometry, CSolver **solver_contain /*--- Roughness heights. ---*/ if (config->GetKind_Turb_Model() == SA) - numerics->SetRoughness(geometry->node[iPoint]->GetRoughnessHeight(),geometry->node[jPoint]->GetRoughnessHeight()); + numerics->SetRoughness(geometry->node[iPoint]->GetRoughnessHeight(), + geometry->node[jPoint]->GetRoughnessHeight()); /*--- Compute residual, and Jacobians ---*/ diff --git a/SU2_CFD/src/variables/CTurbSSTVariable.cpp b/SU2_CFD/src/variables/CTurbSSTVariable.cpp index aa961b321f2d..b2ea04777bfc 100644 --- a/SU2_CFD/src/variables/CTurbSSTVariable.cpp +++ b/SU2_CFD/src/variables/CTurbSSTVariable.cpp @@ -6,7 +6,7 @@ * * SU2 Project Website: https://su2code.github.io * - * The SU2 Project is maintained by the SU2 Foundation + * The SU2 Project is maintained by the SU2 Foundation * (http://su2foundation.org) * * Copyright 2012-2020, SU2 Contributors (cf. AUTHORS.md) diff --git a/SU2_CFD/src/variables/CTurbVariable.cpp b/SU2_CFD/src/variables/CTurbVariable.cpp index db6a71c2a1bf..6b1bb162615c 100644 --- a/SU2_CFD/src/variables/CTurbVariable.cpp +++ b/SU2_CFD/src/variables/CTurbVariable.cpp @@ -6,7 +6,7 @@ * * SU2 Project Website: https://su2code.github.io * - * The SU2 Project is maintained by the SU2 Foundation + * The SU2 Project is maintained by the SU2 Foundation * (http://su2foundation.org) * * Copyright 2012-2020, SU2 Contributors (cf. AUTHORS.md) @@ -45,7 +45,7 @@ CTurbVariable::CTurbVariable(unsigned long npoint, unsigned long ndim, unsigned if (config->GetReconstructionGradientRequired()) { Gradient_Aux.resize(nPoint,nVar,nDim,0.0); } - + if (config->GetLeastSquaresRequired()) { Rmatrix.resize(nPoint,nDim,nDim,0.0); } @@ -62,5 +62,5 @@ CTurbVariable::CTurbVariable(unsigned long npoint, unsigned long ndim, unsigned /* Under-relaxation parameter. */ UnderRelaxation.resize(nPoint) = su2double(1.0); LocalCFL.resize(nPoint) = su2double(0.0); - + } From 1e22601704b98451a45afec95dd8309bd696ff26 Mon Sep 17 00:00:00 2001 From: koodlyakshay Date: Thu, 30 Apr 2020 10:41:41 +0200 Subject: [PATCH 15/29] Fix issue with cht test case. --- .travis.yml | 2 +- Common/src/CConfig.cpp | 9 ++++++--- 2 files changed, 7 insertions(+), 4 deletions(-) diff --git a/.travis.yml b/.travis.yml index 99dabda9cc4c..10855686cc03 100644 --- a/.travis.yml +++ b/.travis.yml @@ -18,7 +18,7 @@ compiler: notifications: email: recipients: - - su2code-dev@lists.stanford.edu + - koodlyakshay@gmail.com branches: only: diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index 26cec6c31942..98d757a5b3a3 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -5310,7 +5310,7 @@ void CConfig::SetMarkers(unsigned short val_software) { Marker_CfgFile_KindBC[iMarker_CfgFile] = FLOWLOAD_BOUNDARY; iMarker_CfgFile++; } - + for (iMarker_CfgFile = 0; iMarker_CfgFile < nMarker_CfgFile; iMarker_CfgFile++) { Marker_CfgFile_Monitoring[iMarker_CfgFile] = NO; for (iMarker_Monitoring = 0; iMarker_Monitoring < nMarker_Monitoring; iMarker_Monitoring++) @@ -5438,7 +5438,6 @@ void CConfig::SetMarkers(unsigned short val_software) { if (Marker_CfgFile_TagBound[iMarker_CfgFile] == Marker_PyCustom[iMarker_PyCustom]) Marker_CfgFile_PyCustom[iMarker_CfgFile] = YES; } - } void CConfig::SetOutput(unsigned short val_software, unsigned short val_izone) { @@ -8784,7 +8783,11 @@ su2double CConfig::GetWall_RoughnessHeight(string val_marker) const { } } } - return Roughness_Height[flag]; + + if (flag == -1) + return 0.0; //For cht_wall_interface + else + return Roughness_Height[flag]; } unsigned short CConfig::GetWallFunction_Treatment(string val_marker) const { From b56ff109513a6b15021e262c9be3fa4c06945758 Mon Sep 17 00:00:00 2001 From: koodlyakshay Date: Tue, 19 May 2020 12:12:03 +0200 Subject: [PATCH 16/29] Remove roughness elements in ADT files and add routines to get the correct roughness on every rank. --- Common/include/CConfig.hpp | 16 +++- Common/include/adt_structure.hpp | 11 +-- Common/include/geometry/CGeometry.hpp | 2 + Common/include/geometry/CPhysicalGeometry.hpp | 2 + Common/src/CConfig.cpp | 40 +++++++-- Common/src/adt_structure.cpp | 11 +-- Common/src/fem_geometry_structure.cpp | 5 +- Common/src/fem_wall_distance.cpp | 16 ++-- Common/src/geometry/CPhysicalGeometry.cpp | 87 ++++++++++++++++--- Common/src/geometry_structure_fem_part.cpp | 5 +- SU2_CFD/src/drivers/CDriver.cpp | 4 + 11 files changed, 144 insertions(+), 55 deletions(-) diff --git a/Common/include/CConfig.hpp b/Common/include/CConfig.hpp index 265c887027a4..1f46a1b394bd 100644 --- a/Common/include/CConfig.hpp +++ b/Common/include/CConfig.hpp @@ -308,7 +308,8 @@ class CConfig { su2double *Outlet_Pressure; /*!< \brief Specified back pressures (static) for outlet boundaries. */ su2double *Isothermal_Temperature; /*!< \brief Specified isothermal wall temperatures (static). */ su2double *Heat_Flux; /*!< \brief Specified wall heat fluxes. */ - su2double *Roughness_Height; /*!< \brief Equivalent sand grain roughness for the marker. */ + su2double *Roughness_Height; /*!< \brief Equivalent sand grain roughness for the marker according to config file. */ + su2double *GlobalRoughness_Height; /*!< \brief Equivalent sand grain roughness for markers globally including S-R. */ su2double *Displ_Value; /*!< \brief Specified displacement for displacement boundaries. */ su2double *Load_Value; /*!< \brief Specified force for load boundaries. */ su2double *Damper_Constant; /*!< \brief Specified constant for damper boundaries. */ @@ -389,6 +390,7 @@ class CConfig { *Marker_CfgFile_KindBC; /*!< \brief Global index for boundaries using config file. */ short *Marker_All_SendRecv; /*!< \brief Information about if the boundary is sended (+), received (-). */ short *Marker_All_PerBound; /*!< \brief Global index for periodic bc using the grid information. */ + int *GlobalMarkerStorageDispl; /*!< \brief storage of markers globally. */ unsigned long nExtIter; /*!< \brief Number of external iterations. */ unsigned long ExtIter; /*!< \brief Current external iteration number. */ @@ -2931,6 +2933,12 @@ class CConfig { */ unsigned short GetnMarker_HeatFlux(void) const { return nMarker_HeatFlux; } + /*! + * \brief Get the total number of rough markers. + * \return Total number of heat flux markers. + */ + unsigned short GetnRoughWall(void) const { return nRoughWall; } + /*! * \brief Get the total number of objectives in kind_objective list * \return Total number of objectives in kind_objective list @@ -5092,6 +5100,12 @@ class CConfig { */ su2double GetCoeff_ObjChainRule(unsigned short iVar) const { return Obj_ChainRuleCoeff[iVar]; } + /*--- New roughness mpi comm stuff akshay ----*/ + su2double GetLocalRoughness(int rankID, unsigned short markerID) const; + + void SetGlobalMarkerArray(int *global_displ, int size) ; + void SetGlobalRoughnessArray(su2double *global_rough, int sizeGlobal) ; + /*! * \author H. Kline * \brief Get the flag indicating whether to comput a combined objective. diff --git a/Common/include/adt_structure.hpp b/Common/include/adt_structure.hpp index 7bc389b65490..8fdbbeff508e 100644 --- a/Common/include/adt_structure.hpp +++ b/Common/include/adt_structure.hpp @@ -289,8 +289,6 @@ class CADTElemClass : public CADTBaseClass { of the elements in the ADT. */ vector ranksOfElems; /*!< \brief Vector, which contains the ranks of the elements in the ADT. */ - vector markerRoughness; - #ifdef HAVE_OMP vector >BBoxTargets; /*!< \brief Vector, used to store possible bounding box candidates during the nearest element search. */ @@ -318,7 +316,6 @@ class CADTElemClass : public CADTBaseClass { vector &val_VTKElem, vector &val_markerID, vector &val_elemID, - vector &val_roughheight, const bool globalTree); /*! @@ -360,11 +357,10 @@ class CADTElemClass : public CADTBaseClass { su2double &dist, unsigned short &markerID, unsigned long &elemID, - int &rankID, - su2double &localRoughness) { + int &rankID) { const auto iThread = omp_get_thread_num(); DetermineNearestElement_impl(BBoxTargets[iThread], FrontLeaves[iThread], - FrontLeavesNew[iThread], coor, dist, markerID, elemID, rankID,localRoughness); + FrontLeavesNew[iThread], coor, dist, markerID, elemID, rankID); } private: @@ -392,8 +388,7 @@ class CADTElemClass : public CADTBaseClass { su2double &dist, unsigned short &markerID, unsigned long &elemID, - int &rankID, - su2double &localRoughness) const; + int &rankID) const; /*! * \brief Function, which checks whether or not the given coordinate is diff --git a/Common/include/geometry/CGeometry.hpp b/Common/include/geometry/CGeometry.hpp index 3e3e4392a032..fa16a2a3d1b2 100644 --- a/Common/include/geometry/CGeometry.hpp +++ b/Common/include/geometry/CGeometry.hpp @@ -1692,6 +1692,8 @@ class CGeometry { * \param[in] geometry_container - Geometrical definition of the problem. */ static void ComputeWallDistance(const CConfig * const *config_container, CGeometry ****geometry_container); + + inline virtual void SetGlobalMarkerRoughness(CConfig *config) { } }; diff --git a/Common/include/geometry/CPhysicalGeometry.hpp b/Common/include/geometry/CPhysicalGeometry.hpp index 9bde45df4fd5..ec013389fdf2 100644 --- a/Common/include/geometry/CPhysicalGeometry.hpp +++ b/Common/include/geometry/CPhysicalGeometry.hpp @@ -803,5 +803,7 @@ class CPhysicalGeometry final : public CGeometry { nodes->SetWall_Distance(iPoint, val); } } + + void SetGlobalMarkerRoughness(CConfig *config) override; }; diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index 9c431e154a10..53e43ce010de 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -832,8 +832,9 @@ void CConfig::SetPointersNull(void) { /*--- Boundary Condition settings ---*/ Isothermal_Temperature = nullptr; - Heat_Flux = nullptr; Displ_Value = nullptr; Load_Value = nullptr; - FlowLoad_Value = nullptr; Damper_Constant = nullptr; Wall_Emissivity = nullptr; + Heat_Flux = nullptr; Displ_Value = nullptr; Load_Value = nullptr; + FlowLoad_Value = nullptr; Damper_Constant = nullptr; Wall_Emissivity = nullptr; + Roughness_Height = nullptr; GlobalMarkerStorageDispl = nullptr; GlobalRoughness_Height = nullptr; /*--- Inlet Outlet Boundary Condition settings ---*/ @@ -864,9 +865,7 @@ void CConfig::SetPointersNull(void) { Inlet_FlowDir = nullptr; Inlet_Temperature = nullptr; Inlet_Pressure = nullptr; Inlet_Velocity = nullptr; Inflow_Mach = nullptr; Inflow_Pressure = nullptr; Exhaust_Pressure = nullptr; Outlet_Pressure = nullptr; Isothermal_Temperature= nullptr; - Heat_Flux = nullptr; Displ_Value = nullptr; Load_Value = nullptr; - FlowLoad_Value = nullptr; Roughness_Height = nullptr; - + ElasticityMod = nullptr; PoissonRatio = nullptr; MaterialDensity = nullptr; Load_Dir = nullptr; Load_Dir_Value = nullptr; Load_Dir_Multiplier = nullptr; @@ -5448,6 +5447,7 @@ void CConfig::SetMarkers(unsigned short val_software) { if (Marker_CfgFile_TagBound[iMarker_CfgFile] == Marker_PyCustom[iMarker_PyCustom]) Marker_CfgFile_PyCustom[iMarker_CfgFile] = YES; } + } void CConfig::SetOutput(unsigned short val_software, unsigned short val_izone) { @@ -7624,6 +7624,8 @@ CConfig::~CConfig(void) { delete[] FlowLoad_Value; delete[] Roughness_Height; delete[] Wall_Emissivity; + delete[] GlobalMarkerStorageDispl; + delete[] GlobalRoughness_Height; /*--- related to periodic boundary conditions ---*/ for (iMarker = 0; iMarker < nMarker_PerBound; iMarker++) { @@ -9767,3 +9769,31 @@ void CConfig::SetMultizone(CConfig *driver_config, CConfig **config_container){ Multizone_Residual = true; } + +su2double CConfig::GetLocalRoughness(int rankID, unsigned short markerID) const { + + su2double localRoughness = 0.0; + int index; + + index = GlobalMarkerStorageDispl[rankID]+ markerID; + localRoughness = GlobalRoughness_Height[index]; + +return localRoughness; +} + +void CConfig::SetGlobalMarkerArray(int* global_displ, int size) { + + if (GlobalMarkerStorageDispl == nullptr) GlobalMarkerStorageDispl = new int [size]; + GlobalMarkerStorageDispl[0] = 0; + for (int iRank = 0; iRank < size; iRank++) + GlobalMarkerStorageDispl[iRank] = global_displ[iRank]; + +} + +void CConfig::SetGlobalRoughnessArray(su2double *global_rough, int sizeGlobal) { + + if (GlobalRoughness_Height == nullptr) GlobalRoughness_Height = new su2double [sizeGlobal]; + for (int iMarker = 0; iMarker < sizeGlobal; iMarker++) + GlobalRoughness_Height[iMarker] = global_rough[iMarker]; + +} diff --git a/Common/src/adt_structure.cpp b/Common/src/adt_structure.cpp index 477bbea069b5..192199ef31f6 100644 --- a/Common/src/adt_structure.cpp +++ b/Common/src/adt_structure.cpp @@ -471,7 +471,6 @@ CADTElemClass::CADTElemClass(unsigned short val_nDim, vector &val_VTKElem, vector &val_markerID, vector &val_elemID, - vector &val_roughheight, const bool globalTree) { /* Copy the dimension of the problem into nDim. */ @@ -547,7 +546,6 @@ CADTElemClass::CADTElemClass(unsigned short val_nDim, elemVTK_Type.resize(sizeGlobal); localMarkers.resize(sizeGlobal); - markerRoughness.resize(sizeGlobal); localElemIDs.resize(sizeGlobal); SU2_MPI::Allgatherv(val_VTKElem.data(), sizeLocal, MPI_UNSIGNED_SHORT, elemVTK_Type.data(), @@ -555,9 +553,6 @@ CADTElemClass::CADTElemClass(unsigned short val_nDim, SU2_MPI::Allgatherv(val_markerID.data(), sizeLocal, MPI_UNSIGNED_SHORT, localMarkers.data(), recvCounts.data(), displs.data(), MPI_UNSIGNED_SHORT, MPI_COMM_WORLD); - - SU2_MPI::Allgatherv(val_roughheight.data(), sizeLocal, MPI_DOUBLE, markerRoughness.data(), - recvCounts.data(), displs.data(), MPI_DOUBLE, MPI_COMM_WORLD); SU2_MPI::Allgatherv(val_elemID.data(), sizeLocal, MPI_UNSIGNED_LONG, localElemIDs.data(), recvCounts.data(), displs.data(), MPI_UNSIGNED_LONG, MPI_COMM_WORLD); @@ -600,7 +595,6 @@ CADTElemClass::CADTElemClass(unsigned short val_nDim, elemVTK_Type = val_VTKElem; localMarkers = val_markerID; localElemIDs = val_elemID; - markerRoughness = val_roughheight; ranksOfElems.assign(elemVTK_Type.size(), rank); } @@ -613,7 +607,6 @@ CADTElemClass::CADTElemClass(unsigned short val_nDim, elemVTK_Type = val_VTKElem; localMarkers = val_markerID; localElemIDs = val_elemID; - markerRoughness = val_roughheight; ranksOfElems.assign(elemVTK_Type.size(), MASTER_NODE); @@ -795,8 +788,7 @@ void CADTElemClass::DetermineNearestElement_impl(vector& BBoxT su2double &dist, unsigned short &markerID, unsigned long &elemID, - int &rankID, - su2double &localRoughness) const { + int &rankID) const { AD_BEGIN_PASSIVE @@ -958,7 +950,6 @@ void CADTElemClass::DetermineNearestElement_impl(vector& BBoxT markerID = localMarkers[ii]; elemID = localElemIDs[ii]; rankID = ranksOfElems[ii]; - localRoughness = markerRoughness[ii]; } } diff --git a/Common/src/fem_geometry_structure.cpp b/Common/src/fem_geometry_structure.cpp index 327d8a28fec3..518759f8c3bc 100644 --- a/Common/src/fem_geometry_structure.cpp +++ b/Common/src/fem_geometry_structure.cpp @@ -6220,7 +6220,6 @@ void CMeshFEM_DG::WallFunctionPreprocessing(CConfig *config) { vector subElementIDInParent; vector VTK_TypeElem; vector elemConn; - vector dummyRough; /* Loop over the locally stored volume elements (including halo elements) to create the connectivity of the subelements. */ @@ -6257,7 +6256,6 @@ void CMeshFEM_DG::WallFunctionPreprocessing(CConfig *config) { for(unsigned short j=0; j().swap(parentElement); vector().swap(elemConn); vector().swap(volCoor); - vector().swap(dummyRough); /*--------------------------------------------------------------------------*/ /*--- Step 3. Search for donor elements at the exchange locations in ---*/ diff --git a/Common/src/fem_wall_distance.cpp b/Common/src/fem_wall_distance.cpp index 2b25d23bd9c2..ccc99cd1fad2 100644 --- a/Common/src/fem_wall_distance.cpp +++ b/Common/src/fem_wall_distance.cpp @@ -48,7 +48,6 @@ std::unique_ptr CMeshFEM_DG::ComputeViscousWallADT(const CConfig vector elemIDs; vector VTK_TypeElem; vector markerIDs; - vector dummyRoughs; /* Loop over the boundary markers. */ for(unsigned short iMarker=0; iMarker CMeshFEM_DG::ComputeViscousWallADT(const CConfig unsigned short ii = 0; for(unsigned short j=0; j CMeshFEM_DG::ComputeViscousWallADT(const CConfig /* Build the ADT. */ std::unique_ptr WallADT(new CADTElemClass(nDim, surfaceCoor, surfaceConn, VTK_TypeElem, - markerIDs, elemIDs, dummyRoughs, true)); + markerIDs, elemIDs, true)); return WallADT; @@ -204,9 +202,8 @@ void CMeshFEM_DG::SetWallDistance(const CConfig *config, CADTElemClass *WallADT) unsigned long elemID; int rankID; su2double dist; - su2double dumRough; - WallADT->DetermineNearestElement(coor, dist, markerID, elemID, rankID, dumRough); + WallADT->DetermineNearestElement(coor, dist, markerID, elemID, rankID); volElem[l].wallDistance[i] = dist; } @@ -239,9 +236,8 @@ void CMeshFEM_DG::SetWallDistance(const CConfig *config, CADTElemClass *WallADT) unsigned long elemID; int rankID; su2double dist; - su2double dumRough; - WallADT->DetermineNearestElement(coor, dist, markerID, elemID, rankID, dumRough); + WallADT->DetermineNearestElement(coor, dist, markerID, elemID, rankID); volElem[l].wallDistanceSolDOFs[i] = dist; } @@ -274,9 +270,8 @@ void CMeshFEM_DG::SetWallDistance(const CConfig *config, CADTElemClass *WallADT) unsigned long elemID; int rankID; su2double dist; - su2double dumRough; - WallADT->DetermineNearestElement(coor, dist, markerID, elemID, rankID, dumRough); + WallADT->DetermineNearestElement(coor, dist, markerID, elemID, rankID); matchingFaces[l].wallDistance[i] = dist; } @@ -328,9 +323,8 @@ void CMeshFEM_DG::SetWallDistance(const CConfig *config, CADTElemClass *WallADT) unsigned long elemID; int rankID; su2double dist; - su2double dumRough; - WallADT->DetermineNearestElement(coor, dist, markerID, elemID, rankID, dumRough); + WallADT->DetermineNearestElement(coor, dist, markerID, elemID, rankID); surfElem[l].wallDistance[i] = dist; } diff --git a/Common/src/geometry/CPhysicalGeometry.cpp b/Common/src/geometry/CPhysicalGeometry.cpp index a81775890db2..d9eaec43415a 100644 --- a/Common/src/geometry/CPhysicalGeometry.cpp +++ b/Common/src/geometry/CPhysicalGeometry.cpp @@ -11504,9 +11504,6 @@ std::unique_ptr CPhysicalGeometry::ComputeViscousWallADT(const CC vector elemIDs; vector VTK_TypeElem; vector markerIDs; - vector roughheights; - string Marker_Tag; - su2double roughness; /* Loop over the boundary markers. */ @@ -11516,9 +11513,6 @@ std::unique_ptr CPhysicalGeometry::ComputeViscousWallADT(const CC /* Check for a viscous wall. */ if( config->GetViscous_Wall(iMarker)) { - Marker_Tag = config->GetMarker_All_TagBound(iMarker); - roughness = config->GetWall_RoughnessHeight(Marker_Tag); - /* Loop over the surface elements of this marker. */ for(unsigned long iElem=0; iElem < nElem_Bound[iMarker]; iElem++) { @@ -11538,8 +11532,6 @@ std::unique_ptr CPhysicalGeometry::ComputeViscousWallADT(const CC markerIDs.push_back(iMarker); VTK_TypeElem.push_back(VTK_Type); elemIDs.push_back(iElem); - roughheights.push_back(roughness); - for (unsigned short iNode = 0; iNode < nDOFsPerElem; iNode++) surfaceConn.push_back(bound[iMarker][iElem]->GetNode(iNode)); } @@ -11574,7 +11566,7 @@ std::unique_ptr CPhysicalGeometry::ComputeViscousWallADT(const CC /*--------------------------------------------------------------------------*/ std::unique_ptr WallADT(new CADTElemClass(nDim, surfaceCoor, surfaceConn, VTK_TypeElem, - markerIDs, elemIDs, roughheights, true)); + markerIDs, elemIDs, true)); return WallADT; @@ -11600,13 +11592,84 @@ void CPhysicalGeometry::SetWallDistance(const CConfig *config, CADTElemClass *Wa su2double dist; su2double localRoughness; - WallADT->DetermineNearestElement(nodes->GetCoord(iPoint), dist, markerID, - elemID, rankID, localRoughness); + WallADT->DetermineNearestElement(nodes->GetCoord(iPoint), dist, markerID, elemID, rankID); - /*WallADT->DetermineNearestElement(nodes->GetCoord(iPoint), dist, markerID, elemID, rankID);*/ nodes->SetWall_Distance(iPoint, min(dist,nodes->GetWall_Distance(iPoint))); + localRoughness = config->GetLocalRoughness(rankID, markerID); nodes->SetRoughnessHeight(iPoint, localRoughness); } + } // end SU2_OMP_PARALLEL } + +void CPhysicalGeometry::SetGlobalMarkerRoughness(CConfig *config) { + +unsigned short iMarker; +#ifndef HAVE_MPI + + int *local_displ = new int [1]; + int size = 1; + unsigned short nMarker_All = config->GetnMarker_All(); + local_displ[0] = 0; + config->SetGlobalMarkerArray(local_displ,size) + + su2double *localRough = new su2double [nMarker_All]; + + if (config->GetnRoughWall() > 0) + for (iMarker = 0; iMarker < nMarker_All; iMarker++) + localRough[iMarker] = config->GetWall_RoughnessHeight(config->GetMarker_All_TagBound(iMarker)); + else + for (iMarker = 0; iMarker < nMarker_All; iMarker++) + localRough[iMarker] = 0.0 ; + + config->SetGlobalRoughnessArray(localRough, nMarker_All); + +#else + int rank, size; + unsigned short nMarker_All = config->GetnMarker_All(); + SU2_MPI::Comm_rank(MPI_COMM_WORLD, &rank); + SU2_MPI::Comm_size(MPI_COMM_WORLD, &size); + + int *displs = new int [size]; + int* recvCounts = new int [size];//, displs(size); + int sizeLocal = (int) nMarker_All; // number of local markers + + //Communicate size of local marker array and make an array large enough to hold all data + SU2_MPI::Allgather(&sizeLocal, 1, MPI_INT, recvCounts, 1, + MPI_INT, MPI_COMM_WORLD); + + SU2_MPI::Barrier(MPI_COMM_WORLD); // Not sure if needed + + // displacements based on size on each rank + displs[0] = 0; + for(int i=1; iSetGlobalMarkerArray(displs, size); + + //total size + int sizeGlobal = displs[size-1] + recvCounts[size-1]; + + /*--- Allocate local and global arrays to hold roughness. ---*/ + su2double *localRough = new su2double [nMarker_All]; // local number of markers + su2double *globalRough = new su2double[sizeGlobal]; // all markers including send recieve + + if (config->GetnRoughWall() > 0) + for (iMarker = 0; iMarker < nMarker_All; iMarker++) + localRough[iMarker] = config->GetWall_RoughnessHeight(config->GetMarker_All_TagBound(iMarker)); + else + for (iMarker = 0; iMarker < nMarker_All; iMarker++) + localRough[iMarker] = 0.0 ; + + SU2_MPI::Allgatherv( localRough, sizeLocal, MPI_DOUBLE, globalRough , + recvCounts, displs, MPI_DOUBLE, + MPI_COMM_WORLD); + + SU2_MPI::Barrier(MPI_COMM_WORLD); // Not sure if needed + + /*--- Set the global array of roughness per marker. ---*/ + config->SetGlobalRoughnessArray(globalRough, sizeGlobal); + +#endif +} diff --git a/Common/src/geometry_structure_fem_part.cpp b/Common/src/geometry_structure_fem_part.cpp index e32935416d1c..f2c6c3bdebd8 100644 --- a/Common/src/geometry_structure_fem_part.cpp +++ b/Common/src/geometry_structure_fem_part.cpp @@ -3391,7 +3391,6 @@ void CPhysicalGeometry::DetermineDonorElementsWallFunctions(CConfig *config) { vector subElementIDInParent; vector VTK_TypeElem; vector elemConn; - vector dummyRough; /* Loop over the local volume elements to create the connectivity of the linear sub-elements. */ @@ -3443,7 +3442,6 @@ void CPhysicalGeometry::DetermineDonorElementsWallFunctions(CConfig *config) { parentElement.push_back(elem[l]->GetGlobalElemID()); subElementIDInParent.push_back(jj); VTK_TypeElem.push_back(VTK_Type[i]); - dummyRough.push_back(0.0); for(unsigned short k=0; kGetNode(connSubElems[i][kk]); @@ -3468,7 +3466,7 @@ void CPhysicalGeometry::DetermineDonorElementsWallFunctions(CConfig *config) { /* Build the local ADT. */ CADTElemClass localVolumeADT(nDim, volCoor, elemConn, VTK_TypeElem, - subElementIDInParent, parentElement, dummyRough, false); + subElementIDInParent, parentElement, false); /* Release the memory of the vectors used to build the ADT. To make sure that all the memory is deleted, the swap function is used. */ @@ -3477,7 +3475,6 @@ void CPhysicalGeometry::DetermineDonorElementsWallFunctions(CConfig *config) { vector().swap(parentElement); vector().swap(elemConn); vector().swap(volCoor); - vector().swap(dummyRough); /*--------------------------------------------------------------------------*/ /*--- Step 3. Search for donor elements at the exchange locations in ---*/ diff --git a/SU2_CFD/src/drivers/CDriver.cpp b/SU2_CFD/src/drivers/CDriver.cpp index 897f30705115..c68db37ba52f 100644 --- a/SU2_CFD/src/drivers/CDriver.cpp +++ b/SU2_CFD/src/drivers/CDriver.cpp @@ -853,6 +853,10 @@ void CDriver::Geometrical_Preprocessing_FVM(CConfig *config, CGeometry **&geomet geometry[MESH_0]->ComputeMeshQualityStatistics(config); } + /*--- Store marker list and roughness in a global array ---*/ + + geometry[MESH_0]->SetGlobalMarkerRoughness(config); + geometry[MESH_0]->SetMGLevel(MESH_0); if ((config->GetnMGLevels() != 0) && (rank == MASTER_NODE)) cout << "Setting the multigrid structure." << endl; From 44ee968304c896716b545f9ef59f932b085e68d2 Mon Sep 17 00:00:00 2001 From: koodlyakshay Date: Tue, 19 May 2020 12:38:58 +0200 Subject: [PATCH 17/29] Fix spaces and errors. --- Common/include/adt_structure.hpp | 2 +- Common/include/geometry/CGeometry.hpp | 3 +++ Common/include/geometry/CPhysicalGeometry.hpp | 5 ++++- Common/src/adt_structure.cpp | 20 +++++++++---------- Common/src/fem_wall_distance.cpp | 5 ----- Common/src/geometry/CPhysicalGeometry.cpp | 2 +- 6 files changed, 19 insertions(+), 18 deletions(-) diff --git a/Common/include/adt_structure.hpp b/Common/include/adt_structure.hpp index 8fdbbeff508e..8cfdddfa5fce 100644 --- a/Common/include/adt_structure.hpp +++ b/Common/include/adt_structure.hpp @@ -288,7 +288,7 @@ class CADTElemClass : public CADTBaseClass { vector localElemIDs; /*!< \brief Vector, which contains the local element ID's of the elements in the ADT. */ vector ranksOfElems; /*!< \brief Vector, which contains the ranks - of the elements in the ADT. */ + of the elements in the ADT. */ #ifdef HAVE_OMP vector >BBoxTargets; /*!< \brief Vector, used to store possible bounding box candidates during the nearest element search. */ diff --git a/Common/include/geometry/CGeometry.hpp b/Common/include/geometry/CGeometry.hpp index fa16a2a3d1b2..6833a8e030e4 100644 --- a/Common/include/geometry/CGeometry.hpp +++ b/Common/include/geometry/CGeometry.hpp @@ -1693,6 +1693,9 @@ class CGeometry { */ static void ComputeWallDistance(const CConfig * const *config_container, CGeometry ****geometry_container); + /*! + * \brief Set roughness values for markers in a global array. + */ inline virtual void SetGlobalMarkerRoughness(CConfig *config) { } }; diff --git a/Common/include/geometry/CPhysicalGeometry.hpp b/Common/include/geometry/CPhysicalGeometry.hpp index ec013389fdf2..87b21a446928 100644 --- a/Common/include/geometry/CPhysicalGeometry.hpp +++ b/Common/include/geometry/CPhysicalGeometry.hpp @@ -803,7 +803,10 @@ class CPhysicalGeometry final : public CGeometry { nodes->SetWall_Distance(iPoint, val); } } - + + /*! + * \brief Set roughness values for markers in a global array. + */ void SetGlobalMarkerRoughness(CConfig *config) override; }; diff --git a/Common/src/adt_structure.cpp b/Common/src/adt_structure.cpp index 192199ef31f6..e6f9446ecb15 100644 --- a/Common/src/adt_structure.cpp +++ b/Common/src/adt_structure.cpp @@ -590,11 +590,11 @@ CADTElemClass::CADTElemClass(unsigned short val_nDim, int rank; SU2_MPI::Comm_rank(MPI_COMM_WORLD, &rank); - coorPoints = val_coor; - elemConns = val_connElem; - elemVTK_Type = val_VTKElem; - localMarkers = val_markerID; - localElemIDs = val_elemID; + coorPoints = val_coor; + elemConns = val_connElem; + elemVTK_Type = val_VTKElem; + localMarkers = val_markerID; + localElemIDs = val_elemID; ranksOfElems.assign(elemVTK_Type.size(), rank); } @@ -602,11 +602,11 @@ CADTElemClass::CADTElemClass(unsigned short val_nDim, /*--- Sequential mode. Copy the data from the arguments into the member variables and set the ranks to MASTER_NODE. ---*/ - coorPoints = val_coor; - elemConns = val_connElem; - elemVTK_Type = val_VTKElem; - localMarkers = val_markerID; - localElemIDs = val_elemID; + coorPoints = val_coor; + elemConns = val_connElem; + elemVTK_Type = val_VTKElem; + localMarkers = val_markerID; + localElemIDs = val_elemID; ranksOfElems.assign(elemVTK_Type.size(), MASTER_NODE); diff --git a/Common/src/fem_wall_distance.cpp b/Common/src/fem_wall_distance.cpp index ccc99cd1fad2..9a22e0267c2e 100644 --- a/Common/src/fem_wall_distance.cpp +++ b/Common/src/fem_wall_distance.cpp @@ -174,7 +174,6 @@ void CMeshFEM_DG::SetWallDistance(su2double val){ void CMeshFEM_DG::SetWallDistance(const CConfig *config, CADTElemClass *WallADT){ - /*--------------------------------------------------------------------------*/ /*--- Step 3: Determine the wall distance of the integration points of ---*/ /*--- locally owned volume elements. ---*/ @@ -202,7 +201,6 @@ void CMeshFEM_DG::SetWallDistance(const CConfig *config, CADTElemClass *WallADT) unsigned long elemID; int rankID; su2double dist; - WallADT->DetermineNearestElement(coor, dist, markerID, elemID, rankID); volElem[l].wallDistance[i] = dist; @@ -236,7 +234,6 @@ void CMeshFEM_DG::SetWallDistance(const CConfig *config, CADTElemClass *WallADT) unsigned long elemID; int rankID; su2double dist; - WallADT->DetermineNearestElement(coor, dist, markerID, elemID, rankID); volElem[l].wallDistanceSolDOFs[i] = dist; @@ -270,7 +267,6 @@ void CMeshFEM_DG::SetWallDistance(const CConfig *config, CADTElemClass *WallADT) unsigned long elemID; int rankID; su2double dist; - WallADT->DetermineNearestElement(coor, dist, markerID, elemID, rankID); matchingFaces[l].wallDistance[i] = dist; @@ -323,7 +319,6 @@ void CMeshFEM_DG::SetWallDistance(const CConfig *config, CADTElemClass *WallADT) unsigned long elemID; int rankID; su2double dist; - WallADT->DetermineNearestElement(coor, dist, markerID, elemID, rankID); surfElem[l].wallDistance[i] = dist; diff --git a/Common/src/geometry/CPhysicalGeometry.cpp b/Common/src/geometry/CPhysicalGeometry.cpp index d9eaec43415a..8e555ee5b8f7 100644 --- a/Common/src/geometry/CPhysicalGeometry.cpp +++ b/Common/src/geometry/CPhysicalGeometry.cpp @@ -11612,7 +11612,7 @@ unsigned short iMarker; int size = 1; unsigned short nMarker_All = config->GetnMarker_All(); local_displ[0] = 0; - config->SetGlobalMarkerArray(local_displ,size) + config->SetGlobalMarkerArray(local_displ,size); su2double *localRough = new su2double [nMarker_All]; From e73cba1f30bcb4134e0b1bed4d9a22d38d93bc0b Mon Sep 17 00:00:00 2001 From: koodlyakshay Date: Tue, 19 May 2020 15:11:47 +0200 Subject: [PATCH 18/29] Clean up comments and spaces. --- Common/include/CConfig.hpp | 19 +++++++++++++++++-- Common/include/adt_structure.hpp | 1 - 2 files changed, 17 insertions(+), 3 deletions(-) diff --git a/Common/include/CConfig.hpp b/Common/include/CConfig.hpp index 1f46a1b394bd..08ac3cfdccdd 100644 --- a/Common/include/CConfig.hpp +++ b/Common/include/CConfig.hpp @@ -5100,10 +5100,25 @@ class CConfig { */ su2double GetCoeff_ObjChainRule(unsigned short iVar) const { return Obj_ChainRuleCoeff[iVar]; } - /*--- New roughness mpi comm stuff akshay ----*/ + /*! + * \brief Get roughness height of a marker. + * \param [in] rankID rank where the marker is located. + * \param [in] ID of the marker on rank + */ su2double GetLocalRoughness(int rankID, unsigned short markerID) const; - + + /*! + * \brief Set the displacements for each rank in global array. + * \param [in] global displacemts in order of rank. + * \param [in] total number of ranks + */ void SetGlobalMarkerArray(int *global_displ, int size) ; + + /*! + * \brief Set the roughness for each marker in global array. + * \param [in] global list of roughness for every marker. + * \param [in] total number of marker including send recieves. + */ void SetGlobalRoughnessArray(su2double *global_rough, int sizeGlobal) ; /*! diff --git a/Common/include/adt_structure.hpp b/Common/include/adt_structure.hpp index 8cfdddfa5fce..1e946263046e 100644 --- a/Common/include/adt_structure.hpp +++ b/Common/include/adt_structure.hpp @@ -295,7 +295,6 @@ class CADTElemClass : public CADTBaseClass { #else array,1> BBoxTargets; #endif - public: /*! * \brief Constructor of the class. From bb90518a6a7aad803b5d3ece2331e059fa91a876 Mon Sep 17 00:00:00 2001 From: koodlyakshay Date: Thu, 11 Jun 2020 16:48:01 +0200 Subject: [PATCH 19/29] Simplify input method. --- Common/include/CConfig.hpp | 5 +- Common/src/CConfig.cpp | 120 ++++++++++++++-------- Common/src/geometry/CPhysicalGeometry.cpp | 4 - SU2_CFD/src/solvers/CTurbSASolver.cpp | 1 - 4 files changed, 80 insertions(+), 50 deletions(-) diff --git a/Common/include/CConfig.hpp b/Common/include/CConfig.hpp index 08ac3cfdccdd..20fd57df74e9 100644 --- a/Common/include/CConfig.hpp +++ b/Common/include/CConfig.hpp @@ -254,6 +254,7 @@ class CConfig { *Marker_Outlet, /*!< \brief Outlet flow markers. */ *Marker_Isothermal, /*!< \brief Isothermal wall markers. */ *Marker_HeatFlux, /*!< \brief Constant heat flux wall markers. */ + *Marker_RoughWall, /*!< \brief Constant heat flux wall markers. */ *Marker_EngineInflow, /*!< \brief Engine Inflow flow markers. */ *Marker_EngineExhaust, /*!< \brief Engine Exhaust flow markers. */ *Marker_Clamped, /*!< \brief Clamped markers. */ @@ -927,7 +928,7 @@ class CConfig { nMarkerPitching_Phase, /*!< \brief Number of values provided for pitching phase offset of marker. */ nMarkerPlunging_Omega, /*!< \brief Number of values provided for angular frequency of marker. */ nMarkerPlunging_Ampl, /*!< \brief Number of values provided for plunging amplitude of marker. */ - nRoughWall; /*!< \brief Number of rough walls. */ + nRough_Wall; /*!< \brief Number of rough walls. */ su2double *Omega_HB; /*!< \brief Frequency for Harmonic Balance Operator (in rad/s). */ unsigned short nOmega_HB, /*!< \brief Number of frequencies in Harmonic Balance Operator. */ @@ -2937,7 +2938,7 @@ class CConfig { * \brief Get the total number of rough markers. * \return Total number of heat flux markers. */ - unsigned short GetnRoughWall(void) const { return nRoughWall; } + unsigned short GetnRoughWall(void) const { return nRough_Wall; } /*! * \brief Get the total number of objectives in kind_objective list diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index 53e43ce010de..eb482307b574 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -819,7 +819,7 @@ void CConfig::SetPointersNull(void) { Marker_Inlet = nullptr; Marker_Outlet = nullptr; Marker_Supersonic_Inlet = nullptr; Marker_Supersonic_Outlet= nullptr; Marker_Isothermal = nullptr; Marker_HeatFlux = nullptr; Marker_EngineInflow = nullptr; - Marker_Load = nullptr; Marker_Disp_Dir = nullptr; + Marker_Load = nullptr; Marker_Disp_Dir = nullptr; Marker_RoughWall = nullptr; Marker_EngineExhaust = nullptr; Marker_Displacement = nullptr; Marker_Load = nullptr; Marker_Load_Dir = nullptr; Marker_Load_Sine = nullptr; Marker_Clamped = nullptr; Marker_FlowLoad = nullptr; Marker_Internal = nullptr; @@ -1519,10 +1519,13 @@ void CConfig::SetConfig_Options() { Format: ( Heat flux marker, wall heat flux (static), ... ) \ingroup Config*/ addStringDoubleListOption("MARKER_HEATFLUX", nMarker_HeatFlux, Marker_HeatFlux, Heat_Flux); /*!\brief WALL_TYPE \n DESCRIPTION: List of wall types - SMOOTH or ROUGH (All walls are SMOOTH by default). List length must match number of wall markers or shouldn't appear at all. \ingroup Config*/ - addEnumListOption("WALL_TYPE", nWall_Types, Kind_Wall, WallType_Map); + //addEnumListOption("WALL_TYPE", nWall_Types, Kind_Wall, WallType_Map); /*!\brief Wall roughness height default value = 0 (smooth). \ingroup Config */ default_roughness[0] = 0.0; //default_roughness[1] = 0.0; - addDoubleListOption("WALL_ROUGHNESS", nRoughWall, Roughness_Height); + //addDoubleListOption("WALL_ROUGHNESS", nRoughWall, Roughness_Height); + /*!\brief MARKER_ROUGHWALL \n DESCRIPTION: Specified roughness heights at wall boundary marker(s) + Format: ( Wall marker, roughness_height (static), ... ) \ingroup Config*/ + addStringDoubleListOption("MARKER_ROUGHWALL", nRough_Wall, Marker_RoughWall, Roughness_Height); /*!\brief MARKER_ENGINE_INFLOW \n DESCRIPTION: Engine inflow boundary marker(s) Format: ( nacelle inflow marker, fan face Mach, ... ) \ingroup Config*/ addStringDoubleListOption("MARKER_ENGINE_INFLOW", nMarker_EngineInflow, Marker_EngineInflow, EngineInflow_Target); @@ -4645,46 +4648,63 @@ void CConfig::SetPostprocessing(unsigned short val_software, unsigned short val_ } } - /*--- Check that if the wall roughness array are compatible and set deafult values if needed. ---*/ + /*--- Check that if the wall roughness array are compatible and set deafult values if needed. ---*/ + if ((nMarker_HeatFlux > 0) || (nMarker_Isothermal > 0) || (nMarker_CHTInterface > 0)) { + /*--- If no roughness is specified all walls are assumed to be smooth. ---*/ + if (nRough_Wall == 0) { + nRough_Wall = nMarker_HeatFlux + nMarker_Isothermal + nMarker_CHTInterface; + Roughness_Height = new su2double [nRough_Wall]; + Kind_Wall = new unsigned short [nRough_Wall]; + for (unsigned short iMarker = 0; iMarker < nMarker_HeatFlux; iMarker++) { + Roughness_Height[iMarker] = 0.0; + Kind_Wall[iMarker] = SMOOTH; + } + for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) { + Roughness_Height[nMarker_HeatFlux + iMarker] = 0.0; + Kind_Wall[nMarker_HeatFlux + iMarker] = SMOOTH; + } + for (iMarker = 0; iMarker < nMarker_CHTInterface; iMarker++) { + Roughness_Height[nMarker_HeatFlux + nMarker_Isothermal + iMarker] = 0.0; + Kind_Wall[nMarker_HeatFlux + nMarker_Isothermal + iMarker] = SMOOTH; + } + } else { /*--- Check name of the marker and assign the corresponding roughness. ---*/ + /*--- Store roughness heights in a temp array. ---*/ + vector temp_rough; + for (iMarker = 0; iMarker < nRough_Wall; iMarker++) + temp_rough.push_back(Roughness_Height[iMarker]); + + /*--- Reallocate the roughness arrays in case not all walls are rough. ---*/ + delete Roughness_Height; + delete Kind_Wall; + unsigned short nWall = nMarker_HeatFlux + nMarker_Isothermal + nMarker_CHTInterface; + unsigned short jMarker; + Roughness_Height = new su2double [nWall]; + Kind_Wall = new unsigned short [nWall]; + /*--- Initialize everything to smooth. ---*/ + for (iMarker = 0; iMarker < nWall; iMarker++) { + Roughness_Height[iMarker] = 0.0; + Kind_Wall[iMarker] = SMOOTH; + } - if ((nMarker_HeatFlux > 0) || (nMarker_Isothermal > 0)) { - /*--- If the config option is not declared, then assume smooth walls and set roughness to zero. ---*/ - if (nWall_Types != nMarker_HeatFlux + nMarker_Isothermal) { - /*--- If this option is not used at all, assume all walls are smooth and set roughness height to zero. ---*/ - if (nWall_Types == 0) { - nWall_Types = nMarker_HeatFlux + nMarker_Isothermal; - Roughness_Height = new su2double [nWall_Types]; - Kind_Wall = new unsigned short [nWall_Types]; - for (unsigned short iMarker = 0; iMarker < nMarker_HeatFlux; iMarker++) { - Roughness_Height[iMarker] = 0.0; - Kind_Wall[iMarker] = SMOOTH; - } - for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) { - Roughness_Height[nMarker_HeatFlux + iMarker] = 0.0; - Kind_Wall[nMarker_HeatFlux + iMarker] = SMOOTH; - } - } - else - SU2_MPI::Error("Mismatch in number of wall markers. Specify type of wall for all markers (even if smooth).", CURRENT_FUNCTION); - } else { - if (nRoughWall != nWall_Types) - SU2_MPI::Error("Mismatch in number of roughness height definition and wall type definition.", CURRENT_FUNCTION); - else { - vector temp_rough; - vector temp_kindrough; - for (unsigned short iMarker = 0; iMarker < nWall_Types; iMarker++) { - temp_rough.push_back(Roughness_Height[iMarker]); - temp_kindrough.push_back(Kind_Wall[iMarker]); - } - for (iMarker = 0; iMarker < nMarker_HeatFlux; iMarker++) { - Roughness_Height[iMarker] = temp_rough[iMarker]; - Kind_Wall[iMarker] = temp_kindrough[iMarker]; - } - for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) { - Roughness_Height[nMarker_HeatFlux + iMarker] = temp_rough[nMarker_HeatFlux + iMarker]; - Kind_Wall[nMarker_HeatFlux + iMarker] = temp_kindrough[nMarker_HeatFlux + iMarker]; - } + /*--- Look through heat flux, isothermal and cht_interface markers and assign proper values. ---*/ + for (iMarker = 0; iMarker < nRough_Wall; iMarker++) { + for (jMarker = 0; jMarker < nMarker_HeatFlux; jMarker++) + if (Marker_HeatFlux[jMarker].compare(Marker_RoughWall[iMarker]) == 0) + Roughness_Height[jMarker] = temp_rough[iMarker]; + + for (jMarker = 0; jMarker < nMarker_Isothermal; jMarker++) + if (Marker_Isothermal[jMarker].compare(Marker_RoughWall[iMarker]) == 0) + Roughness_Height[nMarker_HeatFlux + jMarker] = temp_rough[iMarker]; + + for (jMarker = 0; jMarker < nMarker_CHTInterface; jMarker++) + if (Marker_CHTInterface[jMarker].compare(Marker_RoughWall[iMarker]) == 0) + Roughness_Height[nMarker_HeatFlux + nMarker_Isothermal + jMarker] = temp_rough[iMarker]; } + + /*--- Update kind_wall when a non zero roughness value is specified. ---*/ + for (iMarker = 0; iMarker < nWall; iMarker++) + if (Roughness_Height[iMarker] != 0.0) + Kind_Wall[iMarker] = ROUGH; } } @@ -8754,9 +8774,16 @@ unsigned short CConfig::GetKindWall(string val_marker) const { if (flag == -1) for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) if (Marker_Isothermal[iMarker] == val_marker) { - flag = iMarker; + flag = nMarker_HeatFlux + iMarker; break; } + if (flag == -1) { + for (iMarker = 0; iMarker < nMarker_CHTInterface; iMarker++) + if (Marker_CHTInterface[iMarker] == val_marker) { + flag = nMarker_HeatFlux + nMarker_Isothermal + iMarker; + break; + } + } if (flag == -1) return 1; else return Kind_Wall[flag]; } @@ -8778,7 +8805,14 @@ su2double CConfig::GetWall_RoughnessHeight(string val_marker) const { if (flag == -1) { for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) if (Marker_Isothermal[iMarker] == val_marker) { - flag = iMarker; + flag = nMarker_HeatFlux + iMarker; + break; + } + } + if (flag == -1) { + for (iMarker = 0; iMarker < nMarker_CHTInterface; iMarker++) + if (Marker_CHTInterface[iMarker] == val_marker) { + flag = nMarker_HeatFlux + nMarker_Isothermal + iMarker; break; } } diff --git a/Common/src/geometry/CPhysicalGeometry.cpp b/Common/src/geometry/CPhysicalGeometry.cpp index 8e555ee5b8f7..8c38c0c90429 100644 --- a/Common/src/geometry/CPhysicalGeometry.cpp +++ b/Common/src/geometry/CPhysicalGeometry.cpp @@ -11639,8 +11639,6 @@ unsigned short iMarker; SU2_MPI::Allgather(&sizeLocal, 1, MPI_INT, recvCounts, 1, MPI_INT, MPI_COMM_WORLD); - SU2_MPI::Barrier(MPI_COMM_WORLD); // Not sure if needed - // displacements based on size on each rank displs[0] = 0; for(int i=1; iSetGlobalRoughnessArray(globalRough, sizeGlobal); diff --git a/SU2_CFD/src/solvers/CTurbSASolver.cpp b/SU2_CFD/src/solvers/CTurbSASolver.cpp index 4cc7eddadc3f..b3dc665e5a37 100644 --- a/SU2_CFD/src/solvers/CTurbSASolver.cpp +++ b/SU2_CFD/src/solvers/CTurbSASolver.cpp @@ -599,7 +599,6 @@ void CTurbSASolver::BC_Isothermal_Wall(CGeometry *geometry, CSolver **solver_con /*--- Includes 1 in the diagonal ---*/ Jacobian.DeleteValsRowi(iPoint); - } else { /*--- For rough walls, the boundary condition is given by * (\frac{\partial \nu}{\partial n})_wall = \frac{\nu}{0.03*k_s} From 999a6313ddb12e826b441638c26a415d561c4aa6 Mon Sep 17 00:00:00 2001 From: koodlyakshay Date: Fri, 3 Jul 2020 17:44:36 +0200 Subject: [PATCH 20/29] Move global arrays from config into geometry and update option name. --- Common/include/CConfig.hpp | 25 +----- Common/include/geometry/CPhysicalGeometry.hpp | 17 ++++ Common/src/CConfig.cpp | 82 ++++++++----------- Common/src/geometry/CPhysicalGeometry.cpp | 39 ++++++--- externals/meson | 2 +- externals/ninja | 2 +- 6 files changed, 81 insertions(+), 86 deletions(-) diff --git a/Common/include/CConfig.hpp b/Common/include/CConfig.hpp index 118eb27af156..467e6822d4ca 100644 --- a/Common/include/CConfig.hpp +++ b/Common/include/CConfig.hpp @@ -309,7 +309,6 @@ class CConfig { su2double *Isothermal_Temperature; /*!< \brief Specified isothermal wall temperatures (static). */ su2double *Heat_Flux; /*!< \brief Specified wall heat fluxes. */ su2double *Roughness_Height; /*!< \brief Equivalent sand grain roughness for the marker according to config file. */ - su2double *GlobalRoughness_Height; /*!< \brief Equivalent sand grain roughness for markers globally including S-R. */ su2double *Displ_Value; /*!< \brief Specified displacement for displacement boundaries. */ su2double *Load_Value; /*!< \brief Specified force for load boundaries. */ su2double *Damper_Constant; /*!< \brief Specified constant for damper boundaries. */ @@ -390,7 +389,6 @@ class CConfig { *Marker_CfgFile_KindBC; /*!< \brief Global index for boundaries using config file. */ short *Marker_All_SendRecv; /*!< \brief Information about if the boundary is sended (+), received (-). */ short *Marker_All_PerBound; /*!< \brief Global index for periodic bc using the grid information. */ - int *GlobalMarkerStorageDispl; /*!< \brief storage of markers globally. */ unsigned long nExtIter; /*!< \brief Number of external iterations. */ unsigned long ExtIter; /*!< \brief Current external iteration number. */ @@ -5097,28 +5095,7 @@ class CConfig { * Gradients are w.r.t density, velocity[3], and pressure. when 2D gradient w.r.t. 3rd component of velocity set to 0. */ su2double GetCoeff_ObjChainRule(unsigned short iVar) const { return Obj_ChainRuleCoeff[iVar]; } - - /*! - * \brief Get roughness height of a marker. - * \param [in] rankID rank where the marker is located. - * \param [in] ID of the marker on rank - */ - su2double GetLocalRoughness(int rankID, unsigned short markerID) const; - - /*! - * \brief Set the displacements for each rank in global array. - * \param [in] global displacemts in order of rank. - * \param [in] total number of ranks - */ - void SetGlobalMarkerArray(int *global_displ, int size) ; - - /*! - * \brief Set the roughness for each marker in global array. - * \param [in] global list of roughness for every marker. - * \param [in] total number of marker including send recieves. - */ - void SetGlobalRoughnessArray(su2double *global_rough, int sizeGlobal) ; - + /*! * \author H. Kline * \brief Get the flag indicating whether to comput a combined objective. diff --git a/Common/include/geometry/CPhysicalGeometry.hpp b/Common/include/geometry/CPhysicalGeometry.hpp index f2ff8a12eb88..b5b8252776c1 100644 --- a/Common/include/geometry/CPhysicalGeometry.hpp +++ b/Common/include/geometry/CPhysicalGeometry.hpp @@ -31,6 +31,7 @@ #include "meshreader/CMeshReaderFVM.hpp" #include "../toolboxes/C2DContainer.hpp" + /*! * \class CPhysicalGeometry * \brief Class for reading a defining the primal grid which is read from the grid file in .su2 or .cgns format. @@ -104,6 +105,8 @@ class CPhysicalGeometry final : public CGeometry { unsigned long *Elem_ID_Line_Linear{nullptr}; unsigned long *Elem_ID_BoundTria_Linear{nullptr}; unsigned long *Elem_ID_BoundQuad_Linear{nullptr}; + int *GlobalMarkerStorageDispl{nullptr}; + su2double *GlobalRoughness_Height{nullptr}; public: /*--- This is to suppress Woverloaded-virtual, omitting it has no negative impact. ---*/ @@ -808,5 +811,19 @@ class CPhysicalGeometry final : public CGeometry { * \brief Set roughness values for markers in a global array. */ void SetGlobalMarkerRoughness(CConfig *config) override; + + /*! + * \brief Set the displacements for each rank in global array. + * \param [in] global displacemts in order of rank. + * \param [in] total number of ranks + */ + void SetGlobalMarkerArray(int *global_displ, int size) ; + + /*! + * \brief Set the roughness for each marker in global array. + * \param [in] global list of roughness for every marker. + * \param [in] total number of marker including send recieves. + */ + void SetGlobalRoughnessArray(su2double *global_rough, int sizeGlobal) ; }; diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index 2c977c77901f..bf9c624279ed 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -834,7 +834,7 @@ void CConfig::SetPointersNull(void) { Isothermal_Temperature = nullptr; Heat_Flux = nullptr; Displ_Value = nullptr; Load_Value = nullptr; FlowLoad_Value = nullptr; Damper_Constant = nullptr; Wall_Emissivity = nullptr; - Roughness_Height = nullptr; GlobalMarkerStorageDispl = nullptr; GlobalRoughness_Height = nullptr; + Roughness_Height = nullptr; /*--- Inlet Outlet Boundary Condition settings ---*/ @@ -1497,14 +1497,9 @@ void CConfig::SetConfig_Options() { /*!\brief MARKER_HEATFLUX \n DESCRIPTION: Specified heat flux wall boundary marker(s) Format: ( Heat flux marker, wall heat flux (static), ... ) \ingroup Config*/ addStringDoubleListOption("MARKER_HEATFLUX", nMarker_HeatFlux, Marker_HeatFlux, Heat_Flux); - /*!\brief WALL_TYPE \n DESCRIPTION: List of wall types - SMOOTH or ROUGH (All walls are SMOOTH by default). List length must match number of wall markers or shouldn't appear at all. \ingroup Config*/ - //addEnumListOption("WALL_TYPE", nWall_Types, Kind_Wall, WallType_Map); - /*!\brief Wall roughness height default value = 0 (smooth). \ingroup Config */ - default_roughness[0] = 0.0; //default_roughness[1] = 0.0; - //addDoubleListOption("WALL_ROUGHNESS", nRoughWall, Roughness_Height); - /*!\brief MARKER_ROUGHWALL \n DESCRIPTION: Specified roughness heights at wall boundary marker(s) + /*!\brief WALL_ROUGHNESS \n DESCRIPTION: Specified roughness heights at wall boundary marker(s) Format: ( Wall marker, roughness_height (static), ... ) \ingroup Config*/ - addStringDoubleListOption("MARKER_ROUGHWALL", nRough_Wall, Marker_RoughWall, Roughness_Height); + addStringDoubleListOption("WALL_ROUGHNESS", nRough_Wall, Marker_RoughWall, Roughness_Height); /*!\brief MARKER_ENGINE_INFLOW \n DESCRIPTION: Engine inflow boundary marker(s) Format: ( nacelle inflow marker, fan face Mach, ... ) \ingroup Config*/ addStringDoubleListOption("MARKER_ENGINE_INFLOW", nMarker_EngineInflow, Marker_EngineInflow, EngineInflow_Target); @@ -4614,11 +4609,16 @@ void CConfig::SetPostprocessing(unsigned short val_software, unsigned short val_ /*--- Check that if the wall roughness array are compatible and set deafult values if needed. ---*/ if ((nMarker_HeatFlux > 0) || (nMarker_Isothermal > 0) || (nMarker_CHTInterface > 0)) { - /*--- If no roughness is specified all walls are assumed to be smooth. ---*/ + + /*--- The total number of wall markers. ---*/ + unsigned short nWall = nMarker_HeatFlux + nMarker_Isothermal + nMarker_CHTInterface; + + /*--- If no roughness is specified all walls are assumed to be smooth. ---*/ if (nRough_Wall == 0) { - nRough_Wall = nMarker_HeatFlux + nMarker_Isothermal + nMarker_CHTInterface; - Roughness_Height = new su2double [nRough_Wall]; - Kind_Wall = new unsigned short [nRough_Wall]; + + nRough_Wall = nWall; + Roughness_Height = new su2double [nWall]; + Kind_Wall = new unsigned short [nWall]; for (unsigned short iMarker = 0; iMarker < nMarker_HeatFlux; iMarker++) { Roughness_Height[iMarker] = 0.0; Kind_Wall[iMarker] = SMOOTH; @@ -4631,7 +4631,12 @@ void CConfig::SetPostprocessing(unsigned short val_software, unsigned short val_ Roughness_Height[nMarker_HeatFlux + nMarker_Isothermal + iMarker] = 0.0; Kind_Wall[nMarker_HeatFlux + nMarker_Isothermal + iMarker] = SMOOTH; } - } else { /*--- Check name of the marker and assign the corresponding roughness. ---*/ + + /*--- Check name of the marker and assign the corresponding roughness. ---*/ + } else if (nRough_Wall > nWall) { + SU2_MPI::Error("Mismatch in number of rough walls and solid walls. Number of rough walls cannot be more than solid walls.", CURRENT_FUNCTION); + } else { + /*--- Store roughness heights in a temp array. ---*/ vector temp_rough; for (iMarker = 0; iMarker < nRough_Wall; iMarker++) @@ -4640,10 +4645,10 @@ void CConfig::SetPostprocessing(unsigned short val_software, unsigned short val_ /*--- Reallocate the roughness arrays in case not all walls are rough. ---*/ delete Roughness_Height; delete Kind_Wall; - unsigned short nWall = nMarker_HeatFlux + nMarker_Isothermal + nMarker_CHTInterface; - unsigned short jMarker; Roughness_Height = new su2double [nWall]; Kind_Wall = new unsigned short [nWall]; + unsigned short jMarker, chkRough = 0; + /*--- Initialize everything to smooth. ---*/ for (iMarker = 0; iMarker < nWall; iMarker++) { Roughness_Height[iMarker] = 0.0; @@ -4653,22 +4658,31 @@ void CConfig::SetPostprocessing(unsigned short val_software, unsigned short val_ /*--- Look through heat flux, isothermal and cht_interface markers and assign proper values. ---*/ for (iMarker = 0; iMarker < nRough_Wall; iMarker++) { for (jMarker = 0; jMarker < nMarker_HeatFlux; jMarker++) - if (Marker_HeatFlux[jMarker].compare(Marker_RoughWall[iMarker]) == 0) + if (Marker_HeatFlux[jMarker].compare(Marker_RoughWall[iMarker]) == 0) { Roughness_Height[jMarker] = temp_rough[iMarker]; + chkRough++; + } for (jMarker = 0; jMarker < nMarker_Isothermal; jMarker++) - if (Marker_Isothermal[jMarker].compare(Marker_RoughWall[iMarker]) == 0) + if (Marker_Isothermal[jMarker].compare(Marker_RoughWall[iMarker]) == 0) { Roughness_Height[nMarker_HeatFlux + jMarker] = temp_rough[iMarker]; + chkRough++; + } for (jMarker = 0; jMarker < nMarker_CHTInterface; jMarker++) - if (Marker_CHTInterface[jMarker].compare(Marker_RoughWall[iMarker]) == 0) + if (Marker_CHTInterface[jMarker].compare(Marker_RoughWall[iMarker]) == 0) { Roughness_Height[nMarker_HeatFlux + nMarker_Isothermal + jMarker] = temp_rough[iMarker]; + chkRough++; + } } - /*--- Update kind_wall when a non zero roughness value is specified. ---*/ + /*--- Update kind_wall when a non zero roughness value is specified. Also check if a non solid wall marker was specified as rough. ---*/ for (iMarker = 0; iMarker < nWall; iMarker++) if (Roughness_Height[iMarker] != 0.0) Kind_Wall[iMarker] = ROUGH; + + if (chkRough != nRough_Wall) + SU2_MPI::Error("Only solid walls can be rough.", CURRENT_FUNCTION); } } @@ -7622,8 +7636,6 @@ CConfig::~CConfig(void) { delete[] FlowLoad_Value; delete[] Roughness_Height; delete[] Wall_Emissivity; - delete[] GlobalMarkerStorageDispl; - delete[] GlobalRoughness_Height; /*--- related to periodic boundary conditions ---*/ for (iMarker = 0; iMarker < nMarker_PerBound; iMarker++) { @@ -9719,31 +9731,3 @@ void CConfig::SetMultizone(CConfig *driver_config, CConfig **config_container){ Multizone_Residual = true; } - -su2double CConfig::GetLocalRoughness(int rankID, unsigned short markerID) const { - - su2double localRoughness = 0.0; - int index; - - index = GlobalMarkerStorageDispl[rankID]+ markerID; - localRoughness = GlobalRoughness_Height[index]; - -return localRoughness; -} - -void CConfig::SetGlobalMarkerArray(int* global_displ, int size) { - - if (GlobalMarkerStorageDispl == nullptr) GlobalMarkerStorageDispl = new int [size]; - GlobalMarkerStorageDispl[0] = 0; - for (int iRank = 0; iRank < size; iRank++) - GlobalMarkerStorageDispl[iRank] = global_displ[iRank]; - -} - -void CConfig::SetGlobalRoughnessArray(su2double *global_rough, int sizeGlobal) { - - if (GlobalRoughness_Height == nullptr) GlobalRoughness_Height = new su2double [sizeGlobal]; - for (int iMarker = 0; iMarker < sizeGlobal; iMarker++) - GlobalRoughness_Height[iMarker] = global_rough[iMarker]; - -} diff --git a/Common/src/geometry/CPhysicalGeometry.cpp b/Common/src/geometry/CPhysicalGeometry.cpp index 4ef1ae71194f..dcbf6065d229 100644 --- a/Common/src/geometry/CPhysicalGeometry.cpp +++ b/Common/src/geometry/CPhysicalGeometry.cpp @@ -312,6 +312,9 @@ CPhysicalGeometry::CPhysicalGeometry(CGeometry *geometry, delete [] Elem_ID_BoundTria_Linear; delete [] Elem_ID_BoundQuad_Linear; + delete[] GlobalMarkerStorageDispl; + delete[] GlobalRoughness_Height; + } CPhysicalGeometry::~CPhysicalGeometry(void) { @@ -11365,6 +11368,8 @@ void CPhysicalGeometry::SetWallDistance(const CConfig *config, CADTElemClass *Wa /*--- Step 3: Loop over all interior mesh nodes and compute minimum ---*/ /*--- distance to a solid wall element ---*/ /*--------------------------------------------------------------------------*/ + int rank; + SU2_MPI::Comm_rank(MPI_COMM_WORLD, &rank); SU2_OMP_PARALLEL if (!WallADT->IsEmpty()) { @@ -11378,11 +11383,13 @@ void CPhysicalGeometry::SetWallDistance(const CConfig *config, CADTElemClass *Wa int rankID; su2double dist; su2double localRoughness; + int index; WallADT->DetermineNearestElement(nodes->GetCoord(iPoint), dist, markerID, elemID, rankID); nodes->SetWall_Distance(iPoint, min(dist,nodes->GetWall_Distance(iPoint))); - localRoughness = config->GetLocalRoughness(rankID, markerID); + index = GlobalMarkerStorageDispl[rankID]+ markerID; + localRoughness = GlobalRoughness_Height[index]; nodes->SetRoughnessHeight(iPoint, localRoughness); } @@ -11399,7 +11406,8 @@ unsigned short iMarker; int size = 1; unsigned short nMarker_All = config->GetnMarker_All(); local_displ[0] = 0; - config->SetGlobalMarkerArray(local_displ,size); + if (GlobalMarkerStorageDispl == nullptr) GlobalMarkerStorageDispl = new int [size]; + GlobalMarkerStorageDispl[0] = 0; su2double *localRough = new su2double [nMarker_All]; @@ -11410,7 +11418,9 @@ unsigned short iMarker; for (iMarker = 0; iMarker < nMarker_All; iMarker++) localRough[iMarker] = 0.0 ; - config->SetGlobalRoughnessArray(localRough, nMarker_All); + if (GlobalRoughness_Height == nullptr) GlobalRoughness_Height = new su2double [nMarker_All]; + for (int iMarker = 0; iMarker < nMarker_All; iMarker++) + GlobalRoughness_Height[iMarker] = localRough[iMarker]; #else int rank, size; @@ -11419,7 +11429,7 @@ unsigned short iMarker; SU2_MPI::Comm_size(MPI_COMM_WORLD, &size); int *displs = new int [size]; - int* recvCounts = new int [size];//, displs(size); + int* recvCounts = new int [size]; int sizeLocal = (int) nMarker_All; // number of local markers //Communicate size of local marker array and make an array large enough to hold all data @@ -11431,28 +11441,35 @@ unsigned short iMarker; for(int i=1; iSetGlobalMarkerArray(displs, size); + if (GlobalMarkerStorageDispl == nullptr) GlobalMarkerStorageDispl = new int [size]; + GlobalMarkerStorageDispl[0] = 0; + for (int iRank = 1; iRank < size; iRank++) + GlobalMarkerStorageDispl[iRank] = displs[iRank]; //total size int sizeGlobal = displs[size-1] + recvCounts[size-1]; - /*--- Allocate local and global arrays to hold roughness. ---*/ + /*--- Allocate local and global arrays to hold roughness. ---*/ su2double *localRough = new su2double [nMarker_All]; // local number of markers su2double *globalRough = new su2double[sizeGlobal]; // all markers including send recieve - if (config->GetnRoughWall() > 0) + if (config->GetnRoughWall() > 0) { + for (iMarker = 0; iMarker < nMarker_All; iMarker++) + localRough[iMarker] = config->GetWall_RoughnessHeight(config->GetMarker_All_TagBound(iMarker)); + } + else { for (iMarker = 0; iMarker < nMarker_All; iMarker++) - localRough[iMarker] = config->GetWall_RoughnessHeight(config->GetMarker_All_TagBound(iMarker)); - else - for (iMarker = 0; iMarker < nMarker_All; iMarker++) localRough[iMarker] = 0.0 ; + } SU2_MPI::Allgatherv( localRough, sizeLocal, MPI_DOUBLE, globalRough , recvCounts, displs, MPI_DOUBLE, MPI_COMM_WORLD); /*--- Set the global array of roughness per marker. ---*/ - config->SetGlobalRoughnessArray(globalRough, sizeGlobal); + if (GlobalRoughness_Height == nullptr) GlobalRoughness_Height = new su2double [sizeGlobal]; + for (int iMarker = 0; iMarker < sizeGlobal; iMarker++) + GlobalRoughness_Height[iMarker] = globalRough[iMarker]; #endif } diff --git a/externals/meson b/externals/meson index 0435691e83fb..29ef4478df6d 160000 --- a/externals/meson +++ b/externals/meson @@ -1 +1 @@ -Subproject commit 0435691e83fb7172e2a9635d2eb32d5521089916 +Subproject commit 29ef4478df6d3aaca40c7993f125b29409be1de2 diff --git a/externals/ninja b/externals/ninja index 2d15b04e4112..52649de2c56b 160000 --- a/externals/ninja +++ b/externals/ninja @@ -1 +1 @@ -Subproject commit 2d15b04e411229cb902332957281622119025e77 +Subproject commit 52649de2c56b63f42bc59513d51286531c595b44 From c60cb5f511df42cce75ad5542801da9bd42351a7 Mon Sep 17 00:00:00 2001 From: koodlyakshay Date: Tue, 7 Jul 2020 21:18:45 +0200 Subject: [PATCH 21/29] Clean up and remove unused code. --- Common/include/CConfig.hpp | 7 -- Common/include/geometry/CGeometry.hpp | 2 +- Common/include/geometry/CPhysicalGeometry.hpp | 16 +--- Common/src/CConfig.cpp | 84 ++++++++----------- Common/src/geometry/CPhysicalGeometry.cpp | 32 +++---- SU2_CFD/src/drivers/CDriver.cpp | 4 +- .../src/numerics/turbulent/turb_sources.cpp | 7 -- SU2_CFD/src/solvers/CTurbSASolver.cpp | 6 +- 8 files changed, 55 insertions(+), 103 deletions(-) diff --git a/Common/include/CConfig.hpp b/Common/include/CConfig.hpp index 467e6822d4ca..29c6b75f61e3 100644 --- a/Common/include/CConfig.hpp +++ b/Common/include/CConfig.hpp @@ -6790,13 +6790,6 @@ class CConfig { */ su2double GetWall_RoughnessHeight(string val_marker) const; - /*! - * \brief Set the wall roughness height on a wall boundary (Heatflux or Isothermal). - * \param[in] val_index - Index corresponding to the boundary. - * \return The wall roughness height. - */ - void SetWall_RoughnessHeight(string val_marker, su2double val_roughness); - /*! * \brief Get the target (pressure, massflow, etc) at an engine inflow boundary. * \param[in] val_index - Index corresponding to the engine inflow boundary. diff --git a/Common/include/geometry/CGeometry.hpp b/Common/include/geometry/CGeometry.hpp index 673eff1fea4e..85a0f8e3e011 100644 --- a/Common/include/geometry/CGeometry.hpp +++ b/Common/include/geometry/CGeometry.hpp @@ -1716,7 +1716,7 @@ class CGeometry { /*! * \brief Set roughness values for markers in a global array. */ - inline virtual void SetGlobalMarkerRoughness(CConfig *config) { } + inline virtual void SetGlobalMarkerRoughness(const CConfig* config) { } }; diff --git a/Common/include/geometry/CPhysicalGeometry.hpp b/Common/include/geometry/CPhysicalGeometry.hpp index b5b8252776c1..a90d6dc3be68 100644 --- a/Common/include/geometry/CPhysicalGeometry.hpp +++ b/Common/include/geometry/CPhysicalGeometry.hpp @@ -810,20 +810,6 @@ class CPhysicalGeometry final : public CGeometry { /*! * \brief Set roughness values for markers in a global array. */ - void SetGlobalMarkerRoughness(CConfig *config) override; - - /*! - * \brief Set the displacements for each rank in global array. - * \param [in] global displacemts in order of rank. - * \param [in] total number of ranks - */ - void SetGlobalMarkerArray(int *global_displ, int size) ; - - /*! - * \brief Set the roughness for each marker in global array. - * \param [in] global list of roughness for every marker. - * \param [in] total number of marker including send recieves. - */ - void SetGlobalRoughnessArray(su2double *global_rough, int sizeGlobal) ; + void SetGlobalMarkerRoughness(const CConfig* config) override; }; diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index bf9c624279ed..a86a063bef50 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -4607,13 +4607,13 @@ void CConfig::SetPostprocessing(unsigned short val_software, unsigned short val_ } } - /*--- Check that if the wall roughness array are compatible and set deafult values if needed. ---*/ + /*--- Check that if the wall roughness array are compatible and set deafult values if needed. ---*/ if ((nMarker_HeatFlux > 0) || (nMarker_Isothermal > 0) || (nMarker_CHTInterface > 0)) { - /*--- The total number of wall markers. ---*/ + /*--- The total number of wall markers. ---*/ unsigned short nWall = nMarker_HeatFlux + nMarker_Isothermal + nMarker_CHTInterface; - /*--- If no roughness is specified all walls are assumed to be smooth. ---*/ + /*--- If no roughness is specified all walls are assumed to be smooth. ---*/ if (nRough_Wall == 0) { nRough_Wall = nWall; @@ -4632,15 +4632,15 @@ void CConfig::SetPostprocessing(unsigned short val_software, unsigned short val_ Kind_Wall[nMarker_HeatFlux + nMarker_Isothermal + iMarker] = SMOOTH; } - /*--- Check name of the marker and assign the corresponding roughness. ---*/ + /*--- Check for mismatch in number of rough walls and solid walls. ---*/ } else if (nRough_Wall > nWall) { SU2_MPI::Error("Mismatch in number of rough walls and solid walls. Number of rough walls cannot be more than solid walls.", CURRENT_FUNCTION); - } else { - + /*--- Check name of the marker and assign the corresponding roughness. ---*/ + } else { /*--- Store roughness heights in a temp array. ---*/ vector temp_rough; - for (iMarker = 0; iMarker < nRough_Wall; iMarker++) - temp_rough.push_back(Roughness_Height[iMarker]); + for (iMarker = 0; iMarker < nRough_Wall; iMarker++) + temp_rough.push_back(Roughness_Height[iMarker]); /*--- Reallocate the roughness arrays in case not all walls are rough. ---*/ delete Roughness_Height; @@ -4657,13 +4657,13 @@ void CConfig::SetPostprocessing(unsigned short val_software, unsigned short val_ /*--- Look through heat flux, isothermal and cht_interface markers and assign proper values. ---*/ for (iMarker = 0; iMarker < nRough_Wall; iMarker++) { - for (jMarker = 0; jMarker < nMarker_HeatFlux; jMarker++) + for (jMarker = 0; jMarker < nMarker_HeatFlux; jMarker++) if (Marker_HeatFlux[jMarker].compare(Marker_RoughWall[iMarker]) == 0) { Roughness_Height[jMarker] = temp_rough[iMarker]; chkRough++; } - for (jMarker = 0; jMarker < nMarker_Isothermal; jMarker++) + for (jMarker = 0; jMarker < nMarker_Isothermal; jMarker++) if (Marker_Isothermal[jMarker].compare(Marker_RoughWall[iMarker]) == 0) { Roughness_Height[nMarker_HeatFlux + jMarker] = temp_rough[iMarker]; chkRough++; @@ -4676,11 +4676,12 @@ void CConfig::SetPostprocessing(unsigned short val_software, unsigned short val_ } } - /*--- Update kind_wall when a non zero roughness value is specified. Also check if a non solid wall marker was specified as rough. ---*/ + /*--- Update kind_wall when a non zero roughness value is specified. ---*/ for (iMarker = 0; iMarker < nWall; iMarker++) if (Roughness_Height[iMarker] != 0.0) Kind_Wall[iMarker] = ROUGH; + /*--- Check if a non solid wall marker was specified as rough. ---*/ if (chkRough != nRough_Wall) SU2_MPI::Error("Only solid walls can be rough.", CURRENT_FUNCTION); } @@ -8749,62 +8750,51 @@ su2double CConfig::GetWall_HeatFlux(string val_marker) const { unsigned short CConfig::GetKindWall(string val_marker) const { unsigned short iMarker; short flag = -1; + for (iMarker = 0; iMarker < nMarker_HeatFlux; iMarker++) if (Marker_HeatFlux[iMarker] == val_marker) { flag = iMarker; - break; + return Kind_Wall[flag]; } - if (flag == -1) - for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) - if (Marker_Isothermal[iMarker] == val_marker) { + + for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) + if (Marker_Isothermal[iMarker] == val_marker) { flag = nMarker_HeatFlux + iMarker; - break; + return Kind_Wall[flag]; } - if (flag == -1) { - for (iMarker = 0; iMarker < nMarker_CHTInterface; iMarker++) - if (Marker_CHTInterface[iMarker] == val_marker) { - flag = nMarker_HeatFlux + nMarker_Isothermal + iMarker; - break; - } + + for (iMarker = 0; iMarker < nMarker_CHTInterface; iMarker++) + if (Marker_CHTInterface[iMarker] == val_marker) { + flag = nMarker_HeatFlux + nMarker_Isothermal + iMarker; + return Kind_Wall[flag]; } - if (flag == -1) return 1; - else return Kind_Wall[flag]; -} -void CConfig::SetWall_RoughnessHeight(string val_marker, su2double val_roughness) { - //Roughness_Height[val_marker] = val_roughness; + return SMOOTH; } su2double CConfig::GetWall_RoughnessHeight(string val_marker) const { unsigned short iMarker = 0; short flag = -1; - if (nMarker_HeatFlux > 0 || nMarker_Isothermal > 0) { + if (nMarker_HeatFlux > 0 || nMarker_Isothermal > 0 || nMarker_CHTInterface > 0) { for (iMarker = 0; iMarker < nMarker_HeatFlux; iMarker++) if (Marker_HeatFlux[iMarker] == val_marker) { flag = iMarker; - break; + return Roughness_Height[flag]; + } + for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) + if (Marker_Isothermal[iMarker] == val_marker) { + flag = nMarker_HeatFlux + iMarker; + return Roughness_Height[flag]; + } + for (iMarker = 0; iMarker < nMarker_CHTInterface; iMarker++) + if (Marker_CHTInterface[iMarker] == val_marker) { + flag = nMarker_HeatFlux + nMarker_Isothermal + iMarker; + return Roughness_Height[flag]; } - if (flag == -1) { - for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) - if (Marker_Isothermal[iMarker] == val_marker) { - flag = nMarker_HeatFlux + iMarker; - break; - } - } - if (flag == -1) { - for (iMarker = 0; iMarker < nMarker_CHTInterface; iMarker++) - if (Marker_CHTInterface[iMarker] == val_marker) { - flag = nMarker_HeatFlux + nMarker_Isothermal + iMarker; - break; - } - } } - if (flag == -1) - return 0.0; //For cht_wall_interface - else - return Roughness_Height[flag]; + return 0.0; } unsigned short CConfig::GetWallFunction_Treatment(string val_marker) const { diff --git a/Common/src/geometry/CPhysicalGeometry.cpp b/Common/src/geometry/CPhysicalGeometry.cpp index dcbf6065d229..d0a919d817e5 100644 --- a/Common/src/geometry/CPhysicalGeometry.cpp +++ b/Common/src/geometry/CPhysicalGeometry.cpp @@ -11382,51 +11382,41 @@ void CPhysicalGeometry::SetWallDistance(const CConfig *config, CADTElemClass *Wa unsigned long elemID; int rankID; su2double dist; - su2double localRoughness; - int index; WallADT->DetermineNearestElement(nodes->GetCoord(iPoint), dist, markerID, elemID, rankID); nodes->SetWall_Distance(iPoint, min(dist,nodes->GetWall_Distance(iPoint))); - index = GlobalMarkerStorageDispl[rankID]+ markerID; - localRoughness = GlobalRoughness_Height[index]; - nodes->SetRoughnessHeight(iPoint, localRoughness); + if (config->GetnRoughWall() > 0) { + auto index = GlobalMarkerStorageDispl[rankID]+ markerID; + auto localRoughness = GlobalRoughness_Height[index]; + nodes->SetRoughnessHeight(iPoint, localRoughness); + } } } // end SU2_OMP_PARALLEL } -void CPhysicalGeometry::SetGlobalMarkerRoughness(CConfig *config) { +void CPhysicalGeometry::SetGlobalMarkerRoughness(const CConfig* config) { unsigned short iMarker; #ifndef HAVE_MPI - int *local_displ = new int [1]; int size = 1; unsigned short nMarker_All = config->GetnMarker_All(); - local_displ[0] = 0; if (GlobalMarkerStorageDispl == nullptr) GlobalMarkerStorageDispl = new int [size]; + if (GlobalRoughness_Height == nullptr) GlobalRoughness_Height = new su2double [nMarker_All]; GlobalMarkerStorageDispl[0] = 0; - - su2double *localRough = new su2double [nMarker_All]; if (config->GetnRoughWall() > 0) for (iMarker = 0; iMarker < nMarker_All; iMarker++) - localRough[iMarker] = config->GetWall_RoughnessHeight(config->GetMarker_All_TagBound(iMarker)); - else - for (iMarker = 0; iMarker < nMarker_All; iMarker++) - localRough[iMarker] = 0.0 ; - - if (GlobalRoughness_Height == nullptr) GlobalRoughness_Height = new su2double [nMarker_All]; - for (int iMarker = 0; iMarker < nMarker_All; iMarker++) - GlobalRoughness_Height[iMarker] = localRough[iMarker]; + GlobalRoughness_Height[iMarker] = config->GetWall_RoughnessHeight(config->GetMarker_All_TagBound(iMarker)); + else + for (iMarker = 0; iMarker < nMarker_All; iMarker++) + GlobalRoughness_Height[iMarker] = 0.0 ; #else - int rank, size; unsigned short nMarker_All = config->GetnMarker_All(); - SU2_MPI::Comm_rank(MPI_COMM_WORLD, &rank); - SU2_MPI::Comm_size(MPI_COMM_WORLD, &size); int *displs = new int [size]; int* recvCounts = new int [size]; diff --git a/SU2_CFD/src/drivers/CDriver.cpp b/SU2_CFD/src/drivers/CDriver.cpp index 3ccbf72d6dd1..68be969f419f 100644 --- a/SU2_CFD/src/drivers/CDriver.cpp +++ b/SU2_CFD/src/drivers/CDriver.cpp @@ -859,9 +859,9 @@ void CDriver::Geometrical_Preprocessing_FVM(CConfig *config, CGeometry **&geomet geometry[MESH_0]->ComputeMeshQualityStatistics(config); } - /*--- Store marker list and roughness in a global array ---*/ + /*--- Store marker list and roughness in a global array. ---*/ - geometry[MESH_0]->SetGlobalMarkerRoughness(config); + if (config->GetnRoughWall() > 0) geometry[MESH_0]->SetGlobalMarkerRoughness(config); geometry[MESH_0]->SetMGLevel(MESH_0); if ((config->GetnMGLevels() != 0) && (rank == MASTER_NODE)) diff --git a/SU2_CFD/src/numerics/turbulent/turb_sources.cpp b/SU2_CFD/src/numerics/turbulent/turb_sources.cpp index 074f16ce6bdf..671793e3b726 100644 --- a/SU2_CFD/src/numerics/turbulent/turb_sources.cpp +++ b/SU2_CFD/src/numerics/turbulent/turb_sources.cpp @@ -125,13 +125,6 @@ CNumerics::ResidualType<> CSourcePieceWise_TurbSA::ComputeResidual(const CConfig dist_i_2 = dist_i*dist_i; nu = Laminar_Viscosity_i/Density_i; - /*--- Old values without roughness ---*/ - /*Ji = TurbVar_i[0]/nu; - Ji_2 = Ji*Ji; - Ji_3 = Ji_2*Ji; - fv1 = Ji_3/(Ji_3+cv1_3); - fv2 = 1.0 - Ji/(1.0+Ji*fv1);*/ - /*--- Modified values for roughness ---*/ /*--- Ref: Aupoix, B. and Spalart, P. R., "Extensions of the Spalart-Allmaras Turbulence Model to Account for Wall Roughness," * International Journal of Heat and Fluid Flow, Vol. 24, 2003, pp. 454-462. ---*/ diff --git a/SU2_CFD/src/solvers/CTurbSASolver.cpp b/SU2_CFD/src/solvers/CTurbSASolver.cpp index d7343a8eb252..5575f6d3e1c3 100644 --- a/SU2_CFD/src/solvers/CTurbSASolver.cpp +++ b/SU2_CFD/src/solvers/CTurbSASolver.cpp @@ -300,7 +300,7 @@ void CTurbSASolver::Preprocessing(CGeometry *geometry, CSolver **solver_containe void CTurbSASolver::Postprocessing(CGeometry *geometry, CSolver **solver_container, CConfig *config, unsigned short iMesh) { - const su2double cv1_3 = 7.1*7.1*7.1, cR1 = 0.5; + const su2double cv1_3 = 7.1*7.1*7.1, cR1 = 0.5, rough_const = 0.03; const bool neg_spalart_allmaras = (config->GetKind_Turb_Model() == SA_NEG); @@ -317,10 +317,10 @@ void CTurbSASolver::Postprocessing(CGeometry *geometry, CSolver **solver_contain su2double roughness = geometry->nodes->GetRoughnessHeight(iPoint); su2double dist = geometry->nodes->GetWall_Distance(iPoint); - dist += 0.03*roughness; + dist += rough_const*roughness; su2double Ji = nu_hat/nu ; - if (roughness > 1.0e-10) + if (roughness > EPS) Ji+= cR1*roughness/(dist+EPS); su2double Ji_3 = Ji*Ji*Ji; From eb5c9fe0137f854a57d578ab21fb1cc5ad69178a Mon Sep 17 00:00:00 2001 From: koodlyakshay Date: Tue, 7 Jul 2020 23:24:36 +0200 Subject: [PATCH 22/29] Simplify some routines and clean up. --- Common/include/CConfig.hpp | 14 ++----- Common/include/geometry/CGeometry.hpp | 5 --- Common/include/geometry/CPhysicalGeometry.hpp | 2 +- Common/src/CConfig.cpp | 40 +++++-------------- Common/src/geometry/CPhysicalGeometry.cpp | 26 ++++++------ SU2_CFD/src/drivers/CDriver.cpp | 4 -- SU2_CFD/src/solvers/CIncNSSolver.cpp | 3 +- SU2_CFD/src/solvers/CNSSolver.cpp | 3 +- SU2_CFD/src/solvers/CTurbSASolver.cpp | 12 +++--- SU2_CFD/src/solvers/CTurbSSTSolver.cpp | 6 ++- 10 files changed, 39 insertions(+), 76 deletions(-) diff --git a/Common/include/CConfig.hpp b/Common/include/CConfig.hpp index 29c6b75f61e3..86e489b7662e 100644 --- a/Common/include/CConfig.hpp +++ b/Common/include/CConfig.hpp @@ -5016,12 +5016,6 @@ class CConfig { */ unsigned short GetKind_ActDisk(void) const { return Kind_ActDisk; } - /*! - * \brief Get the kind of wall. - * \return Kind of wall - smooth or rough. - */ - unsigned short GetKindWall(string val_marker) const; - /*! * \brief Set the kind of wall - rough or smooth. */ @@ -6784,11 +6778,11 @@ class CConfig { su2double* GetWallFunction_DoubleInfo(string val_marker); /*! - * \brief Get the wall roughness height on a wall boundary (Heatflux or Isothermal). + * \brief Get the type of wall and roughness height on a wall boundary (Heatflux or Isothermal). * \param[in] val_index - Index corresponding to the boundary. - * \return The wall roughness height. - */ - su2double GetWall_RoughnessHeight(string val_marker) const; + * \return The wall type and roughness height. + */ + pair GetWallRoughnessProperties(string val_marker) const; /*! * \brief Get the target (pressure, massflow, etc) at an engine inflow boundary. diff --git a/Common/include/geometry/CGeometry.hpp b/Common/include/geometry/CGeometry.hpp index 85a0f8e3e011..e8962ae53de3 100644 --- a/Common/include/geometry/CGeometry.hpp +++ b/Common/include/geometry/CGeometry.hpp @@ -1712,11 +1712,6 @@ class CGeometry { * \param[in] geometry_container - Geometrical definition of the problem. */ static void ComputeWallDistance(const CConfig * const *config_container, CGeometry ****geometry_container); - - /*! - * \brief Set roughness values for markers in a global array. - */ - inline virtual void SetGlobalMarkerRoughness(const CConfig* config) { } }; diff --git a/Common/include/geometry/CPhysicalGeometry.hpp b/Common/include/geometry/CPhysicalGeometry.hpp index a90d6dc3be68..5d22a078f3db 100644 --- a/Common/include/geometry/CPhysicalGeometry.hpp +++ b/Common/include/geometry/CPhysicalGeometry.hpp @@ -810,6 +810,6 @@ class CPhysicalGeometry final : public CGeometry { /*! * \brief Set roughness values for markers in a global array. */ - void SetGlobalMarkerRoughness(const CConfig* config) override; + void SetGlobalMarkerRoughness(const CConfig* config); }; diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index ca6541b2292d..804ed25dae5c 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -8747,54 +8747,34 @@ su2double CConfig::GetWall_HeatFlux(string val_marker) const { return Heat_Flux[iMarker_HeatFlux]; } -unsigned short CConfig::GetKindWall(string val_marker) const { - unsigned short iMarker; - short flag = -1; - - for (iMarker = 0; iMarker < nMarker_HeatFlux; iMarker++) - if (Marker_HeatFlux[iMarker] == val_marker) { - flag = iMarker; - return Kind_Wall[flag]; - } - - for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) - if (Marker_Isothermal[iMarker] == val_marker) { - flag = nMarker_HeatFlux + iMarker; - return Kind_Wall[flag]; - } - - for (iMarker = 0; iMarker < nMarker_CHTInterface; iMarker++) - if (Marker_CHTInterface[iMarker] == val_marker) { - flag = nMarker_HeatFlux + nMarker_Isothermal + iMarker; - return Kind_Wall[flag]; - } - - return SMOOTH; -} - -su2double CConfig::GetWall_RoughnessHeight(string val_marker) const { +pair CConfig::GetWallRoughnessProperties(string val_marker) const { unsigned short iMarker = 0; short flag = -1; + pair WallProp; if (nMarker_HeatFlux > 0 || nMarker_Isothermal > 0 || nMarker_CHTInterface > 0) { for (iMarker = 0; iMarker < nMarker_HeatFlux; iMarker++) if (Marker_HeatFlux[iMarker] == val_marker) { flag = iMarker; - return Roughness_Height[flag]; + WallProp = make_pair(Kind_Wall[flag], Roughness_Height[flag]); + return WallProp; } for (iMarker = 0; iMarker < nMarker_Isothermal; iMarker++) if (Marker_Isothermal[iMarker] == val_marker) { flag = nMarker_HeatFlux + iMarker; - return Roughness_Height[flag]; + WallProp = make_pair(Kind_Wall[flag], Roughness_Height[flag]); + return WallProp; } for (iMarker = 0; iMarker < nMarker_CHTInterface; iMarker++) if (Marker_CHTInterface[iMarker] == val_marker) { flag = nMarker_HeatFlux + nMarker_Isothermal + iMarker; - return Roughness_Height[flag]; + WallProp = make_pair(Kind_Wall[flag], Roughness_Height[flag]); + return WallProp; } } - return 0.0; + WallProp = make_pair(SMOOTH, 0.0); + return WallProp; } unsigned short CConfig::GetWallFunction_Treatment(string val_marker) const { diff --git a/Common/src/geometry/CPhysicalGeometry.cpp b/Common/src/geometry/CPhysicalGeometry.cpp index d0a919d817e5..d0dbca621e1f 100644 --- a/Common/src/geometry/CPhysicalGeometry.cpp +++ b/Common/src/geometry/CPhysicalGeometry.cpp @@ -11371,6 +11371,9 @@ void CPhysicalGeometry::SetWallDistance(const CConfig *config, CADTElemClass *Wa int rank; SU2_MPI::Comm_rank(MPI_COMM_WORLD, &rank); + /*--- Store marker list and roughness in a global array. ---*/ + if (config->GetnRoughWall() > 0) SetGlobalMarkerRoughness(config); + SU2_OMP_PARALLEL if (!WallADT->IsEmpty()) { /*--- Solid wall boundary nodes are present. Compute the wall @@ -11408,12 +11411,10 @@ unsigned short iMarker; if (GlobalRoughness_Height == nullptr) GlobalRoughness_Height = new su2double [nMarker_All]; GlobalMarkerStorageDispl[0] = 0; - if (config->GetnRoughWall() > 0) - for (iMarker = 0; iMarker < nMarker_All; iMarker++) - GlobalRoughness_Height[iMarker] = config->GetWall_RoughnessHeight(config->GetMarker_All_TagBound(iMarker)); - else - for (iMarker = 0; iMarker < nMarker_All; iMarker++) - GlobalRoughness_Height[iMarker] = 0.0 ; + for (iMarker = 0; iMarker < nMarker_All; iMarker++) { + wallprop = config->GetWallRoughnessProperties(config->GetMarker_All_TagBound(iMarker)); + GlobalRoughness_Height[iMarker] = wallprop.second; + } #else unsigned short nMarker_All = config->GetnMarker_All(); @@ -11442,14 +11443,11 @@ unsigned short iMarker; /*--- Allocate local and global arrays to hold roughness. ---*/ su2double *localRough = new su2double [nMarker_All]; // local number of markers su2double *globalRough = new su2double[sizeGlobal]; // all markers including send recieve - - if (config->GetnRoughWall() > 0) { - for (iMarker = 0; iMarker < nMarker_All; iMarker++) - localRough[iMarker] = config->GetWall_RoughnessHeight(config->GetMarker_All_TagBound(iMarker)); - } - else { - for (iMarker = 0; iMarker < nMarker_All; iMarker++) - localRough[iMarker] = 0.0 ; + pair wallprop; + + for (iMarker = 0; iMarker < nMarker_All; iMarker++) { + wallprop = config->GetWallRoughnessProperties(config->GetMarker_All_TagBound(iMarker)); + localRough[iMarker] = wallprop.second; } SU2_MPI::Allgatherv( localRough, sizeLocal, MPI_DOUBLE, globalRough , diff --git a/SU2_CFD/src/drivers/CDriver.cpp b/SU2_CFD/src/drivers/CDriver.cpp index 68be969f419f..2dae63977ae5 100644 --- a/SU2_CFD/src/drivers/CDriver.cpp +++ b/SU2_CFD/src/drivers/CDriver.cpp @@ -859,10 +859,6 @@ void CDriver::Geometrical_Preprocessing_FVM(CConfig *config, CGeometry **&geomet geometry[MESH_0]->ComputeMeshQualityStatistics(config); } - /*--- Store marker list and roughness in a global array. ---*/ - - if (config->GetnRoughWall() > 0) geometry[MESH_0]->SetGlobalMarkerRoughness(config); - geometry[MESH_0]->SetMGLevel(MESH_0); if ((config->GetnMGLevels() != 0) && (rank == MASTER_NODE)) cout << "Setting the multigrid structure." << endl; diff --git a/SU2_CFD/src/solvers/CIncNSSolver.cpp b/SU2_CFD/src/solvers/CIncNSSolver.cpp index 09465bcee052..f5fd72caa86c 100644 --- a/SU2_CFD/src/solvers/CIncNSSolver.cpp +++ b/SU2_CFD/src/solvers/CIncNSSolver.cpp @@ -1249,8 +1249,7 @@ void CIncNSSolver::Friction_Forces(CGeometry *geometry, CConfig *config) { /*--- TurbViscosity is usually zero on the walls, except when the SA roughness model is used. * If the SA roughness model is not active, this term will be zero and will not affect the solution---*/ - TurbViscosity = 0.0; - if ((config->GetKindWall(Marker_Tag) == ROUGH ) && (config->GetKind_Turb_Model() == SA)) TurbViscosity = nodes->GetEddyViscosity(iPoint); + TurbViscosity = nodes->GetEddyViscosity(iPoint); Viscosity = nodes->GetLaminarViscosity(iPoint); Viscosity += TurbViscosity; diff --git a/SU2_CFD/src/solvers/CNSSolver.cpp b/SU2_CFD/src/solvers/CNSSolver.cpp index 822cedfdfee9..78b97039de4f 100644 --- a/SU2_CFD/src/solvers/CNSSolver.cpp +++ b/SU2_CFD/src/solvers/CNSSolver.cpp @@ -474,8 +474,7 @@ void CNSSolver::Friction_Forces(CGeometry *geometry, CConfig *config) { /*--- TurbViscosity is usually zero on the walls, except when the SA roughness model is used. * If the SA roughness model is not active, this term will be zero and will not affect the solution---*/ - TurbViscosity = 0.0; - if ((config->GetKindWall(Marker_Tag) == ROUGH ) && (config->GetKind_Turb_Model() == SA)) TurbViscosity = nodes->GetEddyViscosity(iPoint); + TurbViscosity = nodes->GetEddyViscosity(iPoint); Viscosity = nodes->GetLaminarViscosity(iPoint); Viscosity += TurbViscosity; Density = nodes->GetDensity(iPoint); diff --git a/SU2_CFD/src/solvers/CTurbSASolver.cpp b/SU2_CFD/src/solvers/CTurbSASolver.cpp index 5575f6d3e1c3..7e318dd87535 100644 --- a/SU2_CFD/src/solvers/CTurbSASolver.cpp +++ b/SU2_CFD/src/solvers/CTurbSASolver.cpp @@ -320,7 +320,7 @@ void CTurbSASolver::Postprocessing(CGeometry *geometry, CSolver **solver_contain dist += rough_const*roughness; su2double Ji = nu_hat/nu ; - if (roughness > EPS) + if (roughness > 1.0e-10) Ji+= cR1*roughness/(dist+EPS); su2double Ji_3 = Ji*Ji*Ji; @@ -464,12 +464,12 @@ void CTurbSASolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_conta SU2_OMP_BARRIER return; } - + bool rough_wall = false; - string Marker_Tag = config->GetMarker_All_TagBound(val_marker); - - su2double Roughness_Height = config->GetWall_RoughnessHeight(Marker_Tag); - if (config->GetKindWall(Marker_Tag) == ROUGH ) rough_wall = true; + string Marker_Tag = config->GetMarker_All_TagBound(val_marker); + unsigned short WallType; su2double Roughness_Height; + tie(WallType, Roughness_Height) = config->GetWallRoughnessProperties(Marker_Tag); + if (WallType == ROUGH ) rough_wall = true; /*--- The dirichlet condition is used only without wall function, otherwise the convergence is compromised as we are providing nu tilde values for the diff --git a/SU2_CFD/src/solvers/CTurbSSTSolver.cpp b/SU2_CFD/src/solvers/CTurbSSTSolver.cpp index e5f09895a7d2..15acd870b2f1 100644 --- a/SU2_CFD/src/solvers/CTurbSSTSolver.cpp +++ b/SU2_CFD/src/solvers/CTurbSSTSolver.cpp @@ -392,7 +392,9 @@ void CTurbSSTSolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_cont CNumerics *visc_numerics, CConfig *config, unsigned short val_marker) { bool rough_wall = false; string Marker_Tag = config->GetMarker_All_TagBound(val_marker); - if (config->GetKindWall(Marker_Tag) == ROUGH ) rough_wall = true; + unsigned short WallType; su2double Roughness_Height; + tie(WallType, Roughness_Height) = config->GetWallRoughnessProperties(Marker_Tag); + if (WallType == ROUGH ) rough_wall = true; SU2_OMP_FOR_STAT(OMP_MIN_SIZE) for (auto iVertex = 0u; iVertex < geometry->nVertex[val_marker]; iVertex++) { @@ -419,7 +421,7 @@ void CTurbSSTSolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver **solver_cont su2double FrictionVel = sqrt(fabs(WallShearStress)/density); /*--- Compute roughness in wall units. ---*/ - su2double Roughness_Height = config->GetWall_RoughnessHeight(Marker_Tag); + //su2double Roughness_Height = config->GetWall_RoughnessHeight(Marker_Tag); su2double kPlus = FrictionVel*Roughness_Height*density/laminar_viscosity; su2double S_R= 0.0; From c7dbe1990406d65043827e2dcc83322658a9498b Mon Sep 17 00:00:00 2001 From: koodlyakshay Date: Tue, 7 Jul 2020 23:30:25 +0200 Subject: [PATCH 23/29] Fix compile error. --- Common/src/geometry/CPhysicalGeometry.cpp | 1 + 1 file changed, 1 insertion(+) diff --git a/Common/src/geometry/CPhysicalGeometry.cpp b/Common/src/geometry/CPhysicalGeometry.cpp index d0dbca621e1f..e51ab62f23d0 100644 --- a/Common/src/geometry/CPhysicalGeometry.cpp +++ b/Common/src/geometry/CPhysicalGeometry.cpp @@ -11410,6 +11410,7 @@ unsigned short iMarker; if (GlobalMarkerStorageDispl == nullptr) GlobalMarkerStorageDispl = new int [size]; if (GlobalRoughness_Height == nullptr) GlobalRoughness_Height = new su2double [nMarker_All]; GlobalMarkerStorageDispl[0] = 0; + pair wallprop; for (iMarker = 0; iMarker < nMarker_All; iMarker++) { wallprop = config->GetWallRoughnessProperties(config->GetMarker_All_TagBound(iMarker)); From 594ba80076964d59498eff9bccdf4766988a5a6f Mon Sep 17 00:00:00 2001 From: koodlyakshay Date: Wed, 8 Jul 2020 00:24:24 +0200 Subject: [PATCH 24/29] Fix turb viscosity in output routines in flow solver. --- SU2_CFD/src/solvers/CIncNSSolver.cpp | 5 ++++- SU2_CFD/src/solvers/CNSSolver.cpp | 5 ++++- 2 files changed, 8 insertions(+), 2 deletions(-) diff --git a/SU2_CFD/src/solvers/CIncNSSolver.cpp b/SU2_CFD/src/solvers/CIncNSSolver.cpp index f5fd72caa86c..fb928eb329c8 100644 --- a/SU2_CFD/src/solvers/CIncNSSolver.cpp +++ b/SU2_CFD/src/solvers/CIncNSSolver.cpp @@ -1249,7 +1249,10 @@ void CIncNSSolver::Friction_Forces(CGeometry *geometry, CConfig *config) { /*--- TurbViscosity is usually zero on the walls, except when the SA roughness model is used. * If the SA roughness model is not active, this term will be zero and will not affect the solution---*/ - TurbViscosity = nodes->GetEddyViscosity(iPoint); + unsigned short WallType; su2double Roughness_Height; + tie(WallType, Roughness_Height) = config->GetWallRoughnessProperties(Marker_Tag); + TurbViscosity = 0.0; + if (WallType == ROUGH) TurbViscosity = nodes->GetEddyViscosity(iPoint); Viscosity = nodes->GetLaminarViscosity(iPoint); Viscosity += TurbViscosity; diff --git a/SU2_CFD/src/solvers/CNSSolver.cpp b/SU2_CFD/src/solvers/CNSSolver.cpp index 78b97039de4f..5aecd8500f39 100644 --- a/SU2_CFD/src/solvers/CNSSolver.cpp +++ b/SU2_CFD/src/solvers/CNSSolver.cpp @@ -474,7 +474,10 @@ void CNSSolver::Friction_Forces(CGeometry *geometry, CConfig *config) { /*--- TurbViscosity is usually zero on the walls, except when the SA roughness model is used. * If the SA roughness model is not active, this term will be zero and will not affect the solution---*/ - TurbViscosity = nodes->GetEddyViscosity(iPoint); + unsigned short WallType; su2double Roughness_Height; + tie(WallType, Roughness_Height) = config->GetWallRoughnessProperties(Marker_Tag); + TurbViscosity = 0.0; + if (WallType == ROUGH) TurbViscosity = nodes->GetEddyViscosity(iPoint); Viscosity = nodes->GetLaminarViscosity(iPoint); Viscosity += TurbViscosity; Density = nodes->GetDensity(iPoint); From 97c687dd0b2142bc2e8b15de798fc04e6a029f41 Mon Sep 17 00:00:00 2001 From: koodlyakshay Date: Wed, 8 Jul 2020 10:51:11 +0200 Subject: [PATCH 25/29] Small fixes to ease memory handling. --- Common/src/geometry/CPhysicalGeometry.cpp | 32 ++++++----------------- 1 file changed, 8 insertions(+), 24 deletions(-) diff --git a/Common/src/geometry/CPhysicalGeometry.cpp b/Common/src/geometry/CPhysicalGeometry.cpp index f182dbb3f9c5..9c5f3aac1008 100644 --- a/Common/src/geometry/CPhysicalGeometry.cpp +++ b/Common/src/geometry/CPhysicalGeometry.cpp @@ -11276,30 +11276,15 @@ void CPhysicalGeometry::SetWallDistance(const CConfig *config, CADTElemClass *Wa void CPhysicalGeometry::SetGlobalMarkerRoughness(const CConfig* config) { -unsigned short iMarker; -#ifndef HAVE_MPI - - int size = 1; - unsigned short nMarker_All = config->GetnMarker_All(); - if (GlobalMarkerStorageDispl == nullptr) GlobalMarkerStorageDispl = new int [size]; - if (GlobalRoughness_Height == nullptr) GlobalRoughness_Height = new su2double [nMarker_All]; - GlobalMarkerStorageDispl[0] = 0; - pair wallprop; - - for (iMarker = 0; iMarker < nMarker_All; iMarker++) { - wallprop = config->GetWallRoughnessProperties(config->GetMarker_All_TagBound(iMarker)); - GlobalRoughness_Height[iMarker] = wallprop.second; - } - -#else + unsigned short iMarker; unsigned short nMarker_All = config->GetnMarker_All(); - int *displs = new int [size]; - int* recvCounts = new int [size]; + vector displs(size); + vector recvCounts(size); int sizeLocal = (int) nMarker_All; // number of local markers //Communicate size of local marker array and make an array large enough to hold all data - SU2_MPI::Allgather(&sizeLocal, 1, MPI_INT, recvCounts, 1, + SU2_MPI::Allgather(&sizeLocal, 1, MPI_INT, recvCounts.data(), 1, MPI_INT, MPI_COMM_WORLD); // displacements based on size on each rank @@ -11316,8 +11301,8 @@ unsigned short iMarker; int sizeGlobal = displs[size-1] + recvCounts[size-1]; /*--- Allocate local and global arrays to hold roughness. ---*/ - su2double *localRough = new su2double [nMarker_All]; // local number of markers - su2double *globalRough = new su2double[sizeGlobal]; // all markers including send recieve + vector localRough(nMarker_All); // local number of markers + vector globalRough(sizeGlobal); // all markers including send recieve pair wallprop; for (iMarker = 0; iMarker < nMarker_All; iMarker++) { @@ -11325,8 +11310,8 @@ unsigned short iMarker; localRough[iMarker] = wallprop.second; } - SU2_MPI::Allgatherv( localRough, sizeLocal, MPI_DOUBLE, globalRough , - recvCounts, displs, MPI_DOUBLE, + SU2_MPI::Allgatherv( localRough.data(), sizeLocal, MPI_DOUBLE, globalRough.data() , + recvCounts.data(), displs.data(), MPI_DOUBLE, MPI_COMM_WORLD); /*--- Set the global array of roughness per marker. ---*/ @@ -11334,5 +11319,4 @@ unsigned short iMarker; for (int iMarker = 0; iMarker < sizeGlobal; iMarker++) GlobalRoughness_Height[iMarker] = globalRough[iMarker]; -#endif } From a46aff3cb2b6cd378a8845a06613a01685abc77d Mon Sep 17 00:00:00 2001 From: koodlyakshay Date: Wed, 8 Jul 2020 11:16:58 +0200 Subject: [PATCH 26/29] Change travis, config template and add delete statements for free variables. --- .travis.yml | 2 +- Common/src/geometry/CPhysicalGeometry.cpp | 19 ++++++++++++------- config_template.cfg | 12 +++--------- 3 files changed, 16 insertions(+), 17 deletions(-) diff --git a/.travis.yml b/.travis.yml index c93a3c04dec6..284818490731 100644 --- a/.travis.yml +++ b/.travis.yml @@ -18,7 +18,7 @@ compiler: notifications: email: recipients: - - koodlyakshay@gmail.com + - su2code-dev@lists.stanford.edu branches: only: diff --git a/Common/src/geometry/CPhysicalGeometry.cpp b/Common/src/geometry/CPhysicalGeometry.cpp index 9c5f3aac1008..3fbe52414dab 100644 --- a/Common/src/geometry/CPhysicalGeometry.cpp +++ b/Common/src/geometry/CPhysicalGeometry.cpp @@ -11279,12 +11279,12 @@ void CPhysicalGeometry::SetGlobalMarkerRoughness(const CConfig* config) { unsigned short iMarker; unsigned short nMarker_All = config->GetnMarker_All(); - vector displs(size); - vector recvCounts(size); + int *displs = new int [size]; + int* recvCounts = new int [size]; int sizeLocal = (int) nMarker_All; // number of local markers //Communicate size of local marker array and make an array large enough to hold all data - SU2_MPI::Allgather(&sizeLocal, 1, MPI_INT, recvCounts.data(), 1, + SU2_MPI::Allgather(&sizeLocal, 1, MPI_INT, recvCounts, 1, MPI_INT, MPI_COMM_WORLD); // displacements based on size on each rank @@ -11301,8 +11301,8 @@ void CPhysicalGeometry::SetGlobalMarkerRoughness(const CConfig* config) { int sizeGlobal = displs[size-1] + recvCounts[size-1]; /*--- Allocate local and global arrays to hold roughness. ---*/ - vector localRough(nMarker_All); // local number of markers - vector globalRough(sizeGlobal); // all markers including send recieve + su2double *localRough = new su2double [nMarker_All]; // local number of markers + su2double *globalRough = new su2double[sizeGlobal]; // all markers including send recieve pair wallprop; for (iMarker = 0; iMarker < nMarker_All; iMarker++) { @@ -11310,8 +11310,8 @@ void CPhysicalGeometry::SetGlobalMarkerRoughness(const CConfig* config) { localRough[iMarker] = wallprop.second; } - SU2_MPI::Allgatherv( localRough.data(), sizeLocal, MPI_DOUBLE, globalRough.data() , - recvCounts.data(), displs.data(), MPI_DOUBLE, + SU2_MPI::Allgatherv( localRough, sizeLocal, MPI_DOUBLE, globalRough , + recvCounts, displs, MPI_DOUBLE, MPI_COMM_WORLD); /*--- Set the global array of roughness per marker. ---*/ @@ -11319,4 +11319,9 @@ void CPhysicalGeometry::SetGlobalMarkerRoughness(const CConfig* config) { for (int iMarker = 0; iMarker < sizeGlobal; iMarker++) GlobalRoughness_Height[iMarker] = globalRough[iMarker]; + /*--- Deallocate local variables. ---*/ + delete [] displs; + delete [] recvCounts; + delete [] localRough; + delete [] globalRough; } diff --git a/config_template.cfg b/config_template.cfg index 708cd18fec07..19037982bfbd 100644 --- a/config_template.cfg +++ b/config_template.cfg @@ -809,16 +809,10 @@ GILES_EXTRA_RELAXFACTOR= ( 0.05, 0.05) % YES Non reflectivity activated, NO the Giles BC behaves as a normal 1D characteristic-based BC SPATIAL_FOURIER= NO % ------------------------ WALL ROUGHNESS DEFINITION --------------------------% -% Type of wall (SMOOTH or ROUGH). By default all walls are assumed to be smooth. -% If multiple walls are present, list type for each of them (ex- SMOOTH, ROUGH,..) -% in the order they are listed in MARKER_HEATFLUX followed by the MARKER_ISOTHERMAL -% definition. If all walls are smooth this option can be removed. -%WALL_TYPE= SMOOTH - % The equivalent sand grain roughness height (k_s) on each of the wall. This must be in m. -% This is a list of doubles each element corresponding to the MARKER defined in WALL_TYPE. -% For SMOOTH walls, set it to 0.0. -%WALL_ROUGHNESS = 0.0 +% This is a list of (string, double) each element corresponding to the MARKER defined in WALL_TYPE. +%WALL_ROUGHNESS = (wall1, ks1, wall2, ks2) +%WALL_ROUGHNESS = (wall1, ks1, wall2, 0.0) is also allowed % ------------------------ SURFACES IDENTIFICATION ----------------------------% % % Marker(s) of the surface in the surface flow solution file From 62133e65b7557b99ecb358c19178624634c2fffe Mon Sep 17 00:00:00 2001 From: Pedro Gomes Date: Wed, 8 Jul 2020 10:55:05 +0100 Subject: [PATCH 27/29] fix bug in signature of Allgatherv in base mpi wrapper --- Common/include/CConfig.hpp | 4 +- Common/include/geometry/CPhysicalGeometry.hpp | 5 +- Common/include/mpi_structure.hpp | 2 +- Common/include/mpi_structure.inl | 2 +- Common/src/CConfig.cpp | 10 +-- Common/src/geometry/CPhysicalGeometry.cpp | 71 +++++++------------ 6 files changed, 36 insertions(+), 58 deletions(-) diff --git a/Common/include/CConfig.hpp b/Common/include/CConfig.hpp index 86e489b7662e..f4ab25e0dfc7 100644 --- a/Common/include/CConfig.hpp +++ b/Common/include/CConfig.hpp @@ -5089,7 +5089,7 @@ class CConfig { * Gradients are w.r.t density, velocity[3], and pressure. when 2D gradient w.r.t. 3rd component of velocity set to 0. */ su2double GetCoeff_ObjChainRule(unsigned short iVar) const { return Obj_ChainRuleCoeff[iVar]; } - + /*! * \author H. Kline * \brief Get the flag indicating whether to comput a combined objective. @@ -6781,7 +6781,7 @@ class CConfig { * \brief Get the type of wall and roughness height on a wall boundary (Heatflux or Isothermal). * \param[in] val_index - Index corresponding to the boundary. * \return The wall type and roughness height. - */ + */ pair GetWallRoughnessProperties(string val_marker) const; /*! diff --git a/Common/include/geometry/CPhysicalGeometry.hpp b/Common/include/geometry/CPhysicalGeometry.hpp index 45e90e970ba8..ec80887f3fdd 100644 --- a/Common/include/geometry/CPhysicalGeometry.hpp +++ b/Common/include/geometry/CPhysicalGeometry.hpp @@ -105,8 +105,9 @@ class CPhysicalGeometry final : public CGeometry { unsigned long *Elem_ID_Line_Linear{nullptr}; unsigned long *Elem_ID_BoundTria_Linear{nullptr}; unsigned long *Elem_ID_BoundQuad_Linear{nullptr}; - int *GlobalMarkerStorageDispl{nullptr}; - su2double *GlobalRoughness_Height{nullptr}; + + vector GlobalMarkerStorageDispl; + vector GlobalRoughness_Height; public: /*--- This is to suppress Woverloaded-virtual, omitting it has no negative impact. ---*/ diff --git a/Common/include/mpi_structure.hpp b/Common/include/mpi_structure.hpp index 5fc6e5cf0315..f78622f6ffcd 100644 --- a/Common/include/mpi_structure.hpp +++ b/Common/include/mpi_structure.hpp @@ -417,7 +417,7 @@ class CBaseMPIWrapper { void *recvbuf, int recvcnt, Datatype recvtype, Comm comm); static void Allgatherv(void *sendbuf, int sendcnt, Datatype sendtype, - void *recvbuf, int recvcnt, int *displs, Datatype recvtype, Comm comm); + void *recvbuf, int *recvcnt, int *displs, Datatype recvtype, Comm comm); static void Sendrecv(void *sendbuf, int sendcnt, Datatype sendtype, int dest, int sendtag, void *recvbuf, int recvcnt, diff --git a/Common/include/mpi_structure.inl b/Common/include/mpi_structure.inl index 9e7ce2d08426..c6cbf27b1259 100644 --- a/Common/include/mpi_structure.inl +++ b/Common/include/mpi_structure.inl @@ -612,7 +612,7 @@ inline void CBaseMPIWrapper::Scatter(void *sendbuf, int sendcnt, Datatype sendty } inline void CBaseMPIWrapper::Allgatherv(void *sendbuf, int sendcnt, Datatype sendtype, - void *recvbuf, int recvcnt, int *displs, Datatype recvtype, Comm comm){ + void *recvbuf, int *recvcnt, int *displs, Datatype recvtype, Comm comm){ CopyData(sendbuf, recvbuf, sendcnt, sendtype); } diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index 804ed25dae5c..a1d8d3fa4994 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -865,7 +865,7 @@ void CConfig::SetPointersNull(void) { Inlet_FlowDir = nullptr; Inlet_Temperature = nullptr; Inlet_Pressure = nullptr; Inlet_Velocity = nullptr; Inflow_Mach = nullptr; Inflow_Pressure = nullptr; Exhaust_Pressure = nullptr; Outlet_Pressure = nullptr; Isothermal_Temperature= nullptr; - + ElasticityMod = nullptr; PoissonRatio = nullptr; MaterialDensity = nullptr; Load_Dir = nullptr; Load_Dir_Value = nullptr; Load_Dir_Multiplier = nullptr; @@ -4677,8 +4677,8 @@ void CConfig::SetPostprocessing(unsigned short val_software, unsigned short val_ } /*--- Update kind_wall when a non zero roughness value is specified. ---*/ - for (iMarker = 0; iMarker < nWall; iMarker++) - if (Roughness_Height[iMarker] != 0.0) + for (iMarker = 0; iMarker < nWall; iMarker++) + if (Roughness_Height[iMarker] != 0.0) Kind_Wall[iMarker] = ROUGH; /*--- Check if a non solid wall marker was specified as rough. ---*/ @@ -5334,7 +5334,7 @@ void CConfig::SetMarkers(unsigned short val_software) { Marker_CfgFile_KindBC[iMarker_CfgFile] = FLOWLOAD_BOUNDARY; iMarker_CfgFile++; } - + for (iMarker_CfgFile = 0; iMarker_CfgFile < nMarker_CfgFile; iMarker_CfgFile++) { Marker_CfgFile_Monitoring[iMarker_CfgFile] = NO; for (iMarker_Monitoring = 0; iMarker_Monitoring < nMarker_Monitoring; iMarker_Monitoring++) @@ -5460,7 +5460,7 @@ void CConfig::SetMarkers(unsigned short val_software) { if (Marker_CfgFile_TagBound[iMarker_CfgFile] == Marker_PyCustom[iMarker_PyCustom]) Marker_CfgFile_PyCustom[iMarker_CfgFile] = YES; } - + } void CConfig::SetOutput(unsigned short val_software, unsigned short val_izone) { diff --git a/Common/src/geometry/CPhysicalGeometry.cpp b/Common/src/geometry/CPhysicalGeometry.cpp index 3fbe52414dab..6324d4ef4a68 100644 --- a/Common/src/geometry/CPhysicalGeometry.cpp +++ b/Common/src/geometry/CPhysicalGeometry.cpp @@ -312,9 +312,6 @@ CPhysicalGeometry::CPhysicalGeometry(CGeometry *geometry, delete [] Elem_ID_BoundTria_Linear; delete [] Elem_ID_BoundQuad_Linear; - delete[] GlobalMarkerStorageDispl; - delete[] GlobalRoughness_Height; - } CPhysicalGeometry::~CPhysicalGeometry(void) { @@ -11242,8 +11239,6 @@ void CPhysicalGeometry::SetWallDistance(const CConfig *config, CADTElemClass *Wa /*--- Step 3: Loop over all interior mesh nodes and compute minimum ---*/ /*--- distance to a solid wall element ---*/ /*--------------------------------------------------------------------------*/ - int rank; - SU2_MPI::Comm_rank(MPI_COMM_WORLD, &rank); /*--- Store marker list and roughness in a global array. ---*/ if (config->GetnRoughWall() > 0) SetGlobalMarkerRoughness(config); @@ -11263,65 +11258,47 @@ void CPhysicalGeometry::SetWallDistance(const CConfig *config, CADTElemClass *Wa WallADT->DetermineNearestElement(nodes->GetCoord(iPoint), dist, markerID, elemID, rankID); nodes->SetWall_Distance(iPoint, min(dist,nodes->GetWall_Distance(iPoint))); + if (config->GetnRoughWall() > 0) { - auto index = GlobalMarkerStorageDispl[rankID]+ markerID; + auto index = GlobalMarkerStorageDispl[rankID] + markerID; auto localRoughness = GlobalRoughness_Height[index]; nodes->SetRoughnessHeight(iPoint, localRoughness); - } + } } - + } // end SU2_OMP_PARALLEL } void CPhysicalGeometry::SetGlobalMarkerRoughness(const CConfig* config) { - unsigned short iMarker; - unsigned short nMarker_All = config->GetnMarker_All(); - - int *displs = new int [size]; - int* recvCounts = new int [size]; - int sizeLocal = (int) nMarker_All; // number of local markers - - //Communicate size of local marker array and make an array large enough to hold all data - SU2_MPI::Allgather(&sizeLocal, 1, MPI_INT, recvCounts, 1, - MPI_INT, MPI_COMM_WORLD); - - // displacements based on size on each rank - displs[0] = 0; - for(int i=1; iGetnMarker_All(); + + vector recvCounts(size); + auto sizeLocal = static_cast(nMarker_All); // number of local markers + + /*--- Communicate size of local marker array and make an array large enough to hold all data. ---*/ + SU2_MPI::Allgather(&sizeLocal, 1, MPI_INT, recvCounts.data(), 1, MPI_INT, MPI_COMM_WORLD); /*--- Set the global array of displacements, needed to access the correct roughness element. ---*/ - if (GlobalMarkerStorageDispl == nullptr) GlobalMarkerStorageDispl = new int [size]; + GlobalMarkerStorageDispl.resize(size); GlobalMarkerStorageDispl[0] = 0; - for (int iRank = 1; iRank < size; iRank++) - GlobalMarkerStorageDispl[iRank] = displs[iRank]; - - //total size - int sizeGlobal = displs[size-1] + recvCounts[size-1]; + for (int iRank = 1; iRank < size; iRank++) + GlobalMarkerStorageDispl[iRank] = GlobalMarkerStorageDispl[iRank-1] + recvCounts[iRank-1]; + + /*--- Total size ---*/ + const auto sizeGlobal = GlobalMarkerStorageDispl[size-1] + recvCounts[size-1]; /*--- Allocate local and global arrays to hold roughness. ---*/ - su2double *localRough = new su2double [nMarker_All]; // local number of markers - su2double *globalRough = new su2double[sizeGlobal]; // all markers including send recieve - pair wallprop; + vector localRough(nMarker_All); // local number of markers + GlobalRoughness_Height.resize(sizeGlobal); // all markers including send recieve - for (iMarker = 0; iMarker < nMarker_All; iMarker++) { - wallprop = config->GetWallRoughnessProperties(config->GetMarker_All_TagBound(iMarker)); + for (auto iMarker = 0u; iMarker < nMarker_All; iMarker++) { + auto wallprop = config->GetWallRoughnessProperties(config->GetMarker_All_TagBound(iMarker)); localRough[iMarker] = wallprop.second; } - SU2_MPI::Allgatherv( localRough, sizeLocal, MPI_DOUBLE, globalRough , - recvCounts, displs, MPI_DOUBLE, - MPI_COMM_WORLD); - - /*--- Set the global array of roughness per marker. ---*/ - if (GlobalRoughness_Height == nullptr) GlobalRoughness_Height = new su2double [sizeGlobal]; - for (int iMarker = 0; iMarker < sizeGlobal; iMarker++) - GlobalRoughness_Height[iMarker] = globalRough[iMarker]; - - /*--- Deallocate local variables. ---*/ - delete [] displs; - delete [] recvCounts; - delete [] localRough; - delete [] globalRough; + /*--- Finally, gather the roughness of all markers. ---*/ + SU2_MPI::Allgatherv(localRough.data(), sizeLocal, MPI_DOUBLE, GlobalRoughness_Height.data(), + recvCounts.data(), GlobalMarkerStorageDispl.data(), MPI_DOUBLE, MPI_COMM_WORLD); } From 3b0317072e51a8fb7cf8667040cb3fb32ebec0fa Mon Sep 17 00:00:00 2001 From: koodlyakshay Date: Wed, 29 Jul 2020 15:39:23 +0200 Subject: [PATCH 28/29] Add test case config file. --- .../rough_flatplate_incomp.cfg | 291 ++++++++++++++++++ config_template.cfg | 2 +- 2 files changed, 292 insertions(+), 1 deletion(-) create mode 100644 TestCases/incomp_rans/rough_flatplate/rough_flatplate_incomp.cfg diff --git a/TestCases/incomp_rans/rough_flatplate/rough_flatplate_incomp.cfg b/TestCases/incomp_rans/rough_flatplate/rough_flatplate_incomp.cfg new file mode 100644 index 000000000000..0142f508d274 --- /dev/null +++ b/TestCases/incomp_rans/rough_flatplate/rough_flatplate_incomp.cfg @@ -0,0 +1,291 @@ +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% % +% SU2 configuration file % +% Case description: Turbulent flow over rough flat plate with zero % +% pressure gradient % +% Author: Akshay Koodly % +% Date: 2020.07.07 % +% File Version 7.0.x "Blackbird" % +% % +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +% ------------- DIRECT, ADJOINT, AND LINEARIZED PROBLEM DEFINITION ------------% +% +% Physical governing equations (EULER, NAVIER_STOKES, +% WAVE_EQUATION, HEAT_EQUATION, FEM_ELASTICITY, +% POISSON_EQUATION) +SOLVER= INC_RANS +% +% If Navier-Stokes, kind of turbulent model (NONE, SA) +KIND_TURB_MODEL= SA + +% Mathematical problem (DIRECT, CONTINUOUS_ADJOINT) +MATH_PROBLEM= DIRECT +% +% Restart solution (NO, YES) +RESTART_SOL= NO +% +% ---------------- INCOMPRESSIBLE FLOW CONDITION DEFINITION -------------------% +% +% Density model within the incompressible flow solver. +% Options are CONSTANT (default), BOUSSINESQ, or VARIABLE. If VARIABLE, +% an appropriate fluid model must be selected. +INC_DENSITY_MODEL= CONSTANT +% +% Solve the energy equation in the incompressible flow solver +INC_ENERGY_EQUATION = NO +% +% Initial density for incompressible flows (1.2886 kg/m^3 by default) +INC_DENSITY_INIT= 1.32905 +% +% Initial velocity for incompressible flows (1.0,0,0 m/s by default) +INC_VELOCITY_INIT= ( 69.4448, 0.0, 0.0 ) +% +% Initial temperature for incompressible flows that include the +% energy equation (288.15 K by default). Value is ignored if +% INC_ENERGY_EQUATION is false. +INC_TEMPERATURE_INIT= 300.0 +% +% Non-dimensionalization scheme for incompressible flows. Options are +% INITIAL_VALUES (default), REFERENCE_VALUES, or DIMENSIONAL. +% INC_*_REF values are ignored unless REFERENCE_VALUES is chosen. +INC_NONDIM= INITIAL_VALUES +% +% Reference density for incompressible flows (1.0 kg/m^3 by default) +INC_DENSITY_REF= 1.0 +% +% Reference velocity for incompressible flows (1.0 m/s by default) +INC_VELOCITY_REF= 1.0 +% +% Reference temperature for incompressible flows that include the +% energy equation (1.0 K by default) +INC_TEMPERATURE_REF = 1.0 +% +% List of inlet types for incompressible flows. List length must +% match number of inlet markers. Options: VELOCITY_INLET, PRESSURE_INLET. +INC_INLET_TYPE= VELOCITY_INLET + +% --------------------------- VISCOSITY MODEL ---------------------------------% +% +% Viscosity model (SUTHERLAND, CONSTANT_VISCOSITY). +VISCOSITY_MODEL= CONSTANT_VISCOSITY +% +% Molecular Viscosity that would be constant (1.716E-5 by default) +MU_CONSTANT= 1.84592e-05 +% +% Sutherland Viscosity Ref (1.716E-5 default value for AIR SI) +MU_REF= 1.716E-5 +% +% Sutherland Temperature Ref (273.15 K default value for AIR SI) +MU_T_REF= 273.15 +% +% Sutherland constant (110.4 default value for AIR SI) +SUTHERLAND_CONSTANT= 110.4 + +% ---------------------- REFERENCE VALUE DEFINITION ---------------------------% +% +% Reference origin for moment computation +REF_ORIGIN_MOMENT_X = 0.25 +REF_ORIGIN_MOMENT_Y = 0.00 +REF_ORIGIN_MOMENT_Z = 0.00 +% +% Reference length for pitching, rolling, and yawing non-dimensional moment +REF_LENGTH= 1.0 +% +% Reference area for force coefficients (0 implies automatic calculation) +REF_AREA= 2.0 + +% -------------------- BOUNDARY CONDITION DEFINITION --------------------------% +% +% Navier-Stokes wall boundary marker(s) (NONE = no marker) +MARKER_HEATFLUX= ( wall, 0.0 ) +%WALL_ROUGHNESS = (wall, 0.000061) +%WALL_ROUGHNESS = (wall, 0.000123) +WALL_ROUGHNESS = (wall, 0.000246) +%WALL_ROUGHNESS = (wall, 0.000984) +% +% +% Inlet boundary marker(s) (NONE = no marker) +% Format: ( inlet marker, total temperature, total pressure, flow_direction_x, +% flow_direction_y, flow_direction_z, ... ) +MARKER_INLET= ( inlet, 300.0, 69.4448, 1.0, 0.0, 0.0 ) +% +% Outlet boundary marker(s) (NONE = no marker) +% Format: ( outlet marker, back pressure, ... ) +MARKER_OUTLET= ( outlet, 0.0, farfield, 0.0 ) +% +INC_OUTLET_TYPE= PRESSURE_OUTLET,PRESSURE_OUTLET +% +% Symmetry boundary marker(s) (NONE = no marker) +MARKER_SYM= ( symmetry ) +% +% Marker(s) of the surface to be plotted or designed +MARKER_PLOTTING= ( wall ) +% +% Marker(s) of the surface where the functional (Cd, Cl, etc.) will be evaluated +MARKER_MONITORING= ( wall ) + +% ------------- COMMON PARAMETERS DEFINING THE NUMERICAL METHOD ---------------% +% +% Numerical method for spatial gradients (GREEN_GAUSS, LEAST_SQUARES, +% WEIGHTED_LEAST_SQUARES) +NUM_METHOD_GRAD= GREEN_GAUSS +% +% Courant-Friedrichs-Lewy condition of the finest grid +CFL_NUMBER= 100.0 +% +% Adaptive CFL number (NO, YES) +CFL_ADAPT= NO +% +% Parameters of the adaptive CFL number (factor down, factor up, CFL min value, +% CFL max value ) +CFL_ADAPT_PARAM= ( 1.5, 0.5, 1.1, 100.0 ) +% +% Runge-Kutta alpha coefficients +RK_ALPHA_COEFF= ( 0.66667, 0.66667, 1.000000 ) +% +% Number of total iterations +ITER= 1000 + +% ------------------------ LINEAR SOLVER DEFINITION ---------------------------% +% +% Linear solver for implicit formulations (BCGSTAB, FGMRES) +LINEAR_SOLVER= FGMRES +% +% Preconditioner of the Krylov linear solver (JACOBI, LINELET, LU_SGS) +LINEAR_SOLVER_PREC= ILU +% +% Linael solver ILU preconditioner fill-in level (0 by default) +LINEAR_SOLVER_ILU_FILL_IN= 0 +% +% Minimum error of the linear solver for implicit formulations +LINEAR_SOLVER_ERROR= 1E-12 +% +% Max number of iterations of the linear solver for the implicit formulation +LINEAR_SOLVER_ITER= 20 + +% ----------------------- SLOPE LIMITER DEFINITION ----------------------------% +% +% Coefficient for the limiter +VENKAT_LIMITER_COEFF= 0.1 +% +% Coefficient for the sharp edges limiter +ADJ_SHARP_LIMITER_COEFF= 3.0 +% +% Reference coefficient (sensitivity) for detecting sharp edges. +REF_SHARP_EDGES= 3.0 +% +% Remove sharp edges from the sensitivity evaluation (NO, YES) +SENS_REMOVE_SHARP= NO + +% -------------------------- MULTIGRID PARAMETERS -----------------------------% +% +% Multi-Grid Levels (0 = no multi-grid) +MGLEVEL= 0 +% +% Multi-grid cycle (V_CYCLE, W_CYCLE, FULLMG_CYCLE) +MGCYCLE= V_CYCLE +% +% Multi-grid pre-smoothing level +MG_PRE_SMOOTH= ( 1, 1, 1, 1 ) +% +% Multi-grid post-smoothing level +MG_POST_SMOOTH= ( 0, 0, 0, 0 ) +% +% Jacobi implicit smoothing of the correction +MG_CORRECTION_SMOOTH= ( 0, 0, 0, 0 ) +% +% Damping factor for the residual restriction +MG_DAMP_RESTRICTION= 0.8 +% +% Damping factor for the correction prolongation +MG_DAMP_PROLONGATION= 0.8 + +% -------------------- FLOW NUMERICAL METHOD DEFINITION -----------------------% +% +% Convective numerical method (JST, LAX-FRIEDRICH, CUSP, ROE, AUSM, HLLC, +% TURKEL_PREC, MSW) +CONV_NUM_METHOD_FLOW= FDS +% +% Monotonic Upwind Scheme for Conservation Laws (TVD) in the flow equations. +% Required for 2nd order upwind schemes (NO, YES) +MUSCL_FLOW= YES +% +% Slope limiter (NONE, VENKATAKRISHNAN, VENKATAKRISHNAN_WANG, +% BARTH_JESPERSEN, VAN_ALBADA_EDGE) +SLOPE_LIMITER_FLOW= NONE +% +% 2nd and 4th order artificial dissipation coefficients +JST_SENSOR_COEFF= ( 0.5, 0.02 ) +% +% Time discretization (RUNGE-KUTTA_EXPLICIT, EULER_IMPLICIT, EULER_EXPLICIT) +TIME_DISCRE_FLOW= EULER_IMPLICIT + +% -------------------- TURBULENT NUMERICAL METHOD DEFINITION ------------------% +% +% Convective numerical method (SCALAR_UPWIND) +CONV_NUM_METHOD_TURB= SCALAR_UPWIND +% +% Monotonic Upwind Scheme for Conservation Laws (TVD) in the turbulence equations. +% Required for 2nd order upwind schemes (NO, YES) +MUSCL_TURB= NO +% +% Slope limiter (VENKATAKRISHNAN, MINMOD) +SLOPE_LIMITER_TURB= VENKATAKRISHNAN +% +% Time discretization (EULER_IMPLICIT) +TIME_DISCRE_TURB= EULER_IMPLICIT +% --------------------------- CONVERGENCE PARAMETERS --------------------------% +% +% Convergence criteria (CAUCHY, RESIDUAL) +% +CONV_FIELD= RMS_VELOCITY-X +% +% Min value of the residual (log10 of the residual) +CONV_RESIDUAL_MINVAL= -14 +% +% Start convergence criteria at iteration number +CONV_STARTITER= 10 + +% ------------------------- INPUT/OUTPUT INFORMATION --------------------------% +% +% Mesh input file +%MESH_FILENAME= mesh_flatplate_turb_69x49.su2 +MESH_FILENAME= mesh_flatplate_turb_137x97.su2 +%MESH_FILENAME= mesh_flatplate_turb_273x193.su2 +% +% Mesh input file format (SU2, CGNS, NETCDF_ASCII) +MESH_FORMAT= SU2 +% +% Mesh output file +MESH_OUT_FILENAME= mesh_out.su2 +% +% Restart flow input file +SOLUTION_FILENAME= restart_flow +% +% Output file format (PARAVIEW, TECPLOT, SLT) +TABULAR_FORMAT= CSV + +OUTPUT_FILES= PARAVIEW, RESTART, SURFACE_PARAVIEW +% +% Output file convergence history (w/o extension) +CONV_FILENAME= history +% +% Output file restart flow +RESTART_FILENAME= restart_flow +% +% Output file flow (w/o extension) variables +VOLUME_FILENAME= flow +% +% Output file surface flow coefficient (w/o extension) +SURFACE_FILENAME= surface_flow +% +% Writing solution file frequency +OUTPUT_WRT_FREQ= 100 +% +% Writing convergence history frequency +WRT_CON_FREQ= 1 +% +SCREEN_OUTPUT= (WALL_TIME,INNER_ITER, RMS_VELOCITY-X, RMS_NU_TILDE, LIFT,DRAG) +% +WRT_FORCES_BREAKDOWN= YES diff --git a/config_template.cfg b/config_template.cfg index a1a94e2fb2dc..47825d0caede 100644 --- a/config_template.cfg +++ b/config_template.cfg @@ -816,7 +816,7 @@ SPATIAL_FOURIER= NO % The equivalent sand grain roughness height (k_s) on each of the wall. This must be in m. % This is a list of (string, double) each element corresponding to the MARKER defined in WALL_TYPE. %WALL_ROUGHNESS = (wall1, ks1, wall2, ks2) -%WALL_ROUGHNESS = (wall1, ks1, wall2, 0.0) is also allowed +%WALL_ROUGHNESS = (wall1, ks1, wall2, 0.0) %is also allowed % ------------------------ SURFACES IDENTIFICATION ----------------------------% % % Marker(s) of the surface in the surface flow solution file From c518c86862e5bad7cb6b0edfdb88cb49fb259620 Mon Sep 17 00:00:00 2001 From: Pedro Gomes Date: Wed, 29 Jul 2020 15:50:31 +0100 Subject: [PATCH 29/29] revert changes to codi and medi --- externals/codi | 2 +- externals/medi | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/externals/codi b/externals/codi index 501dcf0305df..b8ddede642eb 160000 --- a/externals/codi +++ b/externals/codi @@ -1 +1 @@ -Subproject commit 501dcf0305df147481630f20ce37c2e624fb351f +Subproject commit b8ddede642eb2824a56c75cebddb79ba194c2ada diff --git a/externals/medi b/externals/medi index edde14f9ac40..051753df5bb6 160000 --- a/externals/medi +++ b/externals/medi @@ -1 +1 @@ -Subproject commit edde14f9ac4026b72b1e130f61c0a78e8652afa5 +Subproject commit 051753df5bb6e66c9fb6fbe4f2413c4ba187e81d