Characterization of liquefaction response in terms of triggering and manifestation at sites with soils containing significant amounts of fines is critically important. Saturated fines-dominated soils, that are known as liquefaction-susceptible materials, would represent a range of extreme responses to earthquake excitations from cyclic mobility to cyclic liquefaction causing excessive displacements. In a multi-layered soil profile, the interaction of adjacent layers in porewater pressure redistribution along the site profile plays an important role in liquefaction response. In this study, an advanced Nonlinear Dynamic Analysis (NDA) procedure in predicting the liquefaction mechanism and manifestation at a site with multiple fines-dominated soil layers is presented. The procedure specifically aims at simulating the interactions in a multi-layered system containing both sand-like and clay-like materials. The liquefaction data recorded by Wildlife Liquefaction Array (WLA) in 1987 Superstition Hills Earthquake (M = 6.6) was deployed to validate the NDA procedure. WLA site consists of sandy silt, silty sand, clayey silt, silty clay, and silt deposits in its upper 20 m and has provided invaluable records of porewater pressure and ground motions. Two critical-state-based constitutive models, PM4Silt and PM4Sand, in a coupled dynamic-fluid finite difference platform were employed. Calibration of the constitutive models was conducted using Cyclic Direct Simple Shear (DSS) test single-element modeling. The simulated response of the site from the proposed approach was compared with those of other NDA procedures conducted on this case history in terms of liquefaction parameters such as excess porewater pressure ratio and maximum shear strain. These comparisons revealed that incorporating multi-layer effects in an advanced NDA procedure can enhance the liquefaction behavior predictions at sites with fines-dominated materials.