FDS Source: Issue #14221. Eliminate OpenMP race conditions

Source/divg.f90

@@ -185,6 +185,7 @@ SUBROUTINE DIVERGENCE_PART_1(T,DT,NM)
          WC => WALL(IW)
          IF (WC%BOUNDARY_TYPE==NULL_BOUNDARY) CYCLE WALL_LOOP
          BC => BOUNDARY_COORD(WC%BC_INDEX)
+         IF (WC%THIN .AND. BC%IOR<0) CYCLE WALL_LOOP  ! Avoid OpenMP race condition by processing on one side of thin OBST
          BOUNDARY_TYPE_SELECT: SELECT CASE(WC%BOUNDARY_TYPE)
             CASE DEFAULT
                SELECT CASE(BC%IOR)
@@ -322,6 +323,8 @@ SUBROUTINE DIVERGENCE_PART_1(T,DT,NM)
          ENDIF
          B1%RHO_D_DZDN_F(N) = RHO_D_DZDN
 
+         IF (WC%THIN .AND. BC%IOR<0) CYCLE WALL_LOOP_2  ! Avoid OpenMP race condition by processing only one side of thin OBST
+
          IF (STORE_SPECIES_FLUX) THEN
             IF (CORRECTOR) THEN
                SELECT CASE(BC%IOR)
@@ -529,6 +532,9 @@ SUBROUTINE DIVERGENCE_PART_1(T,DT,NM)
    ELSE
       B1%K_G = KP(BC%IIG,BC%JJG,BC%KKG)
    ENDIF
+   ! Q_LEAK accounts for enthalpy moving through leakage paths
+   DP(BC%IIG,BC%JJG,BC%KKG) = DP(BC%IIG,BC%JJG,BC%KKG) - ( B1%AREA_ADJUST*B1%Q_CON_F*B1%RDN - B1%Q_LEAK )
+   IF (WC%THIN .AND. BC%IOR<0) CYCLE CORRECTION_LOOP  ! Avoid OpenMP race condition by processing on one side of thin OBST
    SELECT CASE(BC%IOR)
       CASE( 1) ; KDTDX(BC%II  ,BC%JJ  ,BC%KK  ) = 0._EB
       CASE(-1) ; KDTDX(BC%II-1,BC%JJ  ,BC%KK  ) = 0._EB
@@ -537,8 +543,6 @@ SUBROUTINE DIVERGENCE_PART_1(T,DT,NM)
       CASE( 3) ; KDTDZ(BC%II  ,BC%JJ  ,BC%KK  ) = 0._EB
       CASE(-3) ; KDTDZ(BC%II  ,BC%JJ  ,BC%KK-1) = 0._EB
    END SELECT
-   ! Q_LEAK accounts for enthalpy moving through leakage paths
-   DP(BC%IIG,BC%JJG,BC%KKG) = DP(BC%IIG,BC%JJG,BC%KKG) - ( B1%AREA_ADJUST*B1%Q_CON_F*B1%RDN - B1%Q_LEAK )
 ENDDO CORRECTION_LOOP
 
 ! Compute (q + del dot k del T) and add to the divergence

Source/func.f90

@@ -1399,6 +1399,7 @@ SUBROUTINE GET_SCALAR_FACE_VALUE(A,U,F,I1,I2,J1,J2,K1,K2,IOR,LIMITER)
       ENDDO
       !$OMP END DO
    CASE(5) ! MP5, Suresh and Huynh (1997)
+      !$OMP SINGLE
       DO K=K1,K2
          DO J=J1,J2
             DO I=I1,I2
@@ -1414,6 +1415,7 @@ SUBROUTINE GET_SCALAR_FACE_VALUE(A,U,F,I1,I2,J1,J2,K1,K2,IOR,LIMITER)
             ENDDO
          ENDDO
       ENDDO
+      !$OMP END SINGLE
 END SELECT
 !$OMP END PARALLEL
 
