Skip to content

Commit 10cc4dc

Browse files
committed
evaluate and accumulate each flux separately
1 parent 8a5ca6c commit 10cc4dc

File tree

1 file changed

+17
-24
lines changed

1 file changed

+17
-24
lines changed

SU2_CFD/src/solvers/CTurbSASolver.cpp

Lines changed: 17 additions & 24 deletions
Original file line numberDiff line numberDiff line change
@@ -305,9 +305,9 @@ void CTurbSASolver::Viscous_Residual(const unsigned long iEdge, const CGeometry*
305305
CFlowVariable* flowNodes = solver_container[FLOW_SOL] ?
306306
su2staticcast_p<CFlowVariable*>(solver_container[FLOW_SOL]->GetNodes()) : nullptr;
307307

308-
/*--- Points in edge ---*/
309-
auto iPoint = geometry->edges->GetNode(iEdge, 0);
310-
auto jPoint = geometry->edges->GetNode(iEdge, 1);
308+
/*--- Points in edge ---*/
309+
auto iPoint = geometry->edges->GetNode(iEdge, 0);
310+
auto jPoint = geometry->edges->GetNode(iEdge, 1);
311311

312312
/*--- Helper function to compute the flux ---*/
313313
auto ComputeFlux = [&](unsigned long point_i, unsigned long point_j, const su2double* normal) {
@@ -323,33 +323,26 @@ void CTurbSASolver::Viscous_Residual(const unsigned long iEdge, const CGeometry*
323323

324324
return numerics->ComputeResidual(config);
325325
};
326-
/*--- Compute fluxes at each node ---*/
327-
const su2double* normal = geometry->edges->GetNormal(iEdge);
326+
327+
SolverSpecificNumerics(iPoint, jPoint);
328328

329+
/*--- Compute fluxes and jacobians i->j ---*/
330+
const su2double* normal = geometry->edges->GetNormal(iEdge);
331+
auto residual_ij = ComputeFlux(iPoint, jPoint, normal);
332+
LinSysRes.SubtractBlock(iPoint, residual_ij);
333+
if (implicit) {
334+
Jacobian.UpdateBlocksSub(iEdge, iPoint, jPoint, residual_ij.jacobian_i, residual_ij.jacobian_j);
335+
}
336+
337+
/*--- Compute fluxes and jacobians j->i ---*/
329338
su2double flipped_normal[3];
330339
for (int iDim=0; iDim<nDim; iDim++)
331340
flipped_normal[iDim] = -normal[iDim];
332341

333-
SolverSpecificNumerics(iPoint, jPoint);
334-
335-
auto residual_ij = ComputeFlux(iPoint, jPoint, normal);
336342
auto residual_ji = ComputeFlux(jPoint, iPoint, flipped_normal);
337-
338-
if (ReducerStrategy) {
339-
EdgeFluxes.SubtractBlock(iEdge, residual_ij);
340-
EdgeFluxes.AddBlock(iEdge, residual_ji);
341-
if (implicit) {
342-
Jacobian.UpdateBlocksSub(iEdge, residual_ij.jacobian_i, residual_ij.jacobian_j);
343-
Jacobian.UpdateBlocks(iEdge, residual_ji.jacobian_i, residual_ji.jacobian_j);
344-
}
345-
}
346-
else {
347-
LinSysRes.SubtractBlock(iPoint, residual_ij);
348-
LinSysRes.AddBlock(jPoint, residual_ji);
349-
if (implicit) {
350-
Jacobian.UpdateBlocksSub(iEdge, iPoint, jPoint, residual_ij.jacobian_i, residual_ij.jacobian_j);
351-
Jacobian.UpdateBlocks(iEdge, iPoint, jPoint, residual_ji.jacobian_i, residual_ji.jacobian_j);
352-
}
343+
LinSysRes.AddBlock(jPoint, residual_ji);
344+
if (implicit) {
345+
Jacobian.UpdateBlocks(iEdge, iPoint, jPoint, residual_ij.jacobian_i, residual_ij.jacobian_j);
353346
}
354347
}
355348

0 commit comments

Comments
 (0)