Small-scale disturbances in the stratigraphy of the NEEM ice core: observations and numerical model simulations

Disturbances on the centimetre scale in the stratigraphy of the North Greenland Eemian Ice Drilling (NEEM) ice core (North Greenland) can be mapped by an optical line scanner as long as the ice has visual layering, such as, for example, cloudy bands. Different focal depths allow, to a certain extent, a three-dimensional view of the structures. In this study we present a detailed analysis of the visible folds, discuss their characteristics and frequency, and present examples of typical fold structures. We also analyse the structures with regard to the deformation boundary conditions under which they formed. The structures evolve from gentle waves at about 1500 m to overturned z folds with increasing depth. Occasionally, the folding causes significant thickening of layers. Their similar fold shape indicates that they are passive features and are probably not initiated by rheology differences between alternating layers. Layering is heavily disturbed and tracing of single layers is no longer possible below a depth of 2160 m. C axes orientation distributions for the corresponding core sections were analysed, where available, in addition to visual stratigraphy. The data show axialplane parallel strings of grains with c axis orientations that deviate from that of the matrix, which shows a single maximum fabric at the depth where the folding occurs. Numerical modelling of crystal viscoplastic deformation and dynamic recrystallisation was used to improve the understanding of the formation of the observed structures during deformation. The modelling reproduces the development of bands of grains with a tilted-lattice orientation relative to the single maximum fabric of the matrix, and also the associated local deformation. We conclude from these results that the observed folding can be explained by formation of these tilted-lattice bands.


Introduction
The NEEM (North Greenland Eemian Ice Drilling) ice core, located at 77 • 27 N 51 • 3.6 W in the northwest of Greenland, was drilled between June 2008 and July 2012.It is located on a topographic ridge, which dips towards the northwest; therefore the surface velocities on the ice divide have a non-negligible component of along-ridge flow of about 6 m a −1 (NEEM community members, 2013).In July 2010 the bedrock was reached at 2537.36 m depth.The site was chosen in order to recover an undisturbed Eemian warmperiod ice layer.However, it was found later that the ice below 2200 m was heavily disturbed and probably folded on a large scale (NEEM community members, 2013).
Visual stratigraphy of the NEEM ice core revealed folding on a small scale too, with fold amplitudes varying from less than 1 cm to a few decimetres (Samyn et al., 2011).These types of folds occur well above the large-scale disturbances reported by the NEEM community members (2013).Similar structures have been found in the lower parts of other deep ice cores (Alley et al., 1997;Thorsteinsson, 1996;Svensson et al., 2005;Faria et al., 2010;Fitzpatrick et al., 2014).Stratigraphy bands are visualised by an indirect light source scattering on surfaces inside the ice, mainly particles and air Published by Copernicus Publications on behalf of the European Geosciences Union.D. Jansen et al.: Small-scale disturbances in the stratigraphy of the NEEM ice core bubbles/hydrates (Svensson et al., 2005).High impurity content is found in ice that originates from snow accumulated during glacial periods.Changing impurity contents between ice from glacial and interglacial periods have been linked to rheological differences (e.g.Paterson, 1991) and may lead to shear localisation in distinct horizontal layers.
Due to their potential influence on the integrity of the climatic record, folds have been subject to modelling studies (e.g.Waddington et al., 2001).Thorsteinsson and Waddington (2002) explored the amplification of small disturbances in the layering of ice cores for isotropic and anisotropic conditions, investigating the potential for the existence of overturned folds near ice sheet centres.Azuma and Goto-Azuma (1996) concluded from their proposed anisotropic flow law formulation that an inclined single maximum fabric could lead to vertical strain even in simple shear and thus influence the stratigraphy.They also suggested that horizontal variations in the inclinations could then cause alternating thickening and thinning of layers, leading to folding or boudinage in the stratigraphy.However, the initial formation of the disturbances is not fully understood.
Here we present a characterisation of the small-scale folding observed in the NEEM ice core.Another feature occasionally observed along with folding in deep ice cores are "fabric stripes" (Alley et al., 1997).They describe bands of deviating grain orientations with respect to the surrounding matrix, which is essentially a single maximum fabric in regions where folding occurs.We discuss possible folding mechanisms and the link to the so-called "Alley stripes" in the crystal fabric of grains.Microstructural modelling with ELLE reproduces similar fabrics and fold structures to the ones we observe in the NEEM ice core.

