Skip to content

Updating SetResidual_DualTime needed for unsteady flows #2647

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
Cristopher-Morales merged 5 commits into from
Dec 17, 2025
+21 −33
Merged

Updating SetResidual_DualTime needed for unsteady flows #2647

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
d2586c0
initial commit
Cristopher-Morales Dec 16, 2025
18da912
Update SU2_CFD/src/solvers/CIncEulerSolver.cpp
Cristopher-Morales Dec 16, 2025
a322172
addressing comment
Dec 16, 2025
f788ed7
updating residuals
Dec 17, 2025
7975a6e
Merge branch 'develop' into feature_fix_V2U
Cristopher-Morales Dec 17, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
48 changes: 18 additions & 30 deletions SU2_CFD/src/solvers/CIncEulerSolver.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2799,18 +2799,17 @@ void CIncEulerSolver::SetResidual_DualTime(CGeometry *geometry, CSolver **solver
su2double U_time_nM1[MAXNVAR], U_time_n[MAXNVAR], U_time_nP1[MAXNVAR];
su2double Volume_nM1, Volume_nP1, TimeStep;
const su2double *Normal = nullptr, *GridVel_i = nullptr, *GridVel_j = nullptr;
su2double Density, Cp;
su2double Density;

const bool implicit = (config->GetKind_TimeIntScheme() == EULER_IMPLICIT);
const bool first_order = (config->GetTime_Marching() == TIME_MARCHING::DT_STEPPING_1ST);
const bool second_order = (config->GetTime_Marching() == TIME_MARCHING::DT_STEPPING_2ND);
const bool energy = config->GetEnergy_Equation();

const int ndim = nDim;
auto V2U = [ndim](su2double Density, su2double Cp, const su2double* V, su2double* U) {
const int nvar = nVar;
auto V2U = [nvar](su2double Density, const su2double* V, su2double* U) {
U[0] = Density;
for (int iDim = 0; iDim < ndim; iDim++) U[iDim+1] = Density*V[iDim+1];
U[ndim+1] = Density*Cp*V[ndim+1];
for (int iVar = 1; iVar < nvar; ++iVar) U[iVar] = Density * V[iVar];
};

/*--- Store the physical time step ---*/
Expand All @@ -2837,16 +2836,15 @@ void CIncEulerSolver::SetResidual_DualTime(CGeometry *geometry, CSolver **solver
V_time_n = nodes->GetSolution_time_n(iPoint);
V_time_nP1 = nodes->GetSolution(iPoint);

/*--- Access the density and Cp at this node (constant for now). ---*/
/*--- Access the density at this node (constant for now). ---*/

Density = nodes->GetDensity(iPoint);
Cp = nodes->GetSpecificHeatCp(iPoint);

/*--- Compute the conservative variable vector for all time levels. ---*/

V2U(Density, Cp, V_time_nM1, U_time_nM1);
V2U(Density, Cp, V_time_n, U_time_n);
V2U(Density, Cp, V_time_nP1, U_time_nP1);
V2U(Density, V_time_nM1, U_time_nM1);
V2U(Density, V_time_n, U_time_n);
V2U(Density, V_time_nP1, U_time_nP1);

/*--- CV volume at time n+1. As we are on a static mesh, the volume
of the CV will remained fixed for all time steps. ---*/
Expand All @@ -2867,13 +2865,9 @@ void CIncEulerSolver::SetResidual_DualTime(CGeometry *geometry, CSolver **solver
/*--- Compute the Jacobian contribution due to the dual time source term. ---*/

if (implicit) {
su2double delta = (second_order? 1.5 : 1.0) * Volume_nP1 * Density / TimeStep;

for (iDim = 0; iDim < nDim; iDim++)
Jacobian.AddVal2Diag(iPoint, iDim+1, delta);
su2double delta = (second_order ? 1.5 : 1.0) * Volume_nP1 * Density / TimeStep;

if (energy) delta *= Cp;
Jacobian.AddVal2Diag(iPoint, nDim+1, delta);
for (iVar = 1; iVar < nVar; ++iVar) Jacobian.AddVal2Diag(iPoint, iVar, delta);
}
}
END_SU2_OMP_FOR
Expand All @@ -2898,8 +2892,7 @@ void CIncEulerSolver::SetResidual_DualTime(CGeometry *geometry, CSolver **solver

V_time_n = nodes->GetSolution_time_n(iPoint);
Density = nodes->GetDensity(iPoint);
Cp = nodes->GetSpecificHeatCp(iPoint);
V2U(Density, Cp, V_time_n, U_time_n);
V2U(Density, V_time_n, U_time_n);

GridVel_i = geometry->nodes->GetGridVel(iPoint);

Expand Down Expand Up @@ -2954,8 +2947,7 @@ void CIncEulerSolver::SetResidual_DualTime(CGeometry *geometry, CSolver **solver

V_time_n = nodes->GetSolution_time_n(iPoint);
Density = nodes->GetDensity(iPoint);
Cp = nodes->GetSpecificHeatCp(iPoint);
V2U(Density, Cp, V_time_n, U_time_n);
V2U(Density, V_time_n, U_time_n);

for (iVar = 0; iVar < nVar-!energy; iVar++)
LinSysRes(iPoint,iVar) += U_time_n[iVar]*Residual_GCL;
Expand All @@ -2981,16 +2973,15 @@ void CIncEulerSolver::SetResidual_DualTime(CGeometry *geometry, CSolver **solver
V_time_n = nodes->GetSolution_time_n(iPoint);
V_time_nP1 = nodes->GetSolution(iPoint);

/*--- Access the density and Cp at this node (constant for now). ---*/
/*--- Access the density at this node (constant for now). ---*/

Density = nodes->GetDensity(iPoint);
Cp = nodes->GetSpecificHeatCp(iPoint);

/*--- Compute the conservative variable vector for all time levels. ---*/

V2U(Density, Cp, V_time_nM1, U_time_nM1);
V2U(Density, Cp, V_time_n, U_time_n);
V2U(Density, Cp, V_time_nP1, U_time_nP1);
V2U(Density, V_time_nM1, U_time_nM1);
V2U(Density, V_time_n, U_time_n);
V2U(Density, V_time_nP1, U_time_nP1);

/*--- CV volume at time n-1 and n+1. In the case of dynamically deforming
grids, the volumes will change. On rigidly transforming grids, the
Expand All @@ -3016,11 +3007,8 @@ void CIncEulerSolver::SetResidual_DualTime(CGeometry *geometry, CSolver **solver
if (implicit) {
su2double delta = (second_order? 1.5 : 1.0) * Volume_nP1 * Density / TimeStep;

for (iDim = 0; iDim < nDim; iDim++)
Jacobian.AddVal2Diag(iPoint, iDim+1, delta);

if (energy) delta *= Cp;
Jacobian.AddVal2Diag(iPoint, nDim+1, delta);
for (iVar = 1; iVar < nVar; ++iVar)
Jacobian.AddVal2Diag(iPoint, iVar, delta);
}
}
END_SU2_OMP_FOR
Expand Down
4 changes: 2 additions & 2 deletions TestCases/parallel_regression_AD.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -316,8 +316,8 @@ def main():
da_unsteadyCHT_cylinder.cfg_dir = "coupled_cht/disc_adj_unsteadyCHT_cylinder"
da_unsteadyCHT_cylinder.cfg_file = "chtMaster.cfg"
da_unsteadyCHT_cylinder.test_iter = 2
da_unsteadyCHT_cylinder.test_vals = [-12.131633, -12.617906, -12.688933, -16.179747, -6.432277, 0.000000, 75.761000, 0.247780]
da_unsteadyCHT_cylinder.test_vals_aarch64 = [-12.131633, -12.617906, -12.688933, -16.179747, -6.432277, 0.000000, 75.761000, 0.247780]
da_unsteadyCHT_cylinder.test_vals = [-8.479629, -9.239920, -9.234868, -15.934511, -13.662012, 0.000000, 89.932000, 0.295190]
da_unsteadyCHT_cylinder.test_vals_aarch64 = [-8.479629, -9.239920, -9.234868, -15.934511, -13.662012, 0.000000, 89.932000, 0.295190]
da_unsteadyCHT_cylinder.unsteady = True
da_unsteadyCHT_cylinder.multizone = True
test_list.append(da_unsteadyCHT_cylinder)
Expand Down
2 changes: 1 addition & 1 deletion TestCases/tutorials.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -49,7 +49,7 @@ def main():
cht_incompressible_unsteady.cfg_dir = "../Tutorials/multiphysics/unsteady_cht/"
cht_incompressible_unsteady.cfg_file = "cht_2d_3cylinders.cfg"
cht_incompressible_unsteady.test_iter = 2
cht_incompressible_unsteady.test_vals = [-2.661443, -2.546289, -0.080399, -0.080399, -0.080399, -12.421980, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 238.240000] #last columns
cht_incompressible_unsteady.test_vals = [-2.661440, -2.534489, -0.080399, -0.080399, -0.080399, -12.421979, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 238.240000] #last columns
cht_incompressible_unsteady.multizone = True
cht_incompressible_unsteady.unsteady = True
test_list.append(cht_incompressible_unsteady)
Expand Down