diff --git a/Common/include/config_structure.hpp b/Common/include/config_structure.hpp index 8bc4e3438695..68cfab7e930d 100644 --- a/Common/include/config_structure.hpp +++ b/Common/include/config_structure.hpp @@ -532,7 +532,8 @@ class CConfig { Kind_Matrix_Coloring, /*!< \brief Type of matrix coloring for sparse Jacobian computation. */ Kind_Solver_Fluid_FSI, /*!< \brief Kind of solver for the fluid in FSI applications. */ Kind_Solver_Struc_FSI, /*!< \brief Kind of solver for the structure in FSI applications. */ - Kind_BGS_RelaxMethod; /*!< \brief Kind of relaxation method for Block Gauss Seidel method in FSI problems. */ + Kind_BGS_RelaxMethod, /*!< \brief Kind of relaxation method for Block Gauss Seidel method in FSI problems. */ + Kind_CHT_Coupling; /*!< \brief Kind of coupling method used at CHT interfaces. */ bool ReconstructionGradientRequired; /*!< \brief Enable or disable a second gradient calculation for upwind reconstruction only. */ bool LeastSquaresRequired; /*!< \brief Enable or disable memory allocation for least-squares gradient methods. */ bool Energy_Equation; /*!< \brief Solve the energy equation for incompressible flows. */ @@ -903,7 +904,6 @@ class CConfig { bool Sine_Load; /*!< \brief option for sine load */ su2double *SineLoad_Coeff; /*!< \brief Stores the load coefficient */ su2double Thermal_Diffusivity; /*!< \brief Thermal diffusivity used in the heat solver. */ - bool CHT_Robin; /*!< \brief Option for boundary condition method at CHT interfaces. */ su2double Cyclic_Pitch, /*!< \brief Cyclic pitch for rotorcraft simulations. */ Collective_Pitch; /*!< \brief Collective pitch for rotorcraft simulations. */ su2double Mach_Motion; /*!< \brief Mach number based on mesh velocity and freestream quantities. */ @@ -8985,10 +8985,10 @@ class CConfig { bool GetWeakly_Coupled_Heat(void); /*! - * \brief Get the boundary condition method for CHT. - * \return YES if Robin BC is used. + * \brief Get the CHT couling method. + * \return Kind of the method. */ - bool GetCHT_Robin(void); + unsigned short GetKind_CHT_Coupling(void) const; /*! * \brief Check if values passed to the BC_HeatFlux-Routine are already integrated. diff --git a/Common/include/config_structure.inl b/Common/include/config_structure.inl index e5ce3843521b..01c99f8d31f7 100644 --- a/Common/include/config_structure.inl +++ b/Common/include/config_structure.inl @@ -1917,7 +1917,7 @@ inline bool CConfig::GetJacobian_Spatial_Discretization_Only(void) {return Jacob inline bool CConfig::GetWeakly_Coupled_Heat(void) { return Weakly_Coupled_Heat; } -inline bool CConfig::GetCHT_Robin(void) { return CHT_Robin; } +inline unsigned short CConfig::GetKind_CHT_Coupling(void) const { return Kind_CHT_Coupling; } inline bool CConfig::GetIntegrated_HeatFlux(void) { return Integrated_HeatFlux; } diff --git a/Common/include/option_structure.hpp b/Common/include/option_structure.hpp index 9c0db07466e4..8f8f45939b8c 100644 --- a/Common/include/option_structure.hpp +++ b/Common/include/option_structure.hpp @@ -1112,6 +1112,23 @@ static const map DVFEA_Map = CCreateMap ("DEAD_WEIGHT", DEAD_WEIGHT) ("ELECTRIC_FIELD", ELECTRIC_FIELD); +/*! + * \brief Kinds of coupling methods at CHT interfaces. + * The first (temperature) part determines the BC method on the fluid side, the second (heatflux) part determines + * the BC method on the solid side of the CHT interface. + */ +enum ENUM_CHT_COUPLING { + DIRECT_TEMPERATURE_NEUMANN_HEATFLUX = 0, + AVERAGED_TEMPERATURE_NEUMANN_HEATFLUX = 1, + DIRECT_TEMPERATURE_ROBIN_HEATFLUX = 2, + AVERAGED_TEMPERATURE_ROBIN_HEATFLUX = 3 +}; +static const map CHT_Coupling_Map = CCreateMap +("DIRECT_TEMPERATURE_NEUMANN_HEATFLUX", DIRECT_TEMPERATURE_NEUMANN_HEATFLUX) +("AVERAGED_TEMPERATURE_NEUMANN_HEATFLUX", AVERAGED_TEMPERATURE_NEUMANN_HEATFLUX) +("DIRECT_TEMPERATURE_ROBIN_HEATFLUX", DIRECT_TEMPERATURE_ROBIN_HEATFLUX) +("AVERAGED_TEMPERATURE_ROBIN_HEATFLUX", AVERAGED_TEMPERATURE_ROBIN_HEATFLUX); + /*! * \brief types Riemann boundary treatments */ diff --git a/Common/src/config_structure.cpp b/Common/src/config_structure.cpp index f2a8727826b7..748b54968c4e 100644 --- a/Common/src/config_structure.cpp +++ b/Common/src/config_structure.cpp @@ -2268,9 +2268,9 @@ void CConfig::SetConfig_Options() { /*!\par CONFIG_CATEGORY: Heat solver \ingroup Config*/ /*--- options related to the heat solver ---*/ - /* DESCRIPTION: Use Robin (default) or Neumann BC at CHT interface. */ + /* DESCRIPTION: CHT interface coupling methods */ /* Options: NO, YES \ingroup Config */ - addBoolOption("CHT_ROBIN", CHT_Robin, true); + addEnumOption("CHT_COUPLING_METHOD", Kind_CHT_Coupling, CHT_Coupling_Map, DIRECT_TEMPERATURE_ROBIN_HEATFLUX); /* DESCRIPTION: Thermal diffusivity constant */ addDoubleOption("THERMAL_DIFFUSIVITY", Thermal_Diffusivity, 1.172E-5); diff --git a/SU2_CFD/src/drivers/CDriver.cpp b/SU2_CFD/src/drivers/CDriver.cpp index 455a95b5c83b..efc728ac9327 100644 --- a/SU2_CFD/src/drivers/CDriver.cpp +++ b/SU2_CFD/src/drivers/CDriver.cpp @@ -3547,7 +3547,7 @@ void CDriver::Interface_Preprocessing(CConfig **config, CSolver***** solver, CGe else if (fluid_donor && heat_target) { nVarTransfer = 0; nVar = 4; - if(config[donorZone]->GetEnergy_Equation()) + if(config[donorZone]->GetEnergy_Equation() || (config[donorZone]->GetKind_Regime() == COMPRESSIBLE)) interface_types[donorZone][targetZone] = CONJUGATE_HEAT_FS; else if (config[donorZone]->GetWeakly_Coupled_Heat()) interface_types[donorZone][targetZone] = CONJUGATE_HEAT_WEAKLY_FS; @@ -3558,7 +3558,7 @@ void CDriver::Interface_Preprocessing(CConfig **config, CSolver***** solver, CGe else if (heat_donor && fluid_target) { nVarTransfer = 0; nVar = 4; - if(config[targetZone]->GetEnergy_Equation()) + if(config[targetZone]->GetEnergy_Equation() || (config[targetZone]->GetKind_Regime() == COMPRESSIBLE)) interface_types[donorZone][targetZone] = CONJUGATE_HEAT_SF; else if (config[targetZone]->GetWeakly_Coupled_Heat()) interface_types[donorZone][targetZone] = CONJUGATE_HEAT_WEAKLY_SF; diff --git a/SU2_CFD/src/interfaces/cht/CConjugateHeatInterface.cpp b/SU2_CFD/src/interfaces/cht/CConjugateHeatInterface.cpp index c33ef4602ca9..835ad3dcef25 100644 --- a/SU2_CFD/src/interfaces/cht/CConjugateHeatInterface.cpp +++ b/SU2_CFD/src/interfaces/cht/CConjugateHeatInterface.cpp @@ -53,19 +53,17 @@ void CConjugateHeatInterface::GetDonor_Variable(CSolver *donor_solution, CGeomet /*--- Check whether the current zone is a solid zone or a fluid zone ---*/ - bool flow = ((donor_config->GetKind_Solver() == NAVIER_STOKES) + bool compressible_flow = ((donor_config->GetKind_Solver() == NAVIER_STOKES) || (donor_config->GetKind_Solver() == RANS) || (donor_config->GetKind_Solver() == DISC_ADJ_NAVIER_STOKES) - || (donor_config->GetKind_Solver() == DISC_ADJ_RANS) - || (donor_config->GetKind_Solver() == INC_NAVIER_STOKES) + || (donor_config->GetKind_Solver() == DISC_ADJ_RANS)); + bool incompressible_flow = ((donor_config->GetKind_Solver() == INC_NAVIER_STOKES) || (donor_config->GetKind_Solver() == INC_RANS) || (donor_config->GetKind_Solver() == DISC_ADJ_INC_NAVIER_STOKES) - || (donor_config->GetKind_Solver() == DISC_ADJ_INC_RANS)); - - bool compressible_flow = (donor_config->GetKind_Regime() == COMPRESSIBLE) && flow; - bool incompressible_flow = (donor_config->GetEnergy_Equation()) && flow; - bool heat_equation = (donor_config->GetKind_Solver() == HEAT_EQUATION_FVM - || donor_config->GetKind_Solver() == DISC_ADJ_HEAT); + || (donor_config->GetKind_Solver() == DISC_ADJ_INC_RANS) + && (donor_config->GetEnergy_Equation())); + bool heat_equation = (donor_config->GetKind_Solver() == HEAT_EQUATION_FVM + || donor_config->GetKind_Solver() == DISC_ADJ_HEAT); Coord = donor_geometry->node[Point_Donor]->GetCoord(); @@ -130,7 +128,8 @@ void CConjugateHeatInterface::GetDonor_Variable(CSolver *donor_solution, CGeomet thermal_conductivityND = Cp*(laminar_viscosity/Prandtl_Lam); heat_flux_density = thermal_conductivityND*dTdn; - if (donor_config->GetCHT_Robin()) { + if ((donor_config->GetKind_CHT_Coupling() == DIRECT_TEMPERATURE_ROBIN_HEATFLUX) || + (donor_config->GetKind_CHT_Coupling() == AVERAGED_TEMPERATURE_ROBIN_HEATFLUX)) { thermal_conductivity = thermal_conductivityND*donor_config->GetViscosity_Ref(); conductivity_over_dist = thermal_conductivity/dist; @@ -143,7 +142,8 @@ void CConjugateHeatInterface::GetDonor_Variable(CSolver *donor_solution, CGeomet thermal_conductivityND = donor_solution->GetNodes()->GetThermalConductivity(iPoint); heat_flux_density = thermal_conductivityND*dTdn; - if (donor_config->GetCHT_Robin()) { + if ((donor_config->GetKind_CHT_Coupling() == DIRECT_TEMPERATURE_ROBIN_HEATFLUX) || + (donor_config->GetKind_CHT_Coupling() == AVERAGED_TEMPERATURE_ROBIN_HEATFLUX)) { switch (donor_config->GetKind_ConductivityModel()) { @@ -167,7 +167,9 @@ void CConjugateHeatInterface::GetDonor_Variable(CSolver *donor_solution, CGeomet thermal_diffusivity = donor_config->GetThermalDiffusivity_Solid(); heat_flux_density = thermal_diffusivity*dTdn; - if (donor_config->GetCHT_Robin()) { + if ((donor_config->GetKind_CHT_Coupling() == DIRECT_TEMPERATURE_ROBIN_HEATFLUX) || + (donor_config->GetKind_CHT_Coupling() == AVERAGED_TEMPERATURE_ROBIN_HEATFLUX)) { + rho_cp_solid = donor_config->GetSpecific_Heat_Cp()*donor_config->GetDensity_Solid(); conductivity_over_dist = thermal_diffusivity*rho_cp_solid/dist; } @@ -180,7 +182,8 @@ void CConjugateHeatInterface::GetDonor_Variable(CSolver *donor_solution, CGeomet /*--- We only need these for the Robin BC option ---*/ - if (donor_config->GetCHT_Robin()) { + if ((donor_config->GetKind_CHT_Coupling() == DIRECT_TEMPERATURE_ROBIN_HEATFLUX) || + (donor_config->GetKind_CHT_Coupling() == AVERAGED_TEMPERATURE_ROBIN_HEATFLUX)) { Donor_Variable[2] = conductivity_over_dist; Donor_Variable[3] = Tnormal*donor_config->GetTemperature_Ref(); @@ -201,7 +204,8 @@ void CConjugateHeatInterface::SetTarget_Variable(CSolver *target_solution, CGeom target_solution->SetConjugateHeatVariable(Marker_Target, Vertex_Target, 1, target_config->GetRelaxation_Factor_CHT(), Target_Variable[1]); - if (target_config->GetCHT_Robin()) { + if ((target_config->GetKind_CHT_Coupling() == DIRECT_TEMPERATURE_ROBIN_HEATFLUX) || + (target_config->GetKind_CHT_Coupling() == AVERAGED_TEMPERATURE_ROBIN_HEATFLUX)) { target_solution->SetConjugateHeatVariable(Marker_Target, Vertex_Target, 2, target_config->GetRelaxation_Factor_CHT(), Target_Variable[2]); diff --git a/SU2_CFD/src/solver_direct_heat.cpp b/SU2_CFD/src/solver_direct_heat.cpp index c11890658dba..83feeeb10058 100644 --- a/SU2_CFD/src/solver_direct_heat.cpp +++ b/SU2_CFD/src/solver_direct_heat.cpp @@ -1113,32 +1113,27 @@ void CHeatSolverFVM::BC_ConjugateHeat_Interface(CGeometry *geometry, CSolver **s if (flow) { - for (iMarker = 0; iMarker < config->GetnMarker_All(); iMarker++) { + for (iVertex = 0; iVertex < geometry->nVertex[val_marker]; iVertex++) { - if (config->GetMarker_All_KindBC(iMarker) == CHT_WALL_INTERFACE) { + iPoint = geometry->vertex[val_marker][iVertex]->GetNode(); - for (iVertex = 0; iVertex < geometry->nVertex[iMarker]; iVertex++) { - iPoint = geometry->vertex[iMarker][iVertex]->GetNode(); - - if (geometry->node[iPoint]->GetDomain()) { + if (geometry->node[iPoint]->GetDomain()) { - Normal = geometry->vertex[iMarker][iVertex]->GetNormal(); - Area = 0.0; - for (iDim = 0; iDim < nDim; iDim++) Area += Normal[iDim]*Normal[iDim]; - Area = sqrt (Area); + Normal = geometry->vertex[val_marker][iVertex]->GetNormal(); + Area = 0.0; + for (iDim = 0; iDim < nDim; iDim++) Area += Normal[iDim]*Normal[iDim]; + Area = sqrt (Area); - T_Conjugate = GetConjugateHeatVariable(iMarker, iVertex, 0)/Temperature_Ref; + T_Conjugate = GetConjugateHeatVariable(val_marker, iVertex, 0)/Temperature_Ref; - nodes->SetSolution_Old(iPoint,&T_Conjugate); - LinSysRes.SetBlock_Zero(iPoint, 0); - nodes->SetRes_TruncErrorZero(iPoint); + nodes->SetSolution_Old(iPoint,&T_Conjugate); + LinSysRes.SetBlock_Zero(iPoint, 0); + nodes->SetRes_TruncErrorZero(iPoint); - if (implicit) { - for (iVar = 0; iVar < nVar; iVar++) { - total_index = iPoint*nVar+iVar; - Jacobian.DeleteValsRowi(total_index); - } - } + if (implicit) { + for (iVar = 0; iVar < nVar; iVar++) { + total_index = iPoint*nVar+iVar; + Jacobian.DeleteValsRowi(total_index); } } } @@ -1146,45 +1141,41 @@ void CHeatSolverFVM::BC_ConjugateHeat_Interface(CGeometry *geometry, CSolver **s } else { - for (iMarker = 0; iMarker < config->GetnMarker_All(); iMarker++) { + for (iVertex = 0; iVertex < geometry->nVertex[val_marker]; iVertex++) { - if (config->GetMarker_All_KindBC(iMarker) == CHT_WALL_INTERFACE) { + iPoint = geometry->vertex[val_marker][iVertex]->GetNode(); - for (iVertex = 0; iVertex < geometry->nVertex[iMarker]; iVertex++) { - iPoint = geometry->vertex[iMarker][iVertex]->GetNode(); - - if (geometry->node[iPoint]->GetDomain()) { + if (geometry->node[iPoint]->GetDomain()) { - Normal = geometry->vertex[iMarker][iVertex]->GetNormal(); - Area = 0.0; - for (iDim = 0; iDim < nDim; iDim++) Area += Normal[iDim]*Normal[iDim]; - Area = sqrt(Area); + Normal = geometry->vertex[val_marker][iVertex]->GetNormal(); + Area = 0.0; + for (iDim = 0; iDim < nDim; iDim++) Area += Normal[iDim]*Normal[iDim]; + Area = sqrt(Area); - thermal_diffusivity = GetConjugateHeatVariable(iMarker, iVertex, 2)/rho_cp_solid; + thermal_diffusivity = GetConjugateHeatVariable(val_marker, iVertex, 2)/rho_cp_solid; - if (config->GetCHT_Robin()) { + if ((config->GetKind_CHT_Coupling() == DIRECT_TEMPERATURE_ROBIN_HEATFLUX) || + (config->GetKind_CHT_Coupling() == AVERAGED_TEMPERATURE_ROBIN_HEATFLUX)) { - Tinterface = nodes->GetSolution(iPoint,0); - Tnormal_Conjugate = GetConjugateHeatVariable(iMarker, iVertex, 3)/Temperature_Ref; + Tinterface = nodes->GetSolution(iPoint,0); + Tnormal_Conjugate = GetConjugateHeatVariable(val_marker, iVertex, 3)/Temperature_Ref; - HeatFluxDensity = thermal_diffusivity*(Tinterface - Tnormal_Conjugate); - HeatFlux = HeatFluxDensity * Area; - } - else { + HeatFluxDensity = thermal_diffusivity*(Tinterface - Tnormal_Conjugate); + HeatFlux = HeatFluxDensity * Area; + } + else { - HeatFluxDensity = GetConjugateHeatVariable(iMarker, iVertex, 1)/config->GetHeat_Flux_Ref(); - HeatFlux = HeatFluxDensity*Area; - } + HeatFluxDensity = GetConjugateHeatVariable(val_marker, iVertex, 1)/config->GetHeat_Flux_Ref(); + HeatFlux = HeatFluxDensity*Area; + } - Res_Visc[0] = -HeatFlux; - LinSysRes.SubtractBlock(iPoint, Res_Visc); + Res_Visc[0] = -HeatFlux; + LinSysRes.SubtractBlock(iPoint, Res_Visc); - if (implicit) { + if (implicit) { - Jacobian_i[0][0] = thermal_diffusivity*Area; - Jacobian.SubtractBlock(iPoint, iPoint, Jacobian_i); - } - } + Jacobian_i[0][0] = thermal_diffusivity*Area; + Jacobian.SubtractBlock(iPoint, iPoint, Jacobian_i); } } } diff --git a/SU2_CFD/src/solver_direct_mean.cpp b/SU2_CFD/src/solver_direct_mean.cpp index 82360d43437d..d96bd40782b8 100644 --- a/SU2_CFD/src/solver_direct_mean.cpp +++ b/SU2_CFD/src/solver_direct_mean.cpp @@ -16003,6 +16003,8 @@ void CNSSolver::BC_ConjugateHeat_Interface(CGeometry *geometry, CSolver **solver su2double Gas_Constant = config->GetGas_ConstantND(); su2double Cp = (Gamma / Gamma_Minus_One) * Gas_Constant; + su2double Temperature_Ref = config->GetTemperature_Ref(); + bool implicit = (config->GetKind_TimeIntScheme_Flow() == EULER_IMPLICIT); /*--- Identify the boundary ---*/ @@ -16092,17 +16094,31 @@ void CNSSolver::BC_ConjugateHeat_Interface(CGeometry *geometry, CSolver **solver /*--- Compute the normal gradient in temperature using Twall ---*/ There = nodes->GetTemperature(Point_Normal); - Tconjugate = GetConjugateHeatVariable(val_marker, iVertex, 0); + Tconjugate = GetConjugateHeatVariable(val_marker, iVertex, 0)/Temperature_Ref; + + if ((config->GetKind_CHT_Coupling() == AVERAGED_TEMPERATURE_NEUMANN_HEATFLUX) || + (config->GetKind_CHT_Coupling() == AVERAGED_TEMPERATURE_ROBIN_HEATFLUX)) { - HF_FactorHere = thermal_conductivity*config->GetViscosity_Ref()/dist_ij; - HF_FactorConjugate = GetConjugateHeatVariable(val_marker, iVertex, 2); + /*--- Compute wall temperature from both temperatures ---*/ - Twall = (There*HF_FactorHere + Tconjugate*HF_FactorConjugate)/(HF_FactorHere + HF_FactorConjugate); + HF_FactorHere = thermal_conductivity*config->GetViscosity_Ref()/dist_ij; + HF_FactorConjugate = GetConjugateHeatVariable(val_marker, iVertex, 2); - // this will be changed soon... - Twall = GetConjugateHeatVariable(val_marker, iVertex, 0); + Twall = (There*HF_FactorHere + Tconjugate*HF_FactorConjugate)/(HF_FactorHere + HF_FactorConjugate); + dTdn = -(There - Twall)/dist_ij; + } + else if ((config->GetKind_CHT_Coupling() == DIRECT_TEMPERATURE_NEUMANN_HEATFLUX) || + (config->GetKind_CHT_Coupling() == DIRECT_TEMPERATURE_ROBIN_HEATFLUX)) { - dTdn = -(There - Twall)/dist_ij; + /*--- (Directly) Set wall temperature to conjugate temperature. ---*/ + + Twall = Tconjugate; + dTdn = -(There - Twall)/dist_ij; + } + else { + + SU2_MPI::Error(string("Unknown CHT coupling method."), CURRENT_FUNCTION); + } /*--- Apply a weak boundary condition for the energy equation. Compute the residual due to the prescribed heat flux. ---*/ @@ -16281,7 +16297,6 @@ void CNSSolver::BC_ConjugateHeat_Interface(CGeometry *geometry, CSolver **solver Jacobian.DeleteValsRowi(total_index); } } - } } } diff --git a/SU2_CFD/src/solver_direct_mean_inc.cpp b/SU2_CFD/src/solver_direct_mean_inc.cpp index 0e2a73da82ac..7bec4988809e 100644 --- a/SU2_CFD/src/solver_direct_mean_inc.cpp +++ b/SU2_CFD/src/solver_direct_mean_inc.cpp @@ -7859,10 +7859,12 @@ void CIncNSSolver::BC_Isothermal_Wall(CGeometry *geometry, CSolver **solver_cont void CIncNSSolver::BC_ConjugateHeat_Interface(CGeometry *geometry, CSolver **solver_container, CNumerics *conv_numerics, CConfig *config, unsigned short val_marker) { unsigned short iVar, jVar, iDim, Wall_Function; - unsigned long iVertex, iPoint, total_index; + unsigned long iVertex, iPoint, total_index, Point_Normal; + + su2double *Coord_i, *Coord_j, dist_ij; + su2double *GridVel, There, Tconjugate, Twall, Temperature_Ref, thermal_conductivity, HF_FactorHere, HF_FactorConjugate; - su2double *GridVel, Tconjugate; - su2double Temperature_Ref = config->GetTemperature_Ref(); + Temperature_Ref = config->GetTemperature_Ref(); bool implicit = (config->GetKind_TimeIntScheme_Flow() == EULER_IMPLICIT); bool energy = config->GetEnergy_Equation(); @@ -7921,12 +7923,48 @@ void CIncNSSolver::BC_ConjugateHeat_Interface(CGeometry *geometry, CSolver **sol Tconjugate = GetConjugateHeatVariable(val_marker, iVertex, 0)/Temperature_Ref; + if ((config->GetKind_CHT_Coupling() == AVERAGED_TEMPERATURE_NEUMANN_HEATFLUX) || + (config->GetKind_CHT_Coupling() == AVERAGED_TEMPERATURE_ROBIN_HEATFLUX)) { + + /*--- Compute closest normal neighbor ---*/ + + Point_Normal = geometry->vertex[val_marker][iVertex]->GetNormal_Neighbor(); + + /*--- Get coordinates of i & nearest normal and compute distance ---*/ + + Coord_i = geometry->node[iPoint]->GetCoord(); + Coord_j = geometry->node[Point_Normal]->GetCoord(); + dist_ij = 0; + for (iDim = 0; iDim < nDim; iDim++) + dist_ij += (Coord_j[iDim]-Coord_i[iDim])*(Coord_j[iDim]-Coord_i[iDim]); + dist_ij = sqrt(dist_ij); + + /*--- Compute wall temperature from both temperatures ---*/ + + thermal_conductivity = nodes->GetThermalConductivity(iPoint); + There = nodes->GetTemperature(Point_Normal); + HF_FactorHere = thermal_conductivity*config->GetViscosity_Ref()/dist_ij; + HF_FactorConjugate = GetConjugateHeatVariable(val_marker, iVertex, 2); + + Twall = (There*HF_FactorHere + Tconjugate*HF_FactorConjugate)/(HF_FactorHere + HF_FactorConjugate); + } + else if ((config->GetKind_CHT_Coupling() == DIRECT_TEMPERATURE_NEUMANN_HEATFLUX) || + (config->GetKind_CHT_Coupling() == DIRECT_TEMPERATURE_ROBIN_HEATFLUX)) { + + /*--- (Directly) Set wall temperature to conjugate temperature. ---*/ + + Twall = Tconjugate; + } + else { + + SU2_MPI::Error(string("Unknown CHT coupling method."), CURRENT_FUNCTION); + } + /*--- Strong imposition of the temperature on the fluid zone. ---*/ LinSysRes.SetBlock_Zero(iPoint, nDim+1); - nodes->SetSolution_Old(iPoint, nDim+1, Tconjugate); + nodes->SetSolution_Old(iPoint, nDim+1, Twall); nodes->SetEnergy_ResTruncError_Zero(iPoint); - } /*--- Enforce the no-slip boundary condition in a strong way by @@ -7942,7 +7980,6 @@ void CIncNSSolver::BC_ConjugateHeat_Interface(CGeometry *geometry, CSolver **sol Jacobian.DeleteValsRowi(total_index); } } - } } } diff --git a/TestCases/coupled_cht/comp_2d/cht_2d_3cylinders.cfg b/TestCases/coupled_cht/comp_2d/cht_2d_3cylinders.cfg new file mode 100644 index 000000000000..5377b2c322a8 --- /dev/null +++ b/TestCases/coupled_cht/comp_2d/cht_2d_3cylinders.cfg @@ -0,0 +1,79 @@ +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% % +% SU2 configuration file % +% Case description: Steady compressible laminar flow around heated cylinders % +% Author: O. Burghardt % +% Institution: Chair for Scientific Computing, TU Kaiserslautern % +% Date: January 6, 2020 % +% File Version 7.0 "Blackbird" % +% % +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +% +% Physical governing equations (EULER, NAVIER_STOKES, +% WAVE_EQUATION, HEAT_EQUATION, FEM_ELASTICITY, +% POISSON_EQUATION) +SOLVER= MULTIPHYSICS +% +% Mathematical problem (DIRECT, CONTINUOUS_ADJOINT, DISCRETE_ADJOINT) +MATH_PROBLEM= DIRECT +% +% +CONFIG_LIST = (flow_cylinder.cfg, solid_cylinder1.cfg, solid_cylinder2.cfg, solid_cylinder3.cfg) +% +% +MARKER_ZONE_INTERFACE= (cylinder_outer1, cylinder_inner1, cylinder_outer2, cylinder_inner2, cylinder_outer3, cylinder_inner3) +% +% +MARKER_CHT_INTERFACE= (cylinder_outer1, cylinder_inner1, cylinder_outer2, cylinder_inner2, cylinder_outer3, cylinder_inner3) +% +% +TIME_DOMAIN = NO +% +% Number of total iterations +OUTER_ITER = 11 +% +% Mesh input file +MESH_FILENAME= mesh_cht_3cyl.su2 +% +% Mesh input file format (SU2, CGNS, NETCDF_ASCII) +MESH_FORMAT= SU2 + +% These are just default parameters so that we can run SU2_DOT_AD, they have no physical meaning for this test case. + +% ----------------------- DESIGN VARIABLE PARAMETERS --------------------------% +% +% Kind of deformation (NO_DEFORMATION, TRANSLATION, ROTATION, SCALE, +% FFD_SETTING, FFD_NACELLE +% FFD_CONTROL_POINT, FFD_CAMBER, FFD_THICKNESS, FFD_TWIST +% FFD_CONTROL_POINT_2D, FFD_CAMBER_2D, FFD_THICKNESS_2D, FFD_TWIST_2D, +% HICKS_HENNE, SURFACE_BUMP) +DV_KIND= HICKS_HENNE +% +% Marker of the surface in which we are going apply the shape deformation +DV_MARKER= ( cylinder_outer1, cylinder_inner1, cylinder_outer2, cylinder_inner2, cylinder_outer3, cylinder_inner3 ) +% +% Parameters of the shape deformation +% - NO_DEFORMATION ( 1.0 ) +% - TRANSLATION ( x_Disp, y_Disp, z_Disp ), as a unit vector +% - ROTATION ( x_Orig, y_Orig, z_Orig, x_End, y_End, z_End ) +% - SCALE ( 1.0 ) +% - ANGLE_OF_ATTACK ( 1.0 ) +% - FFD_SETTING ( 1.0 ) +% - FFD_CONTROL_POINT ( FFD_BoxTag, i_Ind, j_Ind, k_Ind, x_Disp, y_Disp, z_Disp ) +% - FFD_NACELLE ( FFD_BoxTag, rho_Ind, theta_Ind, phi_Ind, rho_Disp, phi_Disp ) +% - FFD_GULL ( FFD_BoxTag, j_Ind ) +% - FFD_ANGLE_OF_ATTACK ( FFD_BoxTag, 1.0 ) +% - FFD_CAMBER ( FFD_BoxTag, i_Ind, j_Ind ) +% - FFD_THICKNESS ( FFD_BoxTag, i_Ind, j_Ind ) +% - FFD_TWIST ( FFD_BoxTag, j_Ind, x_Orig, y_Orig, z_Orig, x_End, y_End, z_End ) +% - FFD_CONTROL_POINT_2D ( FFD_BoxTag, i_Ind, j_Ind, x_Disp, y_Disp ) +% - FFD_CAMBER_2D ( FFD_BoxTag, i_Ind ) +% - FFD_THICKNESS_2D ( FFD_BoxTag, i_Ind ) +% - FFD_TWIST_2D ( FFD_BoxTag, x_Orig, y_Orig ) +% - HICKS_HENNE ( Lower Surface (0)/Upper Surface (1)/Only one Surface (2), x_Loc ) +% - SURFACE_BUMP ( x_Start, x_End, x_Loc ) +DV_PARAM= (0.0, 0.5) +% +% Value of the shape deformation +DV_VALUE= 0.1 diff --git a/TestCases/coupled_cht/comp_2d/flow_cylinder.cfg b/TestCases/coupled_cht/comp_2d/flow_cylinder.cfg new file mode 100644 index 000000000000..0e054d3af89e --- /dev/null +++ b/TestCases/coupled_cht/comp_2d/flow_cylinder.cfg @@ -0,0 +1,207 @@ +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% % +% SU2 configuration file % +% Case description: Steady compressible laminar flow around heated cylinders % +% % +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +% ------------- DIRECT, ADJOINT, AND LINEARIZED PROBLEM DEFINITION ------------% +% +% Physical governing equations (EULER, NAVIER_STOKES, +% WAVE_EQUATION, HEAT_EQUATION, FEM_ELASTICITY, +% POISSON_EQUATION) +SOLVER= NAVIER_STOKES +% +% If Navier-Stokes, kind of turbulent model (NONE, SA) +KIND_TURB_MODEL= NONE +% +% Restart solution (NO, YES) +RESTART_SOL= NO +% +% Objective function in gradient evaluation (DRAG, LIFT, SIDEFORCE, MOMENT_X, +% MOMENT_Y, MOMENT_Z, EFFICIENCY, +% EQUIVALENT_AREA, NEARFIELD_PRESSURE, +% FORCE_X, FORCE_Y, FORCE_Z, THRUST, +% TORQUE, TOTAL_HEATFLUX, +% MAXIMUM_HEATFLUX, INVERSE_DESIGN_PRESSURE, +% INVERSE_DESIGN_HEATFLUX, SURFACE_TOTAL_PRESSURE, +% SURFACE_MASSFLOW, SURFACE_STATIC_PRESSURE, SURFACE_MACH) +% For a weighted sum of objectives: separate by commas, add OBJECTIVE_WEIGHT and MARKER_MONITORING in matching order. +OBJECTIVE_FUNCTION= TOTAL_HEATFLUX +% +% List of weighting values when using more than one OBJECTIVE_FUNCTION. Separate by commas and match with MARKER_MONITORING. +OBJECTIVE_WEIGHT = 1.0 +% +% Read binary restart files (YES, NO) +READ_BINARY_RESTART = YES +% +% Data written to history file +HISTORY_OUTPUT=( ITER, RMS_RES, HEAT ) + +% -------------------- COMPRESSIBLE FREE-STREAM DEFINITION --------------------% +% +% Mach number (non-dimensional, based on the free-stream values) +MACH_NUMBER= 0.2 +% +% Angle of attack (degrees, only for compressible flows) +AOA= 0.0 +% +% Side-slip angle (degrees, only for compressible flows) +SIDESLIP_ANGLE= 0.0 +% +% Free-stream temperature (288.15 K by default) +FREESTREAM_TEMPERATURE= 297.62 +% +% Reynolds number (non-dimensional, based on the free-stream values) +REYNOLDS_NUMBER= 100 +% +% Reynolds length (1 m by default) +REYNOLDS_LENGTH= 1.0 + +% -------------------- BOUNDARY CONDITION DEFINITION --------------------------% +% +% Farfield boundary marker(s) (NONE = no marker) +MARKER_FAR= ( farfield ) +% +% Marker(s) of the surface to be plotted or designed +MARKER_PLOTTING= ( cylinder_outer1, cylinder_outer2, cylinder_outer3 ) +% +% Marker(s) of the surface where the functional (Cd, Cl, etc.) will be evaluated +MARKER_MONITORING= ( cylinder_outer1, cylinder_outer2, cylinder_outer3 ) + +% ------------- COMMON PARAMETERS DEFINING THE NUMERICAL METHOD ---------------% +% +% Numerical method for spatial gradients (GREEN_GAUSS, WEIGHTED_LEAST_SQUARES) +NUM_METHOD_GRAD= GREEN_GAUSS +% +% Courant-Friedrichs-Lewy condition of the finest grid +CFL_NUMBER= 10.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, 10.0, 10000.0 ) +% +% Runge-Kutta alpha coefficients +RK_ALPHA_COEFF= ( 0.66667, 0.66667, 1.000000 ) + +% ------------------------ LINEAR SOLVER DEFINITION ---------------------------% +% +% Linear solver or smoother for implicit formulations (BCGSTAB, FGMRES, SMOOTHER_JACOBI, +% SMOOTHER_ILU, SMOOTHER_LUSGS, +% SMOOTHER_LINELET) +LINEAR_SOLVER= FGMRES +% +% Preconditioner of the Krylov linear solver (ILU, LU_SGS, LINELET, JACOBI) +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-15 +% +% Max number of iterations of the linear solver for the implicit formulation +LINEAR_SOLVER_ITER= 5 + +% -------------------- FLOW NUMERICAL METHOD DEFINITION -----------------------% +% +% Convective numerical method (JST, LAX-FRIEDRICH, CUSP, ROE, AUSM, HLLC, +% TURKEL_PREC, MSW) +CONV_NUM_METHOD_FLOW= JST +% +% 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 +% +% Time discretization (RUNGE-KUTTA_EXPLICIT, EULER_IMPLICIT, EULER_EXPLICIT) +TIME_DISCRE_FLOW= EULER_IMPLICIT + +% --------------------------- CONVERGENCE PARAMETERS --------------------------% +% +% Min value of the residual (log10 of the residual) +CONV_RESIDUAL_MINVAL= -19 +% +% Start convergence criteria at iteration number +CONV_STARTITER= 10 +% +% Number of elements to apply the criteria +CONV_CAUCHY_ELEMS= 100 +% +% Epsilon to control the series convergence +CONV_CAUCHY_EPS= 1E-6 + +% ------------------------- INPUT/OUTPUT INFORMATION --------------------------% +% +% Restart flow input file +SOLUTION_FILENAME= solution_flow.dat +% +% Restart adjoint input file +SOLUTION_ADJ_FILENAME= solution_adj.dat +% +% Output file format (TECPLOT, TECPLOT_BINARY, PARAVIEW, +% FIELDVIEW, FIELDVIEW_BINARY) +TABULAR_FORMAT= TECPLOT +% +% Output file convergence history (w/o extension) +CONV_FILENAME= history +% +% Output file with the forces breakdown +BREAKDOWN_FILENAME= forces_breakdown.dat +% +% Output file restart flow +RESTART_FILENAME= restart_flow.dat +% +% Output file restart adjoint +RESTART_ADJ_FILENAME= restart_adj.dat +% +% Output file flow (w/o extension) variables +VOLUME_FILENAME= flow +% +% Output file adjoint (w/o extension) variables +VOLUME_ADJ_FILENAME= adjoint +% +% Output Objective function +VALUE_OBJFUNC_FILENAME= of_eval.dat +% +% Output objective function gradient (using continuous adjoint) +GRAD_OBJFUNC_FILENAME= of_grad.dat +% +% Output file surface flow coefficient (w/o extension) +SURFACE_FILENAME= surface_flow +% +% Output file surface adjoint coefficient (w/o extension) +SURFACE_ADJ_FILENAME= surface_adjoint +% +% Writing solution file frequency +WRT_SOL_FREQ= 250 +% +% Writing convergence history frequency +WRT_CON_FREQ= 1 + +% ------------------------ GRID DEFORMATION PARAMETERS ------------------------% +% +% Linear solver or smoother for implicit formulations (FGMRES, RESTARTED_FGMRES, BCGSTAB) +DEFORM_LINEAR_SOLVER= FGMRES +% +% Number of smoothing iterations for mesh deformation +DEFORM_LINEAR_SOLVER_ITER= 200 +% +% Number of nonlinear deformation iterations (surface deformation increments) +DEFORM_NONLINEAR_ITER= 1 +% +% Print the residuals during mesh deformation to the console (YES, NO) +DEFORM_CONSOLE_OUTPUT= YES +% +% Type of element stiffness imposed for FEA mesh deformation (INVERSE_VOLUME, +% WALL_DISTANCE, CONSTANT_STIFFNESS) +DEFORM_STIFFNESS_TYPE= INVERSE_VOLUME +% +% Visualize the deformation (NO, YES) +VISUALIZE_VOLUME_DEF= YES diff --git a/TestCases/coupled_cht/comp_2d/solid_cylinder1.cfg b/TestCases/coupled_cht/comp_2d/solid_cylinder1.cfg new file mode 100644 index 000000000000..52bf3556c122 --- /dev/null +++ b/TestCases/coupled_cht/comp_2d/solid_cylinder1.cfg @@ -0,0 +1,186 @@ +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% % +% SU2 configuration file % +% Case description: Steady compressible laminar flow around heated cylinders % +% % +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +% ------------- DIRECT, ADJOINT, AND LINEARIZED PROBLEM DEFINITION ------------% +% +% Physical governing equations (EULER, NAVIER_STOKES, +% WAVE_EQUATION, HEAT_EQUATION, FEM_ELASTICITY, +% POISSON_EQUATION) +SOLVER= HEAT_EQUATION_FVM +% +% Restart solution (NO, YES) +RESTART_SOL= NO +% +% Objective function in gradient evaluation (DRAG, LIFT, SIDEFORCE, MOMENT_X, +% MOMENT_Y, MOMENT_Z, EFFICIENCY, +% EQUIVALENT_AREA, NEARFIELD_PRESSURE, +% FORCE_X, FORCE_Y, FORCE_Z, THRUST, +% TORQUE, TOTAL_HEATFLUX, +% MAXIMUM_HEATFLUX, INVERSE_DESIGN_PRESSURE, +% INVERSE_DESIGN_HEATFLUX, SURFACE_TOTAL_PRESSURE, +% SURFACE_MASSFLOW, SURFACE_STATIC_PRESSURE, SURFACE_MACH) +% For a weighted sum of objectives: separate by commas, add OBJECTIVE_WEIGHT and MARKER_MONITORING in matching order. +OBJECTIVE_FUNCTION= TOTAL_HEATFLUX +% +% List of weighting values when using more than one OBJECTIVE_FUNCTION. Separate by commas and match with MARKER_MONITORING. +OBJECTIVE_WEIGHT = 1.0 +% +% Read binary restart files (YES, NO) +READ_BINARY_RESTART = YES + +% -------------------- BOUNDARY CONDITION DEFINITION --------------------------% +% +MARKER_ISOTHERMAL= ( core1, 350.0 ) +% +% Marker(s) of the surface to be plotted or designed +MARKER_PLOTTING= (cylinder_inner1 ) +% +% Marker(s) of the surface where the functional (Cd, Cl, etc.) will be evaluated +MARKER_MONITORING= ( NONE ) + +% ---------------- (SOLIDS) CONDUCTION CONDITION DEFINITION -------------------% +% +% We should keep the dimensionalization of the coupled flow solver +INC_NONDIM= DIMENSIONAL +% +% Temperature initialization value +SOLID_TEMPERATURE_INIT= 350.0 +% +% Nettis case: hollow cylinder (air w/ 4x the conductivity) +% +% Solid density (kg/m^3) +SOLID_DENSITY= 0.000210322 +% +% Solid specific heat (J/kg*K) +SPECIFIC_HEAT_CP= 1004.703 +% +% Solid thermal conductivity (W/m*K) +SOLID_THERMAL_CONDUCTIVITY= 0.1028 + +% ------------- COMMON PARAMETERS DEFINING THE NUMERICAL METHOD ---------------% +% +% Numerical method for spatial gradients (GREEN_GAUSS, WEIGHTED_LEAST_SQUARES) +NUM_METHOD_GRAD= GREEN_GAUSS +% +% Courant-Friedrichs-Lewy condition of the finest grid +CFL_NUMBER= 10.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, 10.0, 10000.0 ) +% +% Runge-Kutta alpha coefficients +RK_ALPHA_COEFF= ( 0.66667, 0.66667, 1.000000 ) + +% ------------------------ LINEAR SOLVER DEFINITION ---------------------------% +% +% Linear solver or smoother for implicit formulations (BCGSTAB, FGMRES, SMOOTHER_JACOBI, +% SMOOTHER_ILU, SMOOTHER_LUSGS, +% SMOOTHER_LINELET) +LINEAR_SOLVER= FGMRES +% +% Preconditioner of the Krylov linear solver (ILU, LU_SGS, LINELET, JACOBI) +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-15 +% +% Max number of iterations of the linear solver for the implicit formulation +LINEAR_SOLVER_ITER= 5 + +% -------------------- HEAT NUMERICAL METHOD DEFINITION -----------------------% +% +TIME_DISCRE_HEAT= EULER_IMPLICIT + +% --------------------------- CONVERGENCE PARAMETERS --------------------------% +% +% Min value of the residual (log10 of the residual) +CONV_RESIDUAL_MINVAL= -19 +% +% Start convergence criteria at iteration number +CONV_STARTITER= 10 +% +% Number of elements to apply the criteria +CONV_CAUCHY_ELEMS= 100 +% +% Epsilon to control the series convergence +CONV_CAUCHY_EPS= 1E-6 + +% ------------------------- INPUT/OUTPUT INFORMATION --------------------------% +% +% Restart flow input file +SOLUTION_FILENAME= solution_flow.dat +% +% Restart adjoint input file +SOLUTION_ADJ_FILENAME= solution_adj.dat +% +% Output file format (TECPLOT, TECPLOT_BINARY, PARAVIEW, +% FIELDVIEW, FIELDVIEW_BINARY) +TABULAR_FORMAT= TECPLOT +% +% Output file convergence history (w/o extension) +CONV_FILENAME= history +% +% Output file with the forces breakdown +BREAKDOWN_FILENAME= forces_breakdown.dat +% +% Output file restart flow +RESTART_FILENAME= restart_flow.dat +% +% Output file restart adjoint +RESTART_ADJ_FILENAME= restart_adj.dat +% +% Output file flow (w/o extension) variables +VOLUME_FILENAME= flow +% +% Output file adjoint (w/o extension) variables +VOLUME_ADJ_FILENAME= adjoint +% +% Output Objective function +VALUE_OBJFUNC_FILENAME= of_eval.dat +% +% Output objective function gradient (using continuous adjoint) +GRAD_OBJFUNC_FILENAME= of_grad.dat +% +% Output file surface flow coefficient (w/o extension) +SURFACE_FILENAME= surface_flow +% +% Output file surface adjoint coefficient (w/o extension) +SURFACE_ADJ_FILENAME= surface_adjoint +% +% Writing solution file frequency +WRT_SOL_FREQ= 250 +% +% Writing convergence history frequency +WRT_CON_FREQ= 1 + +% ------------------------ GRID DEFORMATION PARAMETERS ------------------------% +% +% Linear solver or smoother for implicit formulations (FGMRES, RESTARTED_FGMRES, BCGSTAB) +DEFORM_LINEAR_SOLVER= FGMRES +% +% Number of smoothing iterations for mesh deformation +DEFORM_LINEAR_SOLVER_ITER= 200 +% +% Number of nonlinear deformation iterations (surface deformation increments) +DEFORM_NONLINEAR_ITER= 1 +% +% Print the residuals during mesh deformation to the console (YES, NO) +DEFORM_CONSOLE_OUTPUT= YES +% +% Type of element stiffness imposed for FEA mesh deformation (INVERSE_VOLUME, +% WALL_DISTANCE, CONSTANT_STIFFNESS) +DEFORM_STIFFNESS_TYPE= INVERSE_VOLUME +% +% Visualize the deformation (NO, YES) +VISUALIZE_VOLUME_DEF= YES diff --git a/TestCases/coupled_cht/comp_2d/solid_cylinder2.cfg b/TestCases/coupled_cht/comp_2d/solid_cylinder2.cfg new file mode 100644 index 000000000000..9ad5951c1be4 --- /dev/null +++ b/TestCases/coupled_cht/comp_2d/solid_cylinder2.cfg @@ -0,0 +1,196 @@ +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% % +% SU2 configuration file % +% Case description: Steady compressible laminar flow around heated cylinders % +% % +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +% ------------- DIRECT, ADJOINT, AND LINEARIZED PROBLEM DEFINITION ------------% +% +% Physical governing equations (EULER, NAVIER_STOKES, +% WAVE_EQUATION, HEAT_EQUATION, FEM_ELASTICITY, +% POISSON_EQUATION) +SOLVER= HEAT_EQUATION_FVM +% +% Restart solution (NO, YES) +RESTART_SOL= NO +% +% Objective function in gradient evaluation (DRAG, LIFT, SIDEFORCE, MOMENT_X, +% MOMENT_Y, MOMENT_Z, EFFICIENCY, +% EQUIVALENT_AREA, NEARFIELD_PRESSURE, +% FORCE_X, FORCE_Y, FORCE_Z, THRUST, +% TORQUE, TOTAL_HEATFLUX, +% MAXIMUM_HEATFLUX, INVERSE_DESIGN_PRESSURE, +% INVERSE_DESIGN_HEATFLUX, SURFACE_TOTAL_PRESSURE, +% SURFACE_MASSFLOW, SURFACE_STATIC_PRESSURE, SURFACE_MACH) +% For a weighted sum of objectives: separate by commas, add OBJECTIVE_WEIGHT and MARKER_MONITORING in matching order. +OBJECTIVE_FUNCTION= TOTAL_HEATFLUX +% +% List of weighting values when using more than one OBJECTIVE_FUNCTION. Separate by commas and match with MARKER_MONITORING. +OBJECTIVE_WEIGHT = 1.0 +% +% Read binary restart files (YES, NO) +READ_BINARY_RESTART = YES + +% -------------------- BOUNDARY CONDITION DEFINITION --------------------------% +% +MARKER_ISOTHERMAL= ( core2, 350.0 ) +% +% Marker(s) of the surface to be plotted or designed +MARKER_PLOTTING= (cylinder_inner2) +% +% Marker(s) of the surface where the functional (Cd, Cl, etc.) will be evaluated +MARKER_MONITORING= ( NONE ) + +% ---------------- (SOLIDS) CONDUCTION CONDITION DEFINITION -------------------% +% +% We should keep the dimensionalization of the coupled flow solver +INC_NONDIM= DIMENSIONAL +% +% Temperature initialization value +SOLID_TEMPERATURE_INIT= 350.0 +% +% Nettis case: hollow cylinder (air w/ 4x the conductivity) +% +% Solid density (kg/m^3) +SOLID_DENSITY= 0.000210322 +% +% Solid specific heat (J/kg*K) +SPECIFIC_HEAT_CP= 1004.703 +% +% Solid thermal conductivity (W/m*K) +SOLID_THERMAL_CONDUCTIVITY= 0.1028 + +% ------------- COMMON PARAMETERS DEFINING THE NUMERICAL METHOD ---------------% +% +% Numerical method for spatial gradients (GREEN_GAUSS, WEIGHTED_LEAST_SQUARES) +NUM_METHOD_GRAD= GREEN_GAUSS +% +% Courant-Friedrichs-Lewy condition of the finest grid +CFL_NUMBER= 10.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, 10.0, 10000.0 ) +% +% Runge-Kutta alpha coefficients +RK_ALPHA_COEFF= ( 0.66667, 0.66667, 1.000000 ) +% +% Objective function in gradient evaluation (DRAG, LIFT, SIDEFORCE, MOMENT_X, +% MOMENT_Y, MOMENT_Z, EFFICIENCY, +% EQUIVALENT_AREA, NEARFIELD_PRESSURE, +% FORCE_X, FORCE_Y, FORCE_Z, THRUST, +% TORQUE, TOTAL_HEATFLUX, +% MAXIMUM_HEATFLUX, INVERSE_DESIGN_PRESSURE, +% INVERSE_DESIGN_HEATFLUX, SURFACE_TOTAL_PRESSURE, +% SURFACE_MASSFLOW, SURFACE_STATIC_PRESSURE, SURFACE_MACH) + +% ------------------------ LINEAR SOLVER DEFINITION ---------------------------% +% +% Linear solver or smoother for implicit formulations (BCGSTAB, FGMRES, SMOOTHER_JACOBI, +% SMOOTHER_ILU, SMOOTHER_LUSGS, +% SMOOTHER_LINELET) +LINEAR_SOLVER= FGMRES +% +% Preconditioner of the Krylov linear solver (ILU, LU_SGS, LINELET, JACOBI) +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-15 +% +% Max number of iterations of the linear solver for the implicit formulation +LINEAR_SOLVER_ITER= 5 + +% -------------------- HEAT NUMERICAL METHOD DEFINITION -----------------------% +% +% Time discretization (EULER_IMPLICIT, EULER_EXPLICIT) +TIME_DISCRE_HEAT= EULER_IMPLICIT + +% --------------------------- CONVERGENCE PARAMETERS --------------------------% +% +% Min value of the residual (log10 of the residual) +CONV_RESIDUAL_MINVAL= -19 +% +% Start convergence criteria at iteration number +CONV_STARTITER= 10 +% +% Number of elements to apply the criteria +CONV_CAUCHY_ELEMS= 100 +% +% Epsilon to control the series convergence +CONV_CAUCHY_EPS= 1E-6 + +% ------------------------- INPUT/OUTPUT INFORMATION --------------------------% +% +% Restart flow input file +SOLUTION_FILENAME= solution_flow.dat +% +% Restart adjoint input file +SOLUTION_ADJ_FILENAME= solution_adj.dat +% +% Output file format (TECPLOT, TECPLOT_BINARY, PARAVIEW, +% FIELDVIEW, FIELDVIEW_BINARY) +TABULAR_FORMAT= TECPLOT +% +% Output file convergence history (w/o extension) +CONV_FILENAME= history +% +% Output file with the forces breakdown +BREAKDOWN_FILENAME= forces_breakdown.dat +% +% Output file restart flow +RESTART_FILENAME= restart_flow.dat +% +% Output file restart adjoint +RESTART_ADJ_FILENAME= restart_adj.dat +% +% Output file flow (w/o extension) variables +VOLUME_FILENAME= flow +% +% Output file adjoint (w/o extension) variables +VOLUME_ADJ_FILENAME= adjoint +% +% Output Objective function +VALUE_OBJFUNC_FILENAME= of_eval.dat +% +% Output objective function gradient (using continuous adjoint) +GRAD_OBJFUNC_FILENAME= of_grad.dat +% +% Output file surface flow coefficient (w/o extension) +SURFACE_FILENAME= surface_flow +% +% Output file surface adjoint coefficient (w/o extension) +SURFACE_ADJ_FILENAME= surface_adjoint +% +% Writing solution file frequency +WRT_SOL_FREQ= 250 +% +% Writing convergence history frequency +WRT_CON_FREQ= 1 + +% ------------------------ GRID DEFORMATION PARAMETERS ------------------------% +% +% Linear solver or smoother for implicit formulations (FGMRES, RESTARTED_FGMRES, BCGSTAB) +DEFORM_LINEAR_SOLVER= FGMRES +% +% Number of smoothing iterations for mesh deformation +DEFORM_LINEAR_SOLVER_ITER= 200 +% +% Number of nonlinear deformation iterations (surface deformation increments) +DEFORM_NONLINEAR_ITER= 1 +% +% Print the residuals during mesh deformation to the console (YES, NO) +DEFORM_CONSOLE_OUTPUT= YES +% +% Type of element stiffness imposed for FEA mesh deformation (INVERSE_VOLUME, +% WALL_DISTANCE, CONSTANT_STIFFNESS) +DEFORM_STIFFNESS_TYPE= INVERSE_VOLUME +% +% Visualize the deformation (NO, YES) +VISUALIZE_VOLUME_DEF= YES diff --git a/TestCases/coupled_cht/comp_2d/solid_cylinder3.cfg b/TestCases/coupled_cht/comp_2d/solid_cylinder3.cfg new file mode 100644 index 000000000000..9166a55eca8f --- /dev/null +++ b/TestCases/coupled_cht/comp_2d/solid_cylinder3.cfg @@ -0,0 +1,196 @@ +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% % +% SU2 configuration file % +% Case description: Steady compressible laminar flow around heated cylinders % +% % +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +% ------------- DIRECT, ADJOINT, AND LINEARIZED PROBLEM DEFINITION ------------% +% +% Physical governing equations (EULER, NAVIER_STOKES, +% WAVE_EQUATION, HEAT_EQUATION, FEM_ELASTICITY, +% POISSON_EQUATION) +SOLVER= HEAT_EQUATION_FVM +% +% Restart solution (NO, YES) +RESTART_SOL= NO +% +% Objective function in gradient evaluation (DRAG, LIFT, SIDEFORCE, MOMENT_X, +% MOMENT_Y, MOMENT_Z, EFFICIENCY, +% EQUIVALENT_AREA, NEARFIELD_PRESSURE, +% FORCE_X, FORCE_Y, FORCE_Z, THRUST, +% TORQUE, TOTAL_HEATFLUX, +% MAXIMUM_HEATFLUX, INVERSE_DESIGN_PRESSURE, +% INVERSE_DESIGN_HEATFLUX, SURFACE_TOTAL_PRESSURE, +% SURFACE_MASSFLOW, SURFACE_STATIC_PRESSURE, SURFACE_MACH) +% For a weighted sum of objectives: separate by commas, add OBJECTIVE_WEIGHT and MARKER_MONITORING in matching order. +OBJECTIVE_FUNCTION= TOTAL_HEATFLUX +% +% List of weighting values when using more than one OBJECTIVE_FUNCTION. Separate by commas and match with MARKER_MONITORING. +OBJECTIVE_WEIGHT = 1.0 +% +% Read binary restart files (YES, NO) +READ_BINARY_RESTART = YES + +% -------------------- BOUNDARY CONDITION DEFINITION --------------------------% +% +MARKER_ISOTHERMAL= ( core3, 350.0 ) +% +% Marker(s) of the surface to be plotted or designed +MARKER_PLOTTING= ( cylinder_inner3 ) +% +% Marker(s) of the surface where the functional (Cd, Cl, etc.) will be evaluated +MARKER_MONITORING= ( NONE ) + +% ---------------- (SOLIDS) CONDUCTION CONDITION DEFINITION -------------------% +% +% We should keep the dimensionalization of the coupled flow solver +INC_NONDIM= DIMENSIONAL +% +% Temperature initialization value +SOLID_TEMPERATURE_INIT= 350.0 +% +% Nettis case: hollow cylinder (air w/ 4x the conductivity) +% +% Solid density (kg/m^3) +SOLID_DENSITY= 0.000210322 +% +% Solid specific heat (J/kg*K) +SPECIFIC_HEAT_CP= 1004.703 +% +% Solid thermal conductivity (W/m*K) +SOLID_THERMAL_CONDUCTIVITY= 0.1028 + +% ------------- COMMON PARAMETERS DEFINING THE NUMERICAL METHOD ---------------% +% +% Numerical method for spatial gradients (GREEN_GAUSS, WEIGHTED_LEAST_SQUARES) +NUM_METHOD_GRAD= GREEN_GAUSS +% +% Courant-Friedrichs-Lewy condition of the finest grid +CFL_NUMBER= 10.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, 10.0, 10000.0 ) +% +% Runge-Kutta alpha coefficients +RK_ALPHA_COEFF= ( 0.66667, 0.66667, 1.000000 ) +% +% Objective function in gradient evaluation (DRAG, LIFT, SIDEFORCE, MOMENT_X, +% MOMENT_Y, MOMENT_Z, EFFICIENCY, +% EQUIVALENT_AREA, NEARFIELD_PRESSURE, +% FORCE_X, FORCE_Y, FORCE_Z, THRUST, +% TORQUE, TOTAL_HEATFLUX, +% MAXIMUM_HEATFLUX, INVERSE_DESIGN_PRESSURE, +% INVERSE_DESIGN_HEATFLUX, SURFACE_TOTAL_PRESSURE, +% SURFACE_MASSFLOW, SURFACE_STATIC_PRESSURE, SURFACE_MACH) + +% ------------------------ LINEAR SOLVER DEFINITION ---------------------------% +% +% Linear solver or smoother for implicit formulations (BCGSTAB, FGMRES, SMOOTHER_JACOBI, +% SMOOTHER_ILU, SMOOTHER_LUSGS, +% SMOOTHER_LINELET) +LINEAR_SOLVER= FGMRES +% +% Preconditioner of the Krylov linear solver (ILU, LU_SGS, LINELET, JACOBI) +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-15 +% +% Max number of iterations of the linear solver for the implicit formulation +LINEAR_SOLVER_ITER= 5 + +% -------------------- HEAT NUMERICAL METHOD DEFINITION -----------------------% +% +% Time discretization (EULER_IMPLICIT, EULER_EXPLICIT) +TIME_DISCRE_HEAT= EULER_IMPLICIT + +% --------------------------- CONVERGENCE PARAMETERS --------------------------% +% +% Min value of the residual (log10 of the residual) +CONV_RESIDUAL_MINVAL= -19 +% +% Start convergence criteria at iteration number +CONV_STARTITER= 10 +% +% Number of elements to apply the criteria +CONV_CAUCHY_ELEMS= 100 +% +% Epsilon to control the series convergence +CONV_CAUCHY_EPS= 1E-6 + +% ------------------------- INPUT/OUTPUT INFORMATION --------------------------% +% +% Restart flow input file +SOLUTION_FILENAME= solution_flow.dat +% +% Restart adjoint input file +SOLUTION_ADJ_FILENAME= solution_adj.dat +% +% Output file format (TECPLOT, TECPLOT_BINARY, PARAVIEW, +% FIELDVIEW, FIELDVIEW_BINARY) +TABULAR_FORMAT= TECPLOT +% +% Output file convergence history (w/o extension) +CONV_FILENAME= history +% +% Output file with the forces breakdown +BREAKDOWN_FILENAME= forces_breakdown.dat +% +% Output file restart flow +RESTART_FILENAME= restart_flow.dat +% +% Output file restart adjoint +RESTART_ADJ_FILENAME= restart_adj.dat +% +% Output file flow (w/o extension) variables +VOLUME_FILENAME= flow +% +% Output file adjoint (w/o extension) variables +VOLUME_ADJ_FILENAME= adjoint +% +% Output Objective function +VALUE_OBJFUNC_FILENAME= of_eval.dat +% +% Output objective function gradient (using continuous adjoint) +GRAD_OBJFUNC_FILENAME= of_grad.dat +% +% Output file surface flow coefficient (w/o extension) +SURFACE_FILENAME= surface_flow +% +% Output file surface adjoint coefficient (w/o extension) +SURFACE_ADJ_FILENAME= surface_adjoint +% +% Writing solution file frequency +WRT_SOL_FREQ= 250 +% +% Writing convergence history frequency +WRT_CON_FREQ= 1 + +% ------------------------ GRID DEFORMATION PARAMETERS ------------------------% +% +% Linear solver or smoother for implicit formulations (FGMRES, RESTARTED_FGMRES, BCGSTAB) +DEFORM_LINEAR_SOLVER= FGMRES +% +% Number of smoothing iterations for mesh deformation +DEFORM_LINEAR_SOLVER_ITER= 200 +% +% Number of nonlinear deformation iterations (surface deformation increments) +DEFORM_NONLINEAR_ITER= 1 +% +% Print the residuals during mesh deformation to the console (YES, NO) +DEFORM_CONSOLE_OUTPUT= YES +% +% Type of element stiffness imposed for FEA mesh deformation (INVERSE_VOLUME, +% WALL_DISTANCE, CONSTANT_STIFFNESS) +DEFORM_STIFFNESS_TYPE= INVERSE_VOLUME +% +% Visualize the deformation (NO, YES) +VISUALIZE_VOLUME_DEF= YES diff --git a/TestCases/coupled_cht/disc_adj_incomp_2d/flow_cylinder.cfg b/TestCases/coupled_cht/disc_adj_incomp_2d/flow_cylinder.cfg index c847420306e6..d9002e435bbb 100644 --- a/TestCases/coupled_cht/disc_adj_incomp_2d/flow_cylinder.cfg +++ b/TestCases/coupled_cht/disc_adj_incomp_2d/flow_cylinder.cfg @@ -27,6 +27,9 @@ RESTART_SOL= NO % INVERSE_DESIGN_HEATFLUX, SURFACE_TOTAL_PRESSURE, % SURFACE_MASSFLOW, SURFACE_STATIC_PRESSURE, SURFACE_MACH) % For a weighted sum of objectives: separate by commas, add OBJECTIVE_WEIGHT and MARKER_MONITORING in matching order. +% +% Gradient for this test case from direct differentiation (to confirm adjoints results, if desired): -0.716778 +% OBJECTIVE_FUNCTION= TOTAL_HEATFLUX % % List of weighting values when using more than one OBJECTIVE_FUNCTION. Separate by commas and match with MARKER_MONITORING. diff --git a/TestCases/coupled_cht/incomp_2d/cht_2d_3cylinders.cfg b/TestCases/coupled_cht/incomp_2d/cht_2d_3cylinders.cfg index 88f0428519c2..d0804e33899b 100644 --- a/TestCases/coupled_cht/incomp_2d/cht_2d_3cylinders.cfg +++ b/TestCases/coupled_cht/incomp_2d/cht_2d_3cylinders.cfg @@ -28,6 +28,9 @@ MARKER_ZONE_INTERFACE= (cylinder_outer1, cylinder_inner1, cylinder_outer2, cylin MARKER_CHT_INTERFACE= (cylinder_outer1, cylinder_inner1, cylinder_outer2, cylinder_inner2, cylinder_outer3, cylinder_inner3) % % +CHT_COUPLING_METHOD= DIRECT_TEMPERATURE_ROBIN_HEATFLUX +% +% TIME_DOMAIN = NO % % Number of total iterations diff --git a/TestCases/parallel_regression.py b/TestCases/parallel_regression.py index 8933726a21a8..2cfe4a9c3255 100644 --- a/TestCases/parallel_regression.py +++ b/TestCases/parallel_regression.py @@ -1084,9 +1084,9 @@ def main(): stat_fsi_restart.tol = 0.00001 test_list.append(stat_fsi_restart) - ########################## - ### Zonal multiphysics ### - ########################## + # ############################### + # ### Conjugate heat transfer ### + # ############################### # CHT incompressible cht_incompressible = TestCase('cht_incompressible') @@ -1100,6 +1100,18 @@ def main(): cht_incompressible.tol = 0.00001 test_list.append(cht_incompressible) + # CHT compressible + cht_compressible = TestCase('cht_compressible') + cht_compressible.cfg_dir = "coupled_cht/comp_2d" + cht_compressible.cfg_file = "cht_2d_3cylinders.cfg" + cht_compressible.test_iter = 10 + cht_compressible.test_vals = [-4.257607, -0.526125, -0.526125, -0.526125] #last 4 columns + cht_compressible.su2_exec = "SU2_CFD" + cht_compressible.timeout = 1600 + cht_compressible.multizone = True + cht_compressible.tol = 0.00001 + test_list.append(cht_compressible) + ########################## ### Python wrapper ### ########################## diff --git a/TestCases/serial_regression.py b/TestCases/serial_regression.py index 332cf3977aa7..3c28aa39894f 100644 --- a/TestCases/serial_regression.py +++ b/TestCases/serial_regression.py @@ -1179,9 +1179,9 @@ def main(): airfoilRBF.tol = 0.00001 test_list.append(airfoilRBF) - # ########################## - # ### Zonal multiphysics ### - # ########################## + # ############################### + # ### Conjugate heat transfer ### + # ############################### # CHT incompressible cht_incompressible = TestCase('cht_incompressible') @@ -1193,6 +1193,18 @@ def main(): cht_incompressible.timeout = 1600 cht_incompressible.multizone = True cht_incompressible.tol = 0.00001 + test_list.append(cht_incompressible) + + # CHT compressible + cht_incompressible = TestCase('cht_compressible') + cht_incompressible.cfg_dir = "coupled_cht/comp_2d" + cht_incompressible.cfg_file = "cht_2d_3cylinders.cfg" + cht_incompressible.test_iter = 10 + cht_incompressible.test_vals = [-4.257607, -0.526125, -0.526125, -0.526125] #last 4 columns + cht_incompressible.su2_exec = "SU2_CFD" + cht_incompressible.timeout = 1600 + cht_incompressible.multizone = True + cht_incompressible.tol = 0.00001 test_list.append(cht_incompressible) ##############################################