Skip to content

Surf layer delamination (fall-off), based on density or temperature #14823

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
shostikk merged 4 commits into from
Jul 3, 2025
Merged

Surf layer delamination (fall-off), based on density or temperature #14823

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
4 commits
Select commit Hold shift + click to select a range
05dae41
FDS Source: wall solver layer delamination (falloff) feature
shostikk Jul 3, 2025
3a7fd1b
Merge remote-tracking branch 'firemodels/master'
shostikk Jul 3, 2025
c5199cc
FDS Source: shorten the comment
shostikk Jul 3, 2025
a435c9e
FDS Source: fix line length
shostikk Jul 3, 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
16 changes: 15 additions & 1 deletion Source/read.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -7809,7 +7809,8 @@ SUBROUTINE READ_SURF(QUICK_READ)
DRAG_COEFFICIENT,MINIMUM_BURNOUT_TIME,DELTA_TMP_MAX,BURN_DURATION,&
REFERENCE_HEAT_FLUX(MAX_QDOTPP_REF),REFERENCE_HEAT_FLUX_TIME_INTERVAL,MINIMUM_SCALING_HEAT_FLUX,&
MAXIMUM_SCALING_HEAT_FLUX,REFERENCE_THICKNESS(MAX_QDOTPP_REF),&
AREA_MULTIPLIER,Z_0,PARTICLE_EXTRACTION_VELOCITY,NEAR_WALL_EDDY_VISCOSITY
AREA_MULTIPLIER,Z_0,PARTICLE_EXTRACTION_VELOCITY,NEAR_WALL_EDDY_VISCOSITY,&
FALLOFF_TMP(MAX_LAYERS),FALLOFF_DENSITY(MAX_LAYERS)
INTEGER :: NPPC,N,IOS,NL,NN,NNN,NNNN,N_LIST,LEAK_PATH(2),DUCT_PATH(2),RGB(3),NR,IL
INTEGER :: N_LAYER_CELLS_MAX(MAX_LAYERS),VEG_LSET_FUEL_INDEX,INDEX_LIST(MAX_MATERIALS**2)
INTEGER :: CHILD_LAYER(MAX_LAYERS,N_MATL),CHILD_SURF(N_MATL)
Expand All @@ -7828,6 +7829,7 @@ SUBROUTINE READ_SURF(QUICK_READ)
DSC_CONVERSION_FACTOR,DT_INSERT,E_COEFFICIENT,&
EMBER_GENERATION_HEIGHT,EMBER_IGNITION_POWER_MEAN,EMBER_IGNITION_POWER_SIGMA,EMBER_TRACKING_RATIO,EMBER_YIELD,&
EMISSIVITY,EMISSIVITY_BACK,EXTERNAL_FLUX,EXTINCTION_TEMPERATURE,&
FALLOFF_DENSITY,FALLOFF_TMP,&
FREE_SLIP,INERT_Q_REF,FYI,GEOMETRY,HEAT_OF_VAPORIZATION,&
HEAT_TRANSFER_COEFFICIENT,HEAT_TRANSFER_COEFFICIENT_BACK,HEAT_TRANSFER_COEFFICIENT_SIGMA,&
HEAT_TRANSFER_MODEL,HORIZONTAL,HRRPUA,VARIABLE_THICKNESS,HT3D,ID,IGNITION_TEMPERATURE,&
Expand Down Expand Up @@ -8761,6 +8763,16 @@ SUBROUTINE READ_SURF(QUICK_READ)
ENDIF
ENDDO

! Set delamination temperature and/or density

SF%FALLOFF_TMP = HUGE_EB
SF%FALLOFF_DENSITY = -1._EB
DO NL = 1,MAX_LAYERS
IF (FALLOFF_TMP(NL) > 0._EB) SF%FALLOFF_TMP(NL) = FALLOFF_TMP(NL) + TMPM
IF (FALLOFF_DENSITY(NL) > 0._EB) SF%FALLOFF_DENSITY(NL) = FALLOFF_DENSITY(NL)
ENDDO


! Adjust min and max temperature and convert C to K