Methods
The data used in this study were obtained by different observational methods, which will be introduced only briefly in the following section.For technical details we recommend checking the original literature cited in the subsections.

Line scan visual stratigraphy
The visual stratigraphy of the NEEM ice core was recorded by means of an automated line scan instrument (see Svensson et al. (2005) for a detailed description of the instrument and data from the North GRIP ice core).Clear ice appears dark when illuminated by an indirect light source.Dust particles or bubbles cause scatter of light and make the ice appear bright in the line scan image.A clear correlation between backscatter and dust content has been found in the North GRIP ice core (Svensson et al., 2005).The method can be applied directly in the field and in the case of the NEEM ice core was applied continuously for the entire core, with a gap between 860 and 1150 m, which corresponds to the brittle zone where the core quality did not allow preparation for the line scanner.For the NEEM ice core the line scan images were recorded with a standardised exposure time and three focal planes within the ice core section with a vertical distance of 1 cm (Kipfstuhl, 2010).This allows, to a certain degree, a three-dimensional mapping of the visible layering in the ice core.The data are stored in high-resolution (118 pixel per centimetre) bitmap images.One drawback of this method is, of course, that it only shows disturbances in the ice if scattering surfaces are present.However, it is possible to even reveal structures at low dust content by means of image processing and filtering.

Automated fabric analyser
The crystal fabric orientation of discrete samples was measured using a G50 Automatic Fabric Analyser (Australian Russell-Head type; see e.g.Peternell et al., 2010, data set: Weikusat andKipfstuhl, 2010).Samples cut from the physical properties part of the NEEM core were cut to 250 µm thin sections to measure c axis crystal lattice orientations by polarised light microscopy, where the thin section is placed between systematically varying crossed polarisers (e.g.Wilson and Russell-Head, 2003).The data coverage is much better than in previous ice cores with continuous sampling of selected core sections (bags) to investigate metre-scale variations in fabric throughout the core.However, due to the timeconsuming preparation of the samples it was not possible to produce a continuous record.

