Skip to content

FDS Source: fix issue with zone breaking in hvac (Issue #16034) #16046

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
drjfloyd merged 1 commit into from
Mar 13, 2026
Merged

FDS Source: fix issue with zone breaking in hvac (Issue #16034) #16046

Changes from all commits
Commits
Show all changes
1 commit
Select commit
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
127 changes: 90 additions & 37 deletions Source/hvac.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1444,8 +1444,10 @@ SUBROUTINE HVAC_CALC(T,DT,FIRST_PASS)
ENDDO
CALL ADJUST_LEAKAGE_AREA
CALL FIND_NETWORKS(CHANGE,T) ! calls determined fixed elements (which calls update fan for fixed fans)
IF (N_NETWORKS == 0) RETURN
IF (HVAC_QFAN .OR. HVAC_MASS_TRANSPORT) CALL FIND_DUCTRUNS(CHANGE)
ELSE FIRST_PASS_IF !Not FIRST_PASS, reset variables to old solution
IF (N_NETWORKS == 0) RETURN
WHERE (.NOT. DUCT%FIXED)
DUCT(:)%VEL(NEW) = DUCT(:)%VEL(OLD)
DUCT(:)%VEL(GUESS) = DUCT(:)%VEL(OLD)
Expand Down Expand Up @@ -1583,7 +1585,7 @@ END SUBROUTINE MATRIX_SOLVE

SUBROUTINE HVAC_UPDATE(NNE,DT)

USE PHYSICAL_FUNCTIONS, ONLY : GET_SPECIFIC_GAS_CONSTANT,GET_ENTHALPY,GET_TEMPERATURE
USE PHYSICAL_FUNCTIONS, ONLY: GET_SPECIFIC_GAS_CONSTANT,GET_ENTHALPY,GET_TEMPERATURE
USE GLOBAL_CONSTANTS, ONLY: N_DUCTS,N_DUCTNODES,N_DUCTRUNS,N_TRACKED_SPECIES,HVAC_MASS_TRANSPORT,I_MAX_TEMP,CPBAR_Z,P_INF,NEW
REAL(EB),INTENT(IN) :: DT
INTEGER, INTENT(IN) :: NNE
Expand Down Expand Up @@ -1855,6 +1857,7 @@ SUBROUTINE HVAC_UPDATE(NNE,DT)
ENDIF
CYCLE DUCT_LOOP_2
ENDIF

IF (N_UPDATED(NN)) THEN
DU%RHO_D = DN%RHO
DU%TMP_D = DN%TMP
Expand Down Expand Up @@ -2539,8 +2542,8 @@ END SUBROUTINE HVAC_BC_IN

SUBROUTINE DETERMINE_FIXED_ELEMENTS(T,CHANGE)

USE MATH_FUNCTIONS,ONLY:EVALUATE_RAMP
USE GLOBAL_CONSTANTS, ONLY: N_DUCTS,NEW,N_DUCTNODES
USE MATH_FUNCTIONS,ONLY: EVALUATE_RAMP
USE GLOBAL_CONSTANTS, ONLY: N_DUCTS,NEW,N_DUCTNODES,CONNECTED_ZONES
INTEGER:: NN,ND, COUNTER
REAL(EB), INTENT(IN):: T
LOGICAL, INTENT(INOUT) :: CHANGE
Expand Down Expand Up @@ -2573,8 +2576,27 @@ SUBROUTINE DETERMINE_FIXED_ELEMENTS(T,CHANGE)

! If leakage and no vent area then set duct area to zero
IF (DU%LEAKAGE) THEN
IF ((NODE_AREA_EX(DU%NODE_INDEX(1)) < TWENTY_EPSILON_EB .AND. .NOT. DUCTNODE(DU%NODE_INDEX(1))%AMBIENT) &
.OR.(NODE_AREA_EX(DU%NODE_INDEX(2)) < TWENTY_EPSILON_EB .AND. .NOT. DUCTNODE(DU%NODE_INDEX(2))%AMBIENT)) DU%AREA=0._EB
ZERO_AREA = .FALSE.
IF (ALLOCATED(CONNECTED_ZONES)) THEN
IF (DUCTNODE(DU%NODE_INDEX(1))%ZONE_INDEX>0 .AND. DUCTNODE(DU%NODE_INDEX(2))%ZONE_INDEX>0) THEN
IF (CONNECTED_ZONES(0,DUCTNODE(DU%NODE_INDEX(1))%ZONE_INDEX)==1 .AND. &
CONNECTED_ZONES(0,DUCTNODE(DU%NODE_INDEX(2))%ZONE_INDEX)==1) ZERO_AREA = .TRUE.
ELSEIF (DUCTNODE(DU%NODE_INDEX(1))%ZONE_INDEX>0 .AND. DUCTNODE(DU%NODE_INDEX(2))%ZONE_INDEX<=0) THEN
IF (CONNECTED_ZONES(0,DUCTNODE(DU%NODE_INDEX(1))%ZONE_INDEX)==1) ZERO_AREA = .TRUE.
ELSEIF (DUCTNODE(DU%NODE_INDEX(1))%ZONE_INDEX<=0 .AND. DUCTNODE(DU%NODE_INDEX(2))%ZONE_INDEX>0) THEN
IF (CONNECTED_ZONES(0,DUCTNODE(DU%NODE_INDEX(2))%ZONE_INDEX)==1) ZERO_AREA = .TRUE.
ENDIF
ENDIF
IF (NODE_AREA_EX(DU%NODE_INDEX(1)) < TWENTY_EPSILON_EB .AND. .NOT. DUCTNODE(DU%NODE_INDEX(1))%AMBIENT) &
ZERO_AREA=.TRUE.
IF(NODE_AREA_EX(DU%NODE_INDEX(2)) < TWENTY_EPSILON_EB .AND. .NOT. DUCTNODE(DU%NODE_INDEX(2))%AMBIENT) &
ZERO_AREA=.TRUE.
IF (ZERO_AREA) THEN
DU%AREA = 0._EB
DU%VEL = 0._EB
ELSE
DU%AREA = DU%AREA_INITIAL
ENDIF
ENDIF

IF (DU%AREA_OLD <TWENTY_EPSILON_EB .NEQV. DU%AREA<TWENTY_EPSILON_EB) CHANGE = .TRUE.
Expand Down Expand Up @@ -2659,14 +2681,14 @@ END SUBROUTINE DETERMINE_FIXED_ELEMENTS
!> \param T Current time (s)

SUBROUTINE FIND_NETWORKS(CHANGE,T)

USE GLOBAL_CONSTANTS, ONLY: N_DUCTS,N_DUCTNODES,N_ZONE,CONNECTED_ZONES,N_NETWORKS,DUCT_NE,DUCTNODE_NE
INTEGER:: NZ,NN,ND,DUCT_COUNTER(N_DUCTS),NODE_COUNTER(N_DUCTNODES),COUNTER,COUNTER2,ZONE_COUNTER(N_ZONE),&
CONNECTED_ZONES2(0:N_ZONE,0:N_ZONE),NETWORK_DCOUNTER(N_DUCTS),NETWORK_NCOUNTER(N_DUCTS),RENUMBER(N_DUCTNODES+N_ZONE)
INTEGER:: NZ,NZ2,NN,ND,DUCT_COUNTER(N_DUCTS),NODE_COUNTER(N_DUCTNODES),COUNTER,ZONE_COUNTER(N_ZONE),&
CONNECTED_ZONES2(0:N_ZONE,0:N_ZONE),NETWORK_DCOUNTER(N_DUCTS),NETWORK_NCOUNTER(N_DUCTS),RENUMBER(0:N_DUCTNODES+N_ZONE)
LOGICAL, INTENT(OUT) :: CHANGE
LOGICAL :: CHANGE2
LOGICAL :: CHANGE2,ZERO_AREA
REAL(EB), INTENT(IN):: T
TYPE(DUCT_TYPE), POINTER :: DU
TYPE(DUCTNODE_TYPE), POINTER :: DN

CHANGE = .FALSE.
IF (.NOT. ALLOCATED(NETWORK)) CHANGE = .TRUE.
Expand All @@ -2678,21 +2700,38 @@ SUBROUTINE FIND_NETWORKS(CHANGE,T)
CALL DETERMINE_FIXED_ELEMENTS(T,CHANGE)

DO NN = 1, N_DUCTNODES
DN => DUCTNODE(NN)
ZERO_AREA = .TRUE.
DO ND = 1,DN%N_DUCTS
IF (DUCT(DN%DUCT_INDEX(ND))%AREA > TWENTY_EPSILON_EB) THEN
ZERO_AREA = .FALSE.
EXIT
ENDIF
ENDDO

IF (ZERO_AREA) CYCLE
NZ = DUCTNODE(NN)%ZONE_INDEX

IF (NZ>=1) THEN
NODE_COUNTER(NN) = FINDLOC(CONNECTED_ZONES2(NZ,:),1,1)-1
ZONE_COUNTER(NZ) = ZONE_COUNTER(NZ) + 1
IF (FINDLOC(CONNECTED_ZONES(NZ,:),1,1)-1>=1) THEN
NODE_COUNTER(NN) = FINDLOC(CONNECTED_ZONES(NZ,:),1,1)-1
ELSE
NODE_COUNTER(NN) = NN+N_ZONE
ENDIF
ELSE
NODE_COUNTER(NN) = NN+N_ZONE
ENDIF
ENDDO

IF (N_ZONE > 0) THEN
DO NZ = 1, N_ZONE
P_ZONE(NZ)%N_DUCTNODES = ZONE_COUNTER(NZ)
IF (.NOT. ALLOCATED(P_ZONE(NZ)%NODE_INDEX)) ALLOCATE(P_ZONE(NZ)%NODE_INDEX(N_DUCTNODES))
P_ZONE(NZ)%NODE_INDEX = 0
COUNTER = 1
DO NN = 1,N_DUCTNODES
IF (NODE_COUNTER(NN) == 0) CYCLE
IF (DUCTNODE(NN)%ZONE_INDEX == NZ) THEN
P_ZONE(NZ)%NODE_INDEX(COUNTER)=NN
COUNTER = COUNTER + 1
Expand All @@ -2702,12 +2741,13 @@ SUBROUTINE FIND_NETWORKS(CHANGE,T)
ENDIF

IF (ALLOCATED(NETWORK)) DEALLOCATE(NETWORK)

CHANGE2 = .TRUE.
DO WHILE (CHANGE2)
CHANGE2 = .FALSE.

DO ND = 1, N_DUCTS
DU => DUCT(ND)
IF (DU%AREA < TWENTY_EPSILON_EB) CYCLE
IF (NODE_COUNTER(DU%NODE_INDEX(1)) /= NODE_COUNTER(DU%NODE_INDEX(2))) THEN
CHANGE2 = .TRUE.
IF (NODE_COUNTER(DU%NODE_INDEX(1)) <=N_ZONE .AND. NODE_COUNTER(DU%NODE_INDEX(2)) <=N_ZONE) THEN
Expand All @@ -2724,67 +2764,80 @@ SUBROUTINE FIND_NETWORKS(CHANGE,T)
ENDDO

IF (N_ZONE > 0) THEN
DO NZ = 1, N_ZONE
COUNTER = FINDLOC(CONNECTED_ZONES2(NZ,:),1,1)-1
DO NN = 1, P_ZONE(NZ)%N_DUCTNODES
COUNTER = MIN(COUNTER,NODE_COUNTER(P_ZONE(NZ)%NODE_INDEX(NN)))
ENDDO
IF (COUNTER < NZ) THEN
COUNTER2 = FINDLOC(CONNECTED_ZONES2(COUNTER,:),1,1)-1
IF (COUNTER2 /= COUNTER) THEN
CHANGE2 = .TRUE.
COUNTER = COUNTER2
DO NZ=1,N_ZONE
IF (P_ZONE(NZ)%N_DUCTNODES==0) CYCLE
COUNTER = NZ
DO NZ2 = 1,NZ-1
IF (CONNECTED_ZONES2(NZ,NZ2)==1) THEN
IF (CONNECTED_ZONES(0,NZ2)==1) CYCLE
COUNTER = NZ2
EXIT
ENDIF
DO NN = 1, P_ZONE(NZ)%N_DUCTNODES
IF (NODE_COUNTER(P_ZONE(NZ)%NODE_INDEX(NN)) /= COUNTER) THEN
ENDDO
IF (COUNTER /=NZ) THEN
DO NN=1,N_DUCTNODES
IF (NODE_COUNTER(NN) == NZ) THEN
NODE_COUNTER(NN) = COUNTER
CHANGE2 = .TRUE.
NODE_COUNTER(P_ZONE(NZ)%NODE_INDEX(NN)) = COUNTER
ENDIF
ENDDO
ENDIF
ENDDO
ENDIF

END DO

RENUMBER = 0
N_NETWORKS = 0
DO NN = 1,N_DUCTNODES
IF (RENUMBER(NODE_COUNTER(NN)) == 0) THEN
N_NETWORKS = N_NETWORKS + 1
RENUMBER(NODE_COUNTER(NN)) = N_NETWORKS
IF (NODE_COUNTER(NN) > 0) THEN
IF (RENUMBER(NODE_COUNTER(NN)) == 0) THEN
N_NETWORKS = N_NETWORKS + 1
RENUMBER(NODE_COUNTER(NN)) = N_NETWORKS
ENDIF
ENDIF
ENDDO

DO NN = 1,N_DUCTNODES
NODE_COUNTER(NN) = RENUMBER(NODE_COUNTER(NN))
ENDDO
DO ND = 1,N_DUCTS
DUCT_COUNTER(ND) = RENUMBER(DUCT_COUNTER(ND))
ENDDO

IF (N_NETWORKS==0) RETURN

ALLOCATE(NETWORK(N_NETWORKS))
NETWORK%N_DUCTS=0
NETWORK%N_DUCTNODES=0
NETWORK_DCOUNTER=0
NETWORK_NCOUNTER=0
COUNTER = 0

DO ND = 1, N_DUCTS
NETWORK(DUCT_COUNTER(ND))%N_DUCTS = NETWORK(DUCT_COUNTER(ND))%N_DUCTS + 1
IF(DUCT_COUNTER(ND) > 0) NETWORK(DUCT_COUNTER(ND))%N_DUCTS = NETWORK(DUCT_COUNTER(ND))%N_DUCTS + 1
ENDDO
DO NN = 1, N_DUCTNODES
NETWORK(NODE_COUNTER(NN))%N_DUCTNODES = NETWORK(NODE_COUNTER(NN))%N_DUCTNODES + 1
IF(NODE_COUNTER(NN) > 0) NETWORK(NODE_COUNTER(NN))%N_DUCTNODES = NETWORK(NODE_COUNTER(NN))%N_DUCTNODES + 1
ENDDO

DO NN = 1, N_NETWORKS
ALLOCATE(NETWORK(NN)%DUCT_INDEX(NETWORK(NN)%N_DUCTS))
NETWORK(NN)%DUCT_INDEX = 0
ALLOCATE(NETWORK(NN)%NODE_INDEX(NETWORK(NN)%N_DUCTNODES))
NETWORK(NN)%NODE_INDEX = 0
ALLOCATE(NETWORK(NN)%MATRIX_INDEX(NETWORK(NN)%N_DUCTS+NETWORK(NN)%N_DUCTNODES))
NETWORK(NN)%MATRIX_INDEX = 0
ENDDO

DO ND = 1, N_DUCTS
IF (DUCT_COUNTER(ND)==0) CYCLE
NETWORK_DCOUNTER(DUCT_COUNTER(ND)) = NETWORK_DCOUNTER(DUCT_COUNTER(ND)) + 1
NETWORK(DUCT_COUNTER(ND))%DUCT_INDEX(NETWORK_DCOUNTER(DUCT_COUNTER(ND))) = ND
DUCT_NE(ND) = NETWORK_DCOUNTER(DUCT_COUNTER(ND))
ENDDO

DO NN = 1, N_DUCTNODES
IF (NODE_COUNTER(NN)==0) CYCLE
NETWORK_NCOUNTER(NODE_COUNTER(NN)) = NETWORK_NCOUNTER(NODE_COUNTER(NN)) + 1
NETWORK(NODE_COUNTER(NN))%NODE_INDEX(NETWORK_NCOUNTER(NODE_COUNTER(NN))) = NN
DUCTNODE_NE(NN) = NETWORK_NCOUNTER(NODE_COUNTER(NN))
Expand Down Expand Up @@ -3151,8 +3204,8 @@ END SUBROUTINE CONVERGENCE_CHECK
SUBROUTINE COLLAPSE_HVAC_BC(T)

! Takes the MPI gathered mesh array of HVAC boundary conditions and updates the DUCTNODE boundary condition values.
USE PHYSICAL_FUNCTIONS, ONLY : GET_SPECIFIC_GAS_CONSTANT,GET_ENTHALPY,GET_TEMPERATURE
USE MATH_FUNCTIONS, ONLY : EVALUATE_RAMP
USE PHYSICAL_FUNCTIONS, ONLY: GET_SPECIFIC_GAS_CONSTANT,GET_ENTHALPY,GET_TEMPERATURE
USE MATH_FUNCTIONS, ONLY: EVALUATE_RAMP
USE GLOBAL_CONSTANTS, ONLY: N_TRACKED_SPECIES,N_DUCTNODES,CONNECTED_ZONES,STRATIFICATION,TMPA,I_RAMP_TMP0_Z,I_RAMP_P0_Z,&
RSUM0,RHOA,P_INF,T_BEGIN,HVAC_PRES_RELAX,N_DUCTS,T_USED
REAL(EB), INTENT(IN) :: T
Expand Down Expand Up @@ -3487,7 +3540,7 @@ END SUBROUTINE LEAKAGE_HVAC

SUBROUTINE FILTER_UPDATE(DT,NODE_INDEX)

USE MATH_FUNCTIONS,ONLY:EVALUATE_RAMP
USE MATH_FUNCTIONS,ONLY: EVALUATE_RAMP
USE GLOBAL_CONSTANTS, ONLY: PREVIOUS
INTEGER,INTENT(IN)::NODE_INDEX
REAL(EB), INTENT(IN) :: DT
Expand Down Expand Up @@ -4439,7 +4492,7 @@ END SUBROUTINE FIND_DUCTRUNS

SUBROUTINE MATRIX_SOLVE_QFAN(DUCTRUN_INDEX,NF)

USE MATH_FUNCTIONS,ONLY : GAUSSJ
USE MATH_FUNCTIONS, ONLY: GAUSSJ
USE GLOBAL_CONSTANTS, ONLY: NEW
INTEGER,INTENT(IN) :: DUCTRUN_INDEX,NF
INTEGER :: IERR
Expand Down Expand Up @@ -4610,7 +4663,7 @@ END SUBROUTINE LHSDUCT_QFAN
!> \param TSI Length of time fan has been operating (s)

REAL(EB) FUNCTION FAN_MAX(FAN_INDEX,TSI)
USE MATH_FUNCTIONS, ONLY : EVALUATE_RAMP
USE MATH_FUNCTIONS, ONLY: EVALUATE_RAMP
INTEGER, INTENT(IN) :: FAN_INDEX
REAL(EB),INTENT(IN) :: TSI

Expand Down Expand Up @@ -4687,8 +4740,8 @@ END SUBROUTINE CONVERGENCE_CHECK_QFAN

SUBROUTINE COIL_UPDATE_QFAN(T,DUCTRUN_INDEX,NF)

USE MATH_FUNCTIONS, ONLY : EVALUATE_RAMP
USE PHYSICAL_FUNCTIONS, ONLY : GET_SPECIFIC_HEAT, GET_AVERAGE_SPECIFIC_HEAT, GET_ENTHALPY
USE MATH_FUNCTIONS, ONLY: EVALUATE_RAMP
USE PHYSICAL_FUNCTIONS, ONLY: GET_SPECIFIC_HEAT, GET_AVERAGE_SPECIFIC_HEAT, GET_ENTHALPY
USE GLOBAL_CONSTANTS, ONLY: N_TRACKED_SPECIES,OLD,NEW
INTEGER, INTENT(IN) :: DUCTRUN_INDEX,NF
REAL(EB), INTENT(IN) :: T
Expand Down
Loading