IF (TMP_FRONT >= -TMPM) TMPMIN = MIN(TMPMIN,TMP_FRONT+TMPM)
Expand Down Expand Up @@ -8961,6 +8973,8 @@ SUBROUTINE SET_SURF_DEFAULTS
EXTERNAL_FLUX = 0._EB
EXTINCTION_TEMPERATURE = -273._EB
INERT_Q_REF = .FALSE.
FALLOFF_TMP = -1._EB
FALLOFF_DENSITY = -1._EB
FREE_SLIP = .FALSE.
FYI = 'null'
GEOMETRY = 'CARTESIAN'
Expand Down
3 changes: 2 additions & 1 deletion Source/type.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -964,7 +964,8 @@ MODULE TYPES
REAL(EB), ALLOCATABLE, DIMENSION(:) :: CELL_SIZE !< Specified constant cell size (m)
REAL(EB), ALLOCATABLE, DIMENSION(:) :: STRETCH_FACTOR
REAL(EB), DIMENSION(MAX_LAYERS) :: LAYER_DENSITY,&
MOISTURE_CONTENT,SURFACE_VOLUME_RATIO,PACKING_RATIO,KAPPA_S=-1._EB,RENODE_DELTA_T
MOISTURE_CONTENT,SURFACE_VOLUME_RATIO,PACKING_RATIO,KAPPA_S=-1._EB,RENODE_DELTA_T
REAL(EB), DIMENSION(MAX_LAYERS) :: FALLOFF_TMP,FALLOFF_DENSITY !< Temp. and density criteria for layer fall-off (delam.)
REAL(EB), DIMENSION(MAX_LAYERS,MAX_MATERIALS) :: DENSITY_ADJUST_FACTOR=1._EB,RHO_S
CHARACTER(LABEL_LENGTH), ALLOCATABLE, DIMENSION(:) :: MATL_NAME
CHARACTER(LABEL_LENGTH), DIMENSION(MAX_LAYERS,MAX_MATERIALS) :: MATL_ID
Expand Down
37 changes: 31 additions & 6 deletions Source/wall.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2043,13 +2043,22 @@ SUBROUTINE SOLID_HEAT_TRANSFER(NM,T,DT_BC,PARTICLE_INDEX,WALL_INDEX,CFACE_INDEX,
B1%N_SUBSTEPS = 1
ONE_D%DELTA_TMP = 0._EB

! Compute number of cells and density vector
NWP = SUM(ONE_D%N_LAYER_CELLS(1:ONE_D%N_LAYERS))
RHO_S = 0._EB
DO I = 1,NWP
DO N=1,ONE_D%N_MATL
RHO_S(I)=RHO_S(I)+ONE_D%MATL_COMP(N)%RHO(I)
ENDDO
ENDDO

! Compute initial thermal properties for Fo number calculation
! CHECK_FO=F by default
! CHECK_FO=T enforces near explicit time step accuracy (Forward-Euler would require 1/2 DT_FO)


CHECK_FO_IF: IF (CHECK_FO) THEN
ONE_D%K_S = 0._EB
RHO_S = 0._EB
ONE_D%RHO_C_S = 0._EB
POINT_LOOP0: DO I=1,NWP
VOLSUM = 0._EB
Expand All @@ -2060,7 +2069,6 @@ SUBROUTINE SOLID_HEAT_TRANSFER(NM,T,DT_BC,PARTICLE_INDEX,WALL_INDEX,CFACE_INDEX,
VOLSUM = VOLSUM + ONE_D%MATL_COMP(N)%RHO(I)/RHO_ADJUSTED(LAYER_INDEX(I),N)
ONE_D%K_S(I) = ONE_D%K_S(I) + ONE_D%MATL_COMP(N)%RHO(I)*ML%K_S(ITMP)/RHO_ADJUSTED(LAYER_INDEX(I),N)
ONE_D%RHO_C_S(I) = ONE_D%RHO_C_S(I) + ONE_D%MATL_COMP(N)%RHO(I)*ML%C_S(ITMP)
RHO_S(I) = RHO_S(I) + ONE_D%MATL_COMP(N)%RHO(I)
ENDDO MATERIAL_LOOP0
IF (SF%PACKING_RATIO(LAYER_INDEX(I))>0._EB) ONE_D%K_S(I) = ONE_D%K_S(I)*SF%PACKING_RATIO(LAYER_INDEX(I))

Expand All @@ -2079,12 +2087,11 @@ SUBROUTINE SOLID_HEAT_TRANSFER(NM,T,DT_BC,PARTICLE_INDEX,WALL_INDEX,CFACE_INDEX,
LAYER_DIVIDE = SF%LAYER_DIVIDE

COMPUTE_GRID: IF (ONE_D%PYROLYSIS_MODEL==PYROLYSIS_PREDICTED .OR. SF%HT_DIM>1 .OR. SF%VARIABLE_THICKNESS) THEN
NWP = SUM(ONE_D%N_LAYER_CELLS(1:ONE_D%N_LAYERS))
CALL GET_WALL_NODE_WEIGHTS(NWP,ONE_D%N_LAYERS,ONE_D%N_LAYER_CELLS(1:ONE_D%N_LAYERS),ONE_D%LAYER_THICKNESS,SF%GEOMETRY, &
ONE_D%X(0:NWP),LAYER_DIVIDE,DX_S(1:NWP),RDX_S(0:NWP+1),RDXN_S(0:NWP),DX_WGT_S(0:NWP),DXF,DXB,&
LAYER_INDEX(0:NWP+1),MF_FRAC(1:NWP),SF%INNER_RADIUS,ONE_D%LAYER_DIVIDE_DEPTH)
ELSE COMPUTE_GRID
NWP = SF%N_CELLS_INI
! NWP = SF%N_CELLS_INI
DXF = SF%DXF
DXB = SF%DXB
DX_S(1:NWP) = SF%DX(1:NWP)
Expand Down Expand Up @@ -2413,6 +2420,25 @@ SUBROUTINE SOLID_HEAT_TRANSFER(NM,T,DT_BC,PARTICLE_INDEX,WALL_INDEX,CFACE_INDEX,

ENDDO POINT_LOOP2

! Delamination layer fall-off

I = 0
DO NL=1,ONE_D%N_LAYERS
IF (ONE_D%LAYER_THICKNESS(NL).GT.ONE_D%MIN_LAYER_THICKNESS(NL) .AND. &
(ONE_D%TMP(I+ONE_D%N_LAYER_CELLS(NL)) > SF%FALLOFF_TMP(NL)) .OR.&
( RHO_S(I+ONE_D%N_LAYER_CELLS(NL)) < SF%FALLOFF_DENSITY(NL))) THEN
REGRID_FACTOR(1:I+ONE_D%N_LAYER_CELLS(NL))=0._EB
ONE_D%RHO_C_S(1:I+ONE_D%N_LAYER_CELLS(NL))=0._EB
Q_S(1:I+ONE_D%N_LAYER_CELLS(NL))=0._EB
DO N=1,ONE_D%N_MATL
ONE_D%MATL_COMP(N)%RHO(1:I+ONE_D%N_LAYER_CELLS(NL))=0._EB
ENDDO
ONE_D%TMP(I+1:I+ONE_D%N_LAYER_CELLS(NL))=ONE_D%TMP(I+ONE_D%N_LAYER_CELLS(NL)+1)
B1%TMP_F = ONE_D%TMP(I+ONE_D%N_LAYER_CELLS(NL)+1)
ENDIF
I = I + ONE_D%N_LAYER_CELLS(NL)
ENDDO

! Compute new coordinates if the solid changes thickness. Save new coordinates in X_S_NEW.
! Remesh layer if any node goes to zero thickness

Expand Down Expand Up @@ -2465,7 +2491,7 @@ SUBROUTINE SOLID_HEAT_TRANSFER(NM,T,DT_BC,PARTICLE_INDEX,WALL_INDEX,CFACE_INDEX,
ENDIF
I = I + ONE_D%N_LAYER_CELLS(NL)
ENDDO

! If any nodes go to zero, apportion Q_S to surrounding nodes.

IF (ANY(CELL_ZERO(1:ONE_D%N_LAYERS)) .AND. NWP > 1) THEN
Expand Down Expand Up @@ -2661,7 +2687,6 @@ SUBROUTINE SOLID_HEAT_TRANSFER(NM,T,DT_BC,PARTICLE_INDEX,WALL_INDEX,CFACE_INDEX,
CALL INTERPOLATE_WALL_ARRAY(N_CELLS,NWP,NWP_NEW,INT_WGT,Q_S(1:N_CELLS))
CALL INTERPOLATE_WALL_ARRAY(N_CELLS,NWP,NWP_NEW,INT_WGT,RHO_H_S(1:N_CELLS))
CALL INTERPOLATE_WALL_ARRAY(N_CELLS,NWP,NWP_NEW,INT_WGT,TMP_S(1:N_CELLS))

DO I=1,NWP_NEW
VOL = (THICKNESS+SF%INNER_RADIUS-X_S_NEW(I-1))**I_GRAD-(THICKNESS+SF%INNER_RADIUS-X_S_NEW(I))**I_GRAD
TMP_S(I) = TMP_S(I) / VOL
Expand Down
Loading