Microstructural modelling with ELLE and full field crystal plasticity
We use two-dimensional numerical modelling to investigate the development of strain localisation in a polycrystalline aggregate.At the moment it is not possible to combine simple shear and pure shear boundary conditions.We chose simple shear boundary conditions as an approximation for the in situ conditions in the lower part of the ice sheet where the bands of deviating grain orientation are observed, and where horizontal shear is dominant.The simulation approach couples a full field method based on the fast Fourier transform (FFT) that simulates viscoplastic deformation, with a front tracking code that simulates dynamic recrystallisation processes (DRX), included within the open-source numerical modelling platform ELLE (http://www.elle.ws;Bons et al., 2008).ELLE has been successfully used to simulate evolution of microstructures during deformation, such as recrystallisation (Piazolo et al., 2010;Roessiger et al., 2011Roessiger et al., , 2014) ) or strain localisation (Jessell et al., 2005;Griera et al., 2011Griera et al., , 2013)).The full-field crystal plasticity (FFT) code (Lebensohn, 2001;Lebensohn et al., 2008;Montagnat et al., 2014a, Llorens et al., 2016) simulates deformation by pure viscoplastic dislocation glide.An experimental run consists of iterative applications of small increments ( γ = 0.04) of simple shear deformation, each followed by a sub-loop of processes simulating dynamic recrystallisation (grain boundary migration and recovery).While grain boundary migration covers the motion of high-angle grain boundaries, recovery achieves a decrease in intra-crystalline heterogeneities by means of local rotation without motion of high-angle boundaries.The recrystallisation sub-loop may be called more than once to simulate the different balance between deformation and recrystallisation as a function of strain rate, since all simulations are performed with the same intrinsic mobility value (M 0 ) and boundary diffusion activation energy (Q) (Roessiger et al., 2014;Llorens et al., 2016).Exchange of data between ELLE and FFT is possible, as both use periodic boundary conditions and the physical space is discretised into a shared regular node mesh.
The ELLE data structure consists of three layers: (1) a network of nodes (boundary nodes or bnodes) that are connected by straight boundary segments that define the highangle grain boundaries that enclose individual ice grains, (2) a set of unconnected nodes (unodes) to map lattice orientations and dislocation densities, used for the FFT calculation, and (3) a passive marker grid utilised to track finite strain.Distances between nodes are kept between 5.5 × 10 −3 and 2.5 × 10 −3 times the unit distance (in a 1x1 bounding box), by removing bnodes when their neighbours are too close or adding bnodes when two nodes are too far apart.The space is discretised in a mesh of 256 × 256 Fourier points, resulting in a unit cell defined by 65 536 discrete nodes.Each unode represents a small area or crystallite with a certain lattice orientation, defined by Euler angles, and a dislocation density value.The ELLE data structure has fully wrapping boundaries.The 10 × 10 cm 2 initial microstructure has 1632 grains, each with a homogeneous lattice orientation, showing a c axis preferred orientation almost perpendicular to the shear plane, in order to simulate an intrinsic anisotropic material.The misorientation between grains was set at < 5 • (i.e. initial noise).Dislocation glide of icesingle crystal was defined by slip on the basal {0001} {11-20}, prismatic {1-100} {11-20} and pyramidal systems {11-22} {11-23}.In these simulations, the ratio A of critical resolved shear stress for non-basal versus basal slip systems was set to A = 20.The same stress exponent (n = 3) is set for all slip systems.The physical values used for recrystallisation are intrinsic mobility M 0 (0.023 m 2 kg −1 s; Nasello et al., 2005), boundary diffusion activation energy Q (40 KJ mol −1 ; Thorsteinsson, 2002), isotropic surface energy γ e (0.065 Jm −2 ; Ketcham and Hobbs, 1969), and temperature, which was set to T = −30 • C. To simulate recovery numerically, a modification of the approach proposed by Borthwick et al. (2013) was used.As the time step for the simulation of recrystallisation is smaller than necessary for the computationally expensive FFT calculation, we modelled 10 DRX steps of t = 6.3 × 10 8 s for each shear strain increment of γ = 0.04, giving a shear strain rate of 6.35 × 10 −12 s −1 .Simulations with other shear strain rates were performed for comparison, but are not presented here.See Llorens et al. (2016) for a complete description of the methods.

Stratigraphy and fold classification
The stratigraphic data were visually inspected for all parts of the ice core containing cloudy bands, in order to categorise disturbances of the visible layers.It has to be noted that this method is only appropriate where sufficient layers are visible, but clear ice may have been deformed as well.Figure 1 shows an overview of the layering structures we find in the NEEM ice core.All images shown are confined to the core sections above the major disturbances in the Eemian ice beginning at a depth of approximately 2200 m.Around this depth the ice is heavily sheared and the layering becomes more and more diffuse.Below that it is no longer possible to see fold structures in the visual stratigraphy data as the Eemian ice is mostly clear.The panels display the scans of entire core sections of about 1.10 m, which were cut into segments of 0.55 m after scanning.The top always represents the upper part of the core segment.Some segments differ in length, as the recovered core pieces are not always exactly 1.10 m long.
www.the-cryosphere.net/10/359/2016/The Cryosphere, 10, 359-370, 2016 Some of the pieces also fractured during the recovery process or during preparation, which is highlighted with red lines in Fig. 1.The images have been partly processed by applying a Gaussian filter to enhance the visibility of the layering, and therefore the grey values are no absolute measure for impurity content or of other parameters that could influence the backscatter within the ice.
The upper part of the NEEM ice core shows little or no disturbances.Figure 1a shows an example from 1430 m depth with perfectly horizontal layers.The layer thickness and opacity vary in the core segment, and single layers have a constant thickness throughout the 10 cm wide core section.A close-up of one of the layers shows no particular structure within it (Fig. 2a).According to the NEEM chronology published in Rasmussen et al. (2013), the annual layer thickness in this depth, which is the upper boundary of the glacial ice, is about 2 cm.Svensson et al. (2005) describe the alternating cloudy and transparent layers of ice in the NorthGRIP ice core as a result of depositional events, which do not necessarily reflect annual cycles.Below a depth of about 1700 m the annual layer thickness decreased to about 1 cm (Rasmussen et al., 2013), which is also reflected in thinner bands visible in the line scan image.In this depth the structure of the layering begins to change, as examples from depths of approximately 1760 and 1867 m show (Fig. 1b and c).Wave-like features with centimetre-scale amplitudes and wavelengths in the order of the core diameter can be observed.In some parts of the core segments these disturbances can be clearly followed through several layers.Figure 2b shows an enlarged section of Fig. 1b, showing a well-developed asymmetric z fold.Its shape indicates sinistral shear and the fold is beginning to overturn.The fold hinge is a sharp feature, which can be followed over several layers.The enlargements in Fig. 2b also show that the cloudy layers themselves appear to be laminated.
For the core sections shown in Fig. 1b and c the layers vary in thickness within the core, as can be clearly seen in Fig. 2b, where the central greyish layer (indicated with green dashed lines) nearly doubles its thickness in the centre of the image due to the folding.This shape is typical for so-called similar folds in geology (Ramsey and Huber, 1987).Figure 1d and e show examples from 1977 and 2098 m depth, where the layering is significantly more disturbed.The vertical scale of the disturbances has risen to a scale of 10 cm (Fig. 2d and e).In between the larger scale folds the layering appears to be more regular again, however the limited width of the core sections limits our interpretation here, as the layers could be overturned folds of which the limbs had become near-horizontal due to a combination of the ongoing shear and vertical thinning.Figure 2c shows a stack of flattened folds, where the doubling of layers may not be immediately obvious to the observer when focusing on the left part of the image.Tracing one boundary between a clear and a cloudy band highlights that the thin layers are probably limbs of an overturned fold (indicated with a blue dashed line).There are also new generations of folds standing out through their well-defined and steeper axial planes and which are not yet overturned (Fig. 2c, on the left, fold axis is indicated with a green dashed line).
At even greater depth the layering becomes less distinct (Fig. 1f-h).In some parts of these sections the layers appear to be undisturbed but inclined, which may indicate that they are part of a larger deformation structure.The now very thin layers still show new generations of folds.

Crystal fabric orientation anomalies connected to folds
In comparison to previous deep ice cores, the amount of data gathered to analyse ice fabric is relatively high.To investigate small-scale variations, entire bags of 55 cm from certain depths were processed.The general evolution of ice fabric with depth in the NEEM ice core was described in Montagnat et al. (2014b).The c axis orientation distribution develops more or less linearly from an isotropic fabric to a single maximum at a depth of about 1400 m, which represents the transition from the Holocene to the last glacial period (Rasmussen et al., 2013).Within the well-developed single maximum fabric we found inclined bands of grains with a deviating c axis orientation.We assume that the bands are planar features, but as the thin sections are vertical cuts through the cylindrical core section, the inclination of the bands is not necessarily equal to the inclination of the planes.Similar bands were described in the GRIP ice core (Thorsteinsson, 1996) and the GISP2 ice core (Alley et al., 1997).In the case of the NEEM ice core, however, significantly more fabric data are available, which enables us to follow these structures through entire core sections.One of the first examples of such a band, shown in Fig. 3a, appears at a depth of 1800 m.The c axis orientation of grains within the bands is tilted anti-clockwise relative to the single maximum, which is indicated by the blue-greenish colours in the colour wheel used to illustrate c axis orientation (inset in Fig. 3a).The grain size does not differ from the average grain size of the sample.The subgrain boundary density in these grains does not differ significantly from the surrounding ones, which indicates that they are most likely not newly nucleated grains (Fig. 4).However, while the subgrain boundaries in grains with vertical c axes are mostly parallel to the basal planes, they are mainly perpendicular to the basal plane in grains within the band, indicating the onset of rotation recrystallisation (Weikusat et al., 2009).
A direct comparison of the fabric data with the line scan images (Fig. 3) reveals that these bands are connected with  disturbances in the layering.The inclination is in agreement with the sense of shear that is derived from the asymmetry of the folded layers.However, layer disturbances are not always visible where fabric anomalies are found.
In Fig. 5a-c the three line scan images available from the different focal depths are plotted next to each other to illustrate the three-dimensional nature of the observed folding in the layering.The shape change of the highlighted layer indicates that the fold axis shifts to the left towards the centre of the core (Fig. 5e).The thin sections for the fabric analysis are prepared from the physical properties sample in the upper part.The line scan measurement is performed on the remaining part of the core with 1 cm between the different focal planes (Fig. 5d).
Figure 6 shows an example from ≈ 1978 m depth where we see finely laminated layers and asymmetric folds that indicate dextral shear.are not very distinct, which is probably due to light diffusion caused by the distorted fine layers.The two distinct bands in the right half of Fig. 6a are relatively steep and exhibit a small tilt in the c axes, while the feature indicated by the central arrow is more flattened and also shows a higher tilt in the c axes.Another example can be found in the Supplement (Fig. S1).Where several bands occur in one core section, their inclination and orientation appears to be consistent throughout (Fig. 7).

Model results
To understand the development of the observed fabric anomalies and the related disturbances in the layering, we simulated the fabric evolution under simple shear with an initially well-developed single maximum orientation distribution.A random noise of < 5 • was added to grain orientations.The setup of the simulation does not fully represent the probable kinematic boundary conditions in the region of the ice core where we observe the structures, which would be a combination of vertical compression and simple shear (Montagnat et al., 2014b).We, however, model these structures in simple shear for simplicity.This approach is reasonable, since there is a significant flow along the ridge of about 6 m a −1 (NEEM Community Members, 2013), which gives the core location the character of a flank with relatively high shear strain rates, rather than a divide.Moreover, the choice of simple shear boundary conditions is also justified by the fact that the bands start to appear in the lower third of the ice sheet, where shear stress becomes the dominant driver for deformation (Montagnat et al., 2014b).However, there might be aspects of the evolution of the folds that cannot be reproduced due to the simple shear approximation.
In the model simulation vertical bands similar to the ones observed in the ice core begin to stand out after a shear strain of γ = 0.6, but start to appear already after small strains (Fig. S3).In the initial state the small deviations of the c axes from the vertical are randomly distributed and do not show alignment (Fig. S3a).During the first deformation steps narrow vertical bands develop with deviations from the single maximum towards the shear direction as well as broader bands with an opposite rotation of the c axes.The rotation of c axes in the narrow bands intensifies during the next steps and the bands begin to tilt due to the continuing shear deformation (Fig. S2).The rotation of the c axes is twice the inclination of the band, which is typical for flexural-slip kink bands (Fig. S4; Dewey, 1965;Tanner, 1989).Figure 8a-c shows the c axis orientations for the sample after shear strains of γ = 1, γ = 2, and γ = 3.The bands seem to develop in different generations, which can be distinguished by their inclination as the new bands are steeper.There are areas between the bands where orientations of c axes rotate anti-clockwise (magenta coloured), but on a larger scale and with less welldefined boundaries (Fig. S3).In later stages of the simulation the oldest bands begin to disintegrate with the grains recrystallising back to a vertical c axis fabric.
Figure 8d-f show the development of a passive marker grid during the simulations.The blue lines were perfectly horizontal at the beginning of the simulation and can be regarded as an analogue to the stratigraphic layering observed in the ice core.It is apparent that the bands with abnormal grain orientation are connected with folding in the layering.At first these disturbances appear as small steps, but they develop into overturned folds with a short and steep limb with progressive deformation.They correspond to the welldeveloped bands in the fabric, and to a long, less inclined limb, representing the area in between the bands.The disturbances in the layering are permanent, and therefore the bands are visible in the passive grid even when they no longer exist in the orientation plot.The evolution of these bands and the corresponding folds is reminiscent to the formation of kink bands or chevron folds (Fleuty, 1964;Dewey, 1965).As these terms may have genetic connotations, referring to either single crystal processes (e.g.Wilson et al., 1986) or macroscopically strong foliated rocks under compression parallel to foliation (e.g.Cobbold et al., 1971), we use the term "tiltedlattice bands" here to describe the bands of grains with lattice orientations that deviate from the dominant lattice preferred orientation (LPO).The rotation rate of the tilted-lattice bands in the model run is controlled by the applied overall perfect simple shear deformation.If there is an additional component of vertical shortening, as in a natural ice sheet, the rotation rate is expected to increase.Figure 8g-i show the equivalent von Mises strain rate for the deformation steps γ = 1, γ = 2, and γ = 3.The strain rate appears to be localised around the margins of the bands where bending strain is the highest, which is most apparent for newer bands with steep inclinations.Fig. 8d-f also highlight strain localisation in horizontal zones, which is not visible in passively deformed horizontal lines.

General discussion of folds
The shape of the observed folds in the NEEM ice core is typical for similar folds, as the layers are thickened in the hinge region and thinned in the fold limbs.Similar folds are passive features, where all layers of the package are deformed in a similar way (Fig. 9a).They form by passive shearing of the layering and can evolve to become overturned z folds or even sheath folds (Quinquis et al., 1978;Bons and Urai, 1996;Alsop and Carreras, 2007).Competence or viscosity contrast between the different layers plays no or only a minor role.In contrast, buckle folds (Fig. 9b) develop when layers have different viscosities.A competence contrast with a ratio of at least about 25 between strong and weak layers is required to develop distinct folds (Llorens et al., 2013a, b).When a stack of strong and weak layers is shortened, the strong layers form folds by bending, which suppresses thickening or thinning of these layers.The weak layers accommodate this, bending by ductile flow into the hinge regions, a process known as flexural flow (Donath and Parker, 1964).Strong and weak layers are thus different in shape (Fig. 9b).The fold shapes observed in the NEEM core, however, appear to be consistent across a stack of several layers (Figs. 1 and 2), which indicates that viscosity contrasts are very low and the folds are formed by passive shearing, although there may be some differences in the flow strength of the ice between the layers due to different impurity content (Paterson, 1991).
Figure 10 gives an overview of the onset of folding (black line) and the evolution of an anisotropic fabric (red line) for several ice cores.Comparison with data from EDML (Faria et al., 2010) and WAIS (Fitzpatrick et al., 2014) in Antarctica and with GRIP (Thorsteinsson, 1996), GISP2 (Alley et al., 1997;Gow et al., 1997) and North GRIP (Svensson et al., 2005;Wang et al., 2006) from Greenland reveals that the onset of visible folding is dependent on the relation between vertical strain rates (shortening) and shear strain rates (Fig. 10).Due to the high vertical strain rates, fold structures are flattened out before they overturn, and are thus no longer visible.This has been theoretically described by Waddington et al. (2001).Thus, the dynamical setting of the borehole location is, in addition to the anisotropy, an essential parameter for the onset of visible folding.An ice core at flanks or on divides with non-negligible flow along the ridge samples ice which experiences more shear strain than an ice core at dome positions.While the GISP2 and North GRIP ice cores are very similar in ice thickness and accumulation rate to GRIP, the onset of folding for the latter is 300 m deeper, which may be due to its dome position and the lower surface velocity (for a comparison of ice core parameters see Faria et al., 2014).In the region of the NEEM ice core there is an even higher along-ridge flow.A comparison of shear strain rates profiles with depth at the NEEM, GRIP, and North GRIP locations can be found in Montagnat et al. (2014b).The depths of the crossing points of the curves for shear and vertical strain rate that are displayed in Montagnat et al. (2014b) correspond approximately to the onset of folding in the cores.
The EDML ice core stands out in the comparison shown in Fig. 10, as the folding begins significantly higher than the establishment of a single maximum fabric.However, Faria et al. (2010) report that a strong girdle fabric has formed in the region of the onset of folding, thus the fabric does show some anisotropy there as well.
The scale of the disturbances found in the layering of the NEEM ice core is very similar to the ones observed at EDML (Faria et al., 2010) and North GRIP (Svensson et al., 2005), for both of which a line scan data set of comparable quality as for the NEEM ice core is available.
The very deep onset of folding in the WAIS ice core might be due to strong basal melting at this site (Fitzpatrick et al., 2014;WAIS Divide Project Members, 2013).

Tilted-lattice bands initiate folding
Ice is a mechanically highly anisotropic mineral, as is polycrystalline ice with a strong single maximum c axes distribution.Our model results and the observations from the NEEM ice core show that this anisotropy and small perturbations thereof are an essential precondition for the development of folds on the centimetre scale.A similar process is common in well-foliated rocks, which therefore exhibit a strong mechanical anisotropy (Hudleston and Treagus, 2010).In that context, there are basically two mechanisms for the formation of kink or chevron folds (Dewey, 1965), which may act in concert.The first is shear localisation, which occurs approximately parallel to the planes that experience the highest net shear stress.The result is a conjugate set of kink bands, originally at a high angle to each other.The second mechanism is a combination of localised bending and flexural slip.There is no thickening or thinning perpendicular to the foliation if the material can only deform by slip parallel to that foliation.A geometric necessity of this type of folding is that the axial plane must be the bisector of the interlimb angle (Fig. 9c) (Frank and Stroh, 1952;Dewey, 1965;Cobbold et al., 1971).With the interlimb angle at the beginning of folding being 180 • , the bands with a slightly tilted-lattice orientation initiate at 90 • to the foliation.We assume that this second mechanism can be applied to ice with strong LPO and explains the observed folding process.
Simple shear parallel to a foliation is a special case where the orientation of tilted-lattice bands formed by both mechanisms is identical: 45 • to the maximum compression, which coincides with parallel and perpendicular to the foliation.Bands at a normal inclination to the foliation have indeed been observed in geological simple shear experiments (Misra and Burgh, 2012;Williams and Price, 1990), and developed in our numerical model.Bands parallel to the foliation are difficult to observe in natural samples, as these would not fold the foliation.However, the numerical model shows shear plane parallel localisation of strain rate as well (Fig. 8g-i).
The tilted-lattice bands rotate passively if there is a layerparallel shear component, which does not necessarily have to be the dominant deformation component.As the band rotates by an angle α to the long limb, the short limb has to rotate by 2α, as is observed in the NEEM core and numerical simulations (Fig. 8d-f).The relation between the tilt of the band and the corresponding tilt of the c axes from the model results is shown in Fig. S4.Alley et al. (1997) suggested that the rotation of the c axes of grains in the stripes lag behind the rotation of the bands themselves, while our model results indicate c axes rotating twice as fast as the axial plane or stripe.Alley et al. (1997) explained the growth in length by the sense of shear at the perimeter of the stripe, which would cause a "spinning" to adjacent grains, causing their c axes to rotate towards the inclination of the band and thus elongate it.In our numerical model we observe that tilted-lattice bands are seeded by individual grains and develop by linking these (Fig. S3).However, their intensification under rotation appears to propagate along their length (Fig. S2).
With progressive rotation of the tilted-lattice bands they become more distinct.Rotation of the short limb occurs by sliding parallel to the basal plane with a sense opposite to the overall shearing direction (Fig. 9c).In the case of chevron folds the features finally "lock up" when the interlimb angle reduces to about 90 • , i.e. when the bands are approximately 45 • to the layering (Dewey, 1965).In the numerical simulations we see that tilted-lattice bands begin to disintegrate at this stage, with recrystallisation and recovery consuming the grains with deviating orientations, and the flow homogenises again (see Fig. 8).However, marker lines, such as the layers in the NEEM core, will still record the tilted-lattice bands, which now continue shearing and develop into passive folds.
In summary, the model results indicate that the evolution of the observed tilted-lattice bands is a consequence of a fabric with a strong anisotropy with superimposed small random disturbances.In this way grains orientated unfavourably for basal glide are rotated by rigid body rotation, as well as shear along the basal plane with a shear sense opposite to the bulk shear strain.Thus, tilting bands appear to be an essential process in ice deformation under shear.Azuma and Goto-Azuma (1996) suggested that horizontal variation in the single maximum direction could explain heterogeneous layer thinning or thickening of initially horizontal layers, eventually leading to folding.The development of tilted-lattice bands is a process providing such variations in the fabric.
A difficulty in comparing the results of the simulation with the observational data is that with fabric measurements we can only capture a two-dimensional section of a threedimensional feature.Assuming that the tilted-lattice bands are planar features, the angle at which the cylinder of the ice core is cut relative to the inclination of the plane determines its appearance in the two-dimensional section.Thus, the inclination of the observed bands in the plane is not sufficient to describe the full orientation of the feature, but instead gives a minimum inclination of the plane.This also has to be taken into account when interpreting the fold structures on the line scan images.
Within one 55 cm section of the ice core (bag) the cutting plane through the core is consistent and so are the samples used to prepare the thin sections for the fabric measurements.Figure 7 shows that within one bag the inclinations of the tilted-lattice bands are consistent as well, strengthening the assumption that they are connected to the local stress environment and to the sense of shear, projected onto the plane of the thin section or line scan image.A consistency between the direction of shear and the shape of z folds in the stratigraphy has also been reported for the GRIP and GISP2 ice cores (Alley et al., 1995).In both examples displayed in Fig. 7, different generations of tilted-lattice bands can be detected, differing in inclination of the bands as the older bands have been subjected to more shear strain since their formation, and in the corresponding shift in c axes orientation, as it is seen in the model results as well.
The connection between bands with an anomalous LPO and stratigraphic disturbances has already been discussed before (Thorsteinsson, 1996;Alley et al., 1997;Samyn et al., 2011).Together with the microstructural model results, the observations can be interpreted with an improved understanding of the underlying process.The model results clearly show that tilted-lattice bands can form in simple shear conditions and that a well-developed anisotropic fabric with small perturbations is required.At the moment it is not clear why the bands sometimes appear dark in the line scan images, but from deeper parts of the core where the crystals are larger in size the line scan images give indication that the backscattering can be affected by the crystal orientation.

Summary and conclusions
The onset of small-scale folding can be observed at the start of the lower third of the NEEM ice core, which is similar to the fold evolution observed in EDML.Below a depth of about 2160 m it is no longer possible to track stratigraphic layers.The shape of the observed structures indicates that they are not buckle folds, which means that they are not originated by a contrast in competence between alternating layers.The amounts of folding as well as the state of disturbance increase with depth.
Folding causes thickening of cloudy bands and can potentially influence the resolution of climate data extracted from the NEEM ice core.Folding causing doubling of layers was observed below a depth of about 2100 m.In some core sections the layering appears to be intact in between larger folds in the line scan data.However, within a core section, only folds and doubling of layers up to the scale of 10 cm can be delineated with certainty.It is therefore difficult to ascertain that the climate signal is not disturbed in regions with parallel layering, as these could potentially be part of larger-scale folds.
Microstructural numerical modelling results indicate that the observed folding is initiated by the formation of bands with a tilted-lattice orientation relative to the bulk LPO.The formation mechanism requires a highly anisotropic material and thus a well-developed single maximum crystal orientation.Local deviations from the single maximum in the direction of shear provide the seeds for tilted-lattice bands and thus folding.Here we have shown that this process is active on the microstructural scale.The possible link between tilted-lattice bands and larger scale folds still has to be investigated, but could be in line with suggestions by Azuma and Goto-Azuma (1996).Grains with inclined basal plane orientations within the tilted-lattice bands are eventually eroded through recrystallisation and recovery.However, the tiltedlattice bands formed folds in material planes that further evolve by passive folding, which is visible in the layering, but not in the c axes patterns.
The Supplement related to this article is available online at doi:10.5194/tc-10-359-2016-supplement.

Figure 1 .
Figure 1.Visual stratigraphy overview.Line scan images of different depths.A Gaussian filter was applied to images shown in (a), (d), (e), (f), (g), and (h) to enhance the visibility of the layers.Red lines indicate fractures.Blue squares and associated figure codes indicate location of enlargements shown in Fig. 2. The 1.10 m line on the right indicates the scaling of the images and is also the typical length of a recovered core section.

Figure 2 .
Figure 2. Close-ups from the overview Fig. 1.(a) Undisturbed layering.(b) Angular z fold consistent throughout layering.The green dashed line indicates a layer discussed in Sect.3.1.(c) Different generation of folds.The dashed lines indicate features discussed in Sect.3.1.(d) Strongly disturbed layer significantly thickened.(e) Strongly disturbed layering with different generations of folds.

Figure 3 .
Figure 3.Comparison of fabric data and visual stratigraphy in detail, bag 3276, approximate depth 1803 m.(a) Fabric data in a vertical section, (b) line scan image in a vertical section, (c) stereo plot of c axes orientations (horizontal plane).

Figure 4 .
Figure 4. (a) Close-up of tilted-lattice band grains at approximately 1803 m depth (bag 3276).Inset shows the colour code for c axes orientation, (b) subgrain structures (blue) visible on LASM (Large Area Scanning Macroscope) data.Black lines indicate grain boundaries; the red outlines highlight the tilted-lattice band grains.

Figure 5 .
Figure 5. Line scan images from the same sample as shown in Fig. 3 from three focal depths with one highlighted layer (a) close to the surface, (b) in the centre of the core section, (c) close to the lower surface.(d) Sketch of the core sections; the upper part represents the physical properties sample, from which the thin sections are prepared.(e) Change of shape of the highlighted layer for the different foci.

Figure 6 .
Figure 6.Comparison of fabric data and visual stratigraphy in detail, bag 3596, approximate depth 1977.8 m.(a) Fabric data in a vertical section, (b) line scan image, in a vertical section, and (c) stereo plot of c axes orientations (horizontal plane).

Figure 7 .
Figure 7.Comparison of entire 55 cm core sections (full bags) of line scan (a, d) and fabric data (b, c).

Figure 8 .
Figure 8. ELLE model results for the simple shear experiment.Panels (a)-(c) show c axes orientations for shear strains of γ = 1, γ = 2, and γ = 3. Panels (d)-(f) show the distortion of the passive grid marker.Panels (g)-(i) show the equivalent von Mises strain rate field.

Figure 9 .
Figure 9. Basic folding mechanisms discussed in the text.(a) Passive folds form by shearing of disturbances in layering, without an active mechanical influence of that layering.Fold geometry is that of similar folds.(b) Buckle folds form by shortening of alternating strong and weak layers, in which the strong layers buckle and weak material flows into fold hinges.Fold geometry is that of parallel folds.(c) Tilted-lattice bands form in case of strong intrinsic anisotropy, but do not require viscosity contrasts between layers.