@@ -79,7 +79,8 @@ CScalarSolver<VariableType>::~CScalarSolver() {
7979
8080template <class VariableType >
8181void CScalarSolver<VariableType>::Upwind_Residual(CGeometry* geometry, CSolver** solver_container,
82- CNumerics** numerics_container, CConfig* config, unsigned short iMesh) {
82+ CNumerics** numerics_container, CConfig* config,
83+ unsigned short iMesh) {
8384 const bool implicit = (config->GetKind_TimeIntScheme () == EULER_IMPLICIT);
8485 const bool muscl = config->GetMUSCL_Turb ();
8586 const bool limiter = (config->GetKind_SlopeLimit_Turb () != NO_LIMITER);
@@ -256,7 +257,7 @@ void CScalarSolver<VariableType>::SumEdgeFluxes(CGeometry* geometry) {
256257
257258template <class VariableType >
258259void CScalarSolver<VariableType>::BC_Periodic(CGeometry* geometry, CSolver** solver_container, CNumerics* numerics,
259- CConfig* config) {
260+ CConfig* config) {
260261 /* --- Complete residuals for periodic boundary conditions. We loop over
261262 the periodic BCs in matching pairs so that, in the event that there are
262263 adjacent periodic markers, the repeated points will have their residuals
@@ -271,7 +272,7 @@ void CScalarSolver<VariableType>::BC_Periodic(CGeometry* geometry, CSolver** sol
271272
272273template <class VariableType >
273274void CScalarSolver<VariableType>::PrepareImplicitIteration(CGeometry* geometry, CSolver** solver_container,
274- CConfig* config) {
275+ CConfig* config) {
275276 const auto flowNodes = solver_container[FLOW_SOL]->GetNodes ();
276277
277278 /* --- Set shared residual variables to 0 and declare
@@ -280,7 +281,6 @@ void CScalarSolver<VariableType>::PrepareImplicitIteration(CGeometry* geometry,
280281 SetResToZero ();
281282
282283 su2double resMax[MAXNVAR] = {0.0 }, resRMS[MAXNVAR] = {0.0 };
283- const su2double* coordMax[MAXNVAR] = {nullptr };
284284 unsigned long idxMax[MAXNVAR] = {0 };
285285
286286 /* --- Build implicit system ---*/
@@ -308,31 +308,19 @@ void CScalarSolver<VariableType>::PrepareImplicitIteration(CGeometry* geometry,
308308 LinSysRes[total_index] = -LinSysRes[total_index];
309309 LinSysSol[total_index] = 0.0 ;
310310
311- su2double Res = fabs (LinSysRes[total_index]);
312- resRMS[iVar] += Res * Res;
313- if (Res > resMax[iVar]) {
314- resMax[iVar] = Res;
315- idxMax[iVar] = iPoint;
316- coordMax[iVar] = geometry->nodes ->GetCoord (iPoint);
317- }
311+ /* --- "Add" residual at (iPoint,iVar) to local residual variables. ---*/
312+ ResidualReductions_PerThread (iPoint, iVar, LinSysRes[total_index], resRMS, resMax, idxMax);
318313 }
319314 }
320315 END_SU2_OMP_FOR
321- SU2_OMP_CRITICAL
322- for (unsigned short iVar = 0 ; iVar < nVar; iVar++) {
323- Residual_RMS[iVar] += resRMS[iVar];
324- AddRes_Max (iVar, resMax[iVar], geometry->nodes ->GetGlobalIndex (idxMax[iVar]), coordMax[iVar]);
325- }
326- END_SU2_OMP_CRITICAL
327- SU2_OMP_BARRIER
328316
329- /* --- Compute the root mean square residual ---*/
330- SetResidual_RMS (geometry, config);
317+ /* --- "Add" residuals from all threads to global residual variables. ---*/
318+ ResidualReductions_FromAllThreads (geometry, config, resRMS, resMax, idxMax );
331319}
332320
333321template <class VariableType >
334322void CScalarSolver<VariableType>::CompleteImplicitIteration(CGeometry* geometry, CSolver** solver_container,
335- CConfig* config) {
323+ CConfig* config) {
336324 const bool compressible = (config->GetKind_Regime () == ENUM_REGIME::COMPRESSIBLE);
337325
338326 const auto flowNodes = solver_container[FLOW_SOL]->GetNodes ();
@@ -380,7 +368,7 @@ void CScalarSolver<VariableType>::CompleteImplicitIteration(CGeometry* geometry,
380368
381369template <class VariableType >
382370void CScalarSolver<VariableType>::ImplicitEuler_Iteration(CGeometry* geometry, CSolver** solver_container,
383- CConfig* config) {
371+ CConfig* config) {
384372 PrepareImplicitIteration (geometry, solver_container, config);
385373
386374 /* --- Solve or smooth the linear system. ---*/
@@ -404,10 +392,41 @@ void CScalarSolver<VariableType>::ImplicitEuler_Iteration(CGeometry* geometry, C
404392 CompleteImplicitIteration (geometry, solver_container, config);
405393}
406394
395+ template <class VariableType >
396+ void CScalarSolver<VariableType>::ExplicitEuler_Iteration(CGeometry* geometry, CSolver** solver_container,
397+ CConfig* config) {
398+ const auto flowNodes = solver_container[FLOW_SOL]->GetNodes ();
399+
400+ /* --- Local residual variables for current thread ---*/
401+ su2double resMax[MAXNVAR] = {0.0 }, resRMS[MAXNVAR] = {0.0 };
402+ unsigned long idxMax[MAXNVAR] = {0 };
403+
404+ SU2_OMP_FOR_STAT (omp_chunk_size)
405+ for (unsigned long iPoint = 0 ; iPoint < nPointDomain; iPoint++) {
406+ const su2double dt = nodes->GetLocalCFL (iPoint) / flowNodes->GetLocalCFL (iPoint) * flowNodes->GetDelta_Time (iPoint);
407+ nodes->SetDelta_Time (iPoint, dt);
408+ const su2double Vol = geometry->nodes ->GetVolume (iPoint) + geometry->nodes ->GetPeriodicVolume (iPoint);
409+
410+ for (auto iVar = 0u ; iVar < nVar; iVar++) {
411+ /* --- "Add" residual at (iPoint,iVar) to local residual variables. ---*/
412+ ResidualReductions_PerThread (iPoint, iVar, LinSysRes (iPoint, iVar), resRMS, resMax, idxMax);
413+ /* --- Explicit Euler step: ---*/
414+ LinSysSol (iPoint, iVar) = -dt / Vol * LinSysRes (iPoint, iVar);
415+ }
416+ }
417+ END_SU2_OMP_FOR
418+
419+ /* --- "Add" residuals from all threads to global residual variables. ---*/
420+ ResidualReductions_FromAllThreads (geometry, config, resRMS, resMax, idxMax);
421+
422+ /* --- Use LinSysSol for solution update. ---*/
423+ CompleteImplicitIteration (geometry, solver_container, config);
424+ }
425+
407426template <class VariableType >
408427void CScalarSolver<VariableType>::SetResidual_DualTime(CGeometry* geometry, CSolver** solver_container, CConfig* config,
409- unsigned short iRKStep, unsigned short iMesh,
410- unsigned short RunTime_EqSystem) {
428+ unsigned short iRKStep, unsigned short iMesh,
429+ unsigned short RunTime_EqSystem) {
411430 const bool implicit = (config->GetKind_TimeIntScheme () == EULER_IMPLICIT);
412431 const bool first_order = (config->GetTime_Marching () == TIME_MARCHING::DT_STEPPING_1ST);
413432 const bool second_order = (config->GetTime_Marching () == TIME_MARCHING::DT_STEPPING_2ND);
0 commit comments