## 1 Introduction  
The Athabasca Oil Sands Region in Alberta, Canada, is home to the third-largest oil reserves after Venezuela and Saudi Arabia. As of 2016, Alberta's oil sands proven reserves were 165.4 billion barrels. Surface mining recovers 20% of this total volume and covers an area of about 4800 km\({}^{2}\)[1]. Fluid fine tailings (FFT) are a by-product of surface bitumen extraction activities, which have accumulated a total volume of 1.3 billion m\({}^{3}\)[2]. The government of Alberta views oil sands operations as temporary land use and requires oil sands operators to return disturbed land to equivalent land capability [3]. Reclamation efforts for such large-scale disturbances require the reconstruction of new landscape-scale ecosystems [4, 5]. Revegetation is key to reclamation in a terrestrial landscape. Therefore, reclamation soil covers are designed to provide an adequate water supply for vegetation over dry summer periods [6, 7]. However, most reclamation soil cover research does not explicitly consider the 3D nature of differential settlement within the context of large-scale landscapes, such as oil sand tailings deposits. Understanding the dynamic nature in which reclamation takes place is fundamental in re-establishing ecosystem functionality. One of the constraints in landscape design is to limit the free water area to less than 10-20% of the reclaimed deposit [8, 9]. However, differential settlement may result in shallow ponded water and wetlands or may deepen or expand the existing wet landscape features. There is a wealth of hydrological studiesthat quantify downslope movement of water within reclaimed upland slopes and temporal variability of soil water distribution [10; 11; 12]. However, little effort seems to be directed at the spatial variability of underlying tailings with respect to geotechnical stability and settlement. The large strain consolidation theory [13] is commonly used to model the non-linear behavior of oil sand tailings consolidation. The spatial variability is of significance when reclamation covers overlay such highly compressible FFT. There is certainly a need to understand not only the settlement with time for these covers but also any potential 3D surface roughness or differential settlement, as it directly affects the contour of the final surface. In addition to reclamation and landscape designing, analyzing the consolidation of the treated tailings is also an essential input to predict the storage capacity of tailings storage facilities. In any tailings storage facility, the three main stages involved in settlements of fine-grained tailings material are flocculation, sedimentation and consolidation [14]. As the consolidation stage dominates the long-term settlement, the focus is typically on analyzing the consolidation stage [15]. A one-dimensional (1D) column consolidation is the basis for most of the current modeling practices and has been extrapolated to the whole deposits' consolidation behavior [16; 17; 18]. Another example of approximating the 3D consolidation based on a 1D large strain model is presented in [19], which studied the 3D analysis of a large strain consolidation model for tailings that are deposited in a mining pit. The goal was to estimate the filling time and tailings storage capacity. The authors performed Hydraulic Consolidation Test (HCT) to obtain the consolidation constitutive relationships, and generated approximate solutions of the 3D analysis. The proposed algorithm discretizes the pit volume in columns and layers to render upper and lower bound solutions based on a 1D large strain model. The quasi-3D models have been developed to estimate the consolidation of tailings deposits, most of them based on a series of 1D columns of concentric cylinders or adjacent columns [20; 21; 22]. There exists a critical assumption in quasi-3D consolidation models; at the same elevation, the sides of the tailings storage facility face the same deformation as the tailings material faces in the center of the deposit.