Opening-mode fractures in layered materials, such as sedimentary rocks, pavement, functionally graded composite materials or surface coating films, often are periodically distributed with spacings scaled to the thickness of the fractured layer. The current general explanation is that when the fracture spacing to layer thickness ratio changes from greater than to less than critical values the normal stress acting perpendicular to the fractures changes from tensile to compressive. This stress state transition is believed to preclude further infilling of fractures, and the critical fracture spacing to layer thickness ratio at this point defines a lower limit, called fracture saturation. To better understand the controls on fracture spacing, we have investigated the problem using a progressive fracture modeling approach that shares many of the natural kinematic features, such as fracture nucleation, fracture infilling and fracture termination. As observed experimentally, our numerical simulations demonstrate that fracture spacing initially decreases as extensional strain increases in the direction perpendicular to the fractures, and at a certain ratio of fracture spacing to layer thickness, no new fractures nucleate (saturated). Beyond this point, the additional strain is accommodated by further opening of existing fractures: the spacing then simply scales with layer thickness, creating fracture saturation. An important observation from our fracture modeling is that saturation may also effectively be achieved by the interface delamination and throughgoing fracturing, which inhibit additional layer-confined fracturing. We believe that these processes may serve as another mechanism to accommodate additional strain for a fracture saturated layer. Because interface debonding stops the transition of stress from the neighboring layers to the embedded central layer, which may preclude further infilling of new fractures, our fracture modeling approach predicts a larger critical length scale of fracture spacing than that predicted by a stress analysis approach based on stress transition theory. Numerical simulations also show that the critical value of the fracture spacing to layer thickness ratio is strongly dependent on the mechanical disorder m the fractured layer. The spacing to thickness ratio decreases with increasing heterogeneity of the mechanical properties.