@@ -3803,6 +3805,8 @@ SUBROUTINE PACK_WALL(NM,OS,WC,SURF_INDEX,RC,IC,LC,UNPACK_IT,COUNT_ONLY,CHECK_BOU
 IC=IC+1 ; IF (.NOT.COUNT_ONLY) CALL EQUATE(OS%INTEGERS(IC),WC%N_INTEGERS,UNPACK_IT)
 IC=IC+1 ; IF (.NOT.COUNT_ONLY) CALL EQUATE(OS%INTEGERS(IC),WC%N_LOGICALS,UNPACK_IT)
 
+LC=LC+1 ; IF (.NOT.COUNT_ONLY) CALL EQUATE(OS%LOGICALS(LC),WC%THIN,UNPACK_IT)
+
 ! Pack or unpack the appropriate derived type variables tied to this wall cell
 
 IF (SF%INCLUDE_BOUNDARY_COORD_TYPE) CALL PACK_BOUNDARY_COORD(NM,IC,RC,OS,BC_INDEX,UNPACK_IT,COUNT_ONLY)

Source/init.f90

@@ -3251,6 +3251,8 @@ SUBROUTINE INIT_WALL_CELL(NM,I,J,K,OBST_INDEX,IW,IOR,SURF_INDEX,IERR,TT)
 WC%OBST_INDEX = OBST_INDEX
 WC%BOUNDARY_TYPE = BOUNDARY_TYPE
 
+IF (.NOT.M%CELL(ICG)%SOLID .AND. .NOT.M%CELL(IC)%SOLID .AND. .NOT.M%CELL(IC)%EXTERIOR) WC%THIN = .TRUE.
+
 IF (IW<=M%N_EXTERNAL_WALL_CELLS) THEN
    EWC => M%EXTERNAL_WALL(IW)
    EWC%NOM     = NOM_FOUND

Source/type.f90

@@ -442,6 +442,8 @@ MODULE TYPES
    INTEGER :: N_INTEGERS=0            !< Number of integers to pack into restart or send/recv buffer
    INTEGER :: N_LOGICALS=0            !< Number of logicals to pack into restart or send/recv buffer
 
+   LOGICAL :: THIN=.FALSE.            !< Indicates if the underlying solid is zero cells thick
+
 END TYPE WALL_TYPE
 
 

Source/wall.f90

@@ -333,9 +333,11 @@ END SUBROUTINE ASSIGN_GHOST_VALUE
 !> \param SF Pointer to SURFACE derived type
 !> \param BC Pointer to BOUNDARY_COORD derived type
 !> \param B1 Pointer to BOUNDARY_PROP1 derived type
-!> \param LP Pointer to LAGRANGIAN_PARTICLE derived type
-!> \param WALL_INDEX Index of wall cell
-!> \param PARTICLE_INDEX Index of particle
+!> \param LP (Optional) Pointer to LAGRANGIAN_PARTICLE derived type
+!> \param TW (Optional) Pointer to THIN_WALL
+!> \param WALL_INDEX (Optional) Index of wall cell
+!> \param PARTICLE_INDEX (Optional) Index of particle
+!> \param THIN_WALL_INDEX (Optional) Index of thin wall cell
 
 SUBROUTINE NEAR_SURFACE_GAS_VARIABLES(T,SF,BC,B1,LP,TW,WALL_INDEX,PARTICLE_INDEX,THIN_WALL_INDEX)
 
@@ -1384,6 +1386,7 @@ SUBROUTINE CALCULATE_ZZ_F(T,DT,WALL_INDEX,CFACE_INDEX,PARTICLE_INDEX)
             ENDIF
          ENDIF
          IF (OTHER_MESH_OBST_INDEX>0) THEN
+            !$OMP CRITICAL
             IF (OBST_INDEX>0) OBSTRUCTION(OBST_INDEX)%MASS = MESHES(EWC%NOM)%OBSTRUCTION(OTHER_MESH_OBST_INDEX)%MASS
             IF (MESHES(EWC%NOM)%OBSTRUCTION(OTHER_MESH_OBST_INDEX)%CONSUMABLE) THEN
                OMESH(EWC%NOM)%N_EXTERNAL_OBST = OMESH(EWC%NOM)%N_EXTERNAL_OBST + 1
@@ -1392,6 +1395,7 @@ SUBROUTINE CALCULATE_ZZ_F(T,DT,WALL_INDEX,CFACE_INDEX,PARTICLE_INDEX)
                OMESH(EWC%NOM)%REAL_SEND_PKG8(LL)   = &
                   (B1%M_DOT_PART_ACTUAL+SUM(B1%M_DOT_G_PP_ACTUAL(1:N_TRACKED_SPECIES)))*DT*B1%AREA
             ENDIF
+            !$OMP END CRITICAL
          ELSE
             !$OMP CRITICAL
             IF (OBST_INDEX>0) OBSTRUCTION(OBST_INDEX)%MASS = OBSTRUCTION(OBST_INDEX)%MASS - &
@@ -3482,7 +3486,7 @@ END SUBROUTINE PYROLYSIS
 !> \param WALL_INDEX_IN Optional wall cell index
 !> \param CFACE_INDEX_IN Optional cface index
 !> \param PARTICLE_INDEX_IN Optional particle index
-!> \param BACK_SIZE Optional flag indicating if the surface is on the back side of the obstruction
+!> \param BACK_SIDE Optional flag indicating if the surface is on the back side of the obstruction
 
 REAL(EB) FUNCTION HEAT_TRANSFER_COEFFICIENT(NM,T,DELTA_N_TMP,SF,WALL_INDEX_IN,CFACE_INDEX_IN,PARTICLE_INDEX_IN,BACK_SIDE)
 
@@ -3740,7 +3744,7 @@ SUBROUTINE HT3D_TEMPERATURE_EXCHANGE(NM)
 ENDDO WALL_LOOP
 !$OMP END DO
 
-!$OMP DO SCHEDULE(GUIDED) PRIVATE(ITW,TW,SF,BC,ONE_D,NWP,THR_D,I,II,IWA,NM2,I_NODE,WC2,BC2,ONE_D2,TW2)
+!$OMP DO SCHEDULE(GUIDED) PRIVATE(ITW,TW,SF,BC,ONE_D,NWP,THR_D,I,II,IWA,NM2,I_NODE,WC2,BC2,ONE_D2,TW2,TMP_1,TMP_NWP)
 THIN_WALL_LOOP: DO ITW=1,M%N_THIN_WALL_CELLS
 
    TW => M%THIN_WALL(ITW)