Physically-based data assimilation

Introduction Conclusions References


Introduction
Data assimilation deals with the optimal combination of observations and a model forecast, or background field, into an analysis field that forms the basis for the next forecast (Daley, 1992).Originating in meteorology, data assimilation is now Correspondence to: G. Levy (gad@nwra.com)used extensively with all operational geophysical models to improve model predictions and performance based on available observations.Most modern data assimilation techniques fall into two main categories: empirical methods, and methods based on statistical estimation theory (Talagard, 1997).Empirical methods, like dynamic relaxation (a.k.a.nudging, Hoke and Anthes, 1976) are most useful with new data and data sets for which error estimates and/or the error covariance structure are not known or available, as is the case for our application.Although implemented differently, it can be shown that all assimilation methods based on statistical estimation theory (hereafter, statistical data assimilation) can be regarded as an extension of optimal interpolation (OI) and are mathematically based on optimization algorithms, most commonly the least square minimization principle.The equivalency of statistical data assimilation methods, including OI, Bayesian data assimilation, Kalman Filtering (KF) and variational techniques, is shown, for example, in Kalnay (2003).Current and future data assimilation systems must cope with some known issues that stem from constraints, limitations, and errors in both the models and the observations assimilated, namely: i.During initialization, if the analyzed field, based on the data, does not match a realizable model state, noise is generated when the model integrates forward in time and this noise can severely impair forecast skill, and may lead to "the rejection problem" (Daley, 1992).
ii. Assimilation systems do not have a simple, quantitative way of representing lower dimensional features contained within a bulk simulation as these features are not G. Levy et al.: Physically-based data assimilation directly defined by either observations or model variables.Assimilation methods that ignore these features often force observed data onto models in a non-physical way.Such features occur frequently in geophysical applications and are associated with discontinuities and other important physical and dynamical processes on multiple scales.
iii.The observed data that are assimilated do not match the variables predicted by the model, requiring parameterizations that are not consistent with the data or the physics.
These issues arise when the validation and assimilation schemes do not maintain the physical principles embodied in the model and are unable to evaluate and assimilate lower dimensional features (e.g., discontinuities) contained within a bulk simulation that are not directly observed or represented by model variables.Under these circumstances, assimilation can lead to the violation of physical principles and the loss of information contained in the lower dimensional features.Conversely, models that resolve such features and the associated physics well, yet imprecisely, are often penalized by traditional schemes, leading to (perceived or real) poor model performance and forecasting skill scores.This loss of information can become deleterious in model improvements when observations are sparse, fuzzy, or irregular.The aforementioned issues have been reported in different applications of data assimilation.For example, in an assimilation experiment with a sea-ice model, Lindsay and Zhang (2006) illustrate that some fields in the model no longer strictly adhere to the physical principles of the model when data assimilation is accomplished through nudging (dynamic relaxation).They also note that similar inconsistencies arise when variables that originate from independent datasets are assimilated independently.Dai et al. (2006), in another seaice model simulation using optimal interpolation, show that the assimilated data degrade the solution at later stages because the underlying physical assumptions in the model are compromised.An ensemble Kalman filter (EnKF) assimilation in a wildfire model by Mandel et al. (2008) resulted in nonphysical states especially far away from the data.When discontinuous processes are modeled the problems of initialization can be exacerbated as shown by Vukicevic and Bao (1998) in a meteorological variational (4DVAR) data assimilation system.In that study, the linearization errors associated with discontinuous convective parameterizations were non-negligible and affected the assimilation results both locally and globally.
This paper considers a new paradigm for assessing model performance and for comparing model results with observational data.The resultant assimilation and validation scheme is compatible with state of the science methods and is capable of handling lower dimensional features in a bulk simulation.This scheme addresses issues (i), (ii), and (iii) above.The underlying principle for our assimilation algorithm is to use data to precondition model variables or material properties that allow internal model physics to generate future model states that agree better with observed data.The new scheme is tested with a sea ice model (Schreyer et al., 2006;described in Sect. 2.3 and Appendix A) and RADARSAT Geophysical Processor System (RGPS; Kwok, 1998) data, in a physically-based nudging context for one assimilation cycle.We employ a fuzzy verification metric of Levy et al. (2008) and standard skill scores of Murphy (1988) to evaluate model performance.We show that implementing the new assimilation scheme introduces model skill (against persistence) several days earlier than in the control run, improves the overall model skill, and delays its drop off at later stages of the simulation.We conclude with some thoughts about extending this method to a full assimilation cycle and to statisticalestimation data assimilation systems.

Data assimilation algorithm
Standard practice in data assimilation is to subtract model estimates of state variables interpolated to observation locations, from observed values (the background field), to produce innovations (or observation increments).Objective analysis of these innovations onto the model grid is then used for initialization of the model.The algorithm for validation and assimilation proposed and tested here assumes an adequate measure for model validation and verification exists, even for lower dimensional features.In it, feature extraction and assessment (see Sects.2.2 and 2.4) take the place of objective analysis for lower dimensional features.We use the fuzzy verification metric, RI, of Levy et al. (2008), Eq. (B1), as a benchmark for assessment.
An underlying principle for our assimilation algorithm is to use data to precondition model variables or material properties that allow internal model physics to generate future model states that agree better with observed data.That is, the model state is updated, as needed, through a physical innovation based on decision criteria from a fuzzy verification (e.g., Appendix B).Similar to four-dimensional (4-D) data assimilation, this procedure produces results that are consistent with internal model physics and dynamics and thus avoids forcing unrealizable model states.Furthermore, in theory, this scheme should work equally well for resolving and assimilating lower dimensional features without the need to transform them into model variables.This is important because in most models the lower dimensional information is not resolved or used.Computationally, this scheme is significantly more efficient than 4-D assimilation and could be regarded as a simplified 4-D assimilation, especially suited for incorporating lower dimensional information.

Lower-dimensional features in sea ice observational data
We choose to test the new assimilation scheme in a sea-ice environment where lower dimensional features in the form of Linear Kinematic Features (LKFs; Kwok, 2001) are abundant and observable.LKFs represent discontinuities in the sea ice (e.g., leads or ridges, Coon et al., 2007).Although not directly measured or resolved by observing systems or most models, LKFs can be related to model resolvable physical or material properties.LKFs persist sufficiently to allow simplified testing of our scheme that adjusts model variables at initialization.A general assimilation scheme would call for the model state to be updated regularly based on the data and their agreement with the model output, rather than just initially.The adjustments can be implemented at set intervals (i.e., in a continuous data assimilation), or only when agreement or skill score falls below a certain threshold.This validation and assimilation technique is schematically depicted in Fig. 1 for both the general case and the ice model implementation tested.The data we use come from the RADARSAT Geophysical Processor System (RGPS), which was developed by the Polar Remote Sensing Group at the Jet Propulsion Laboratory (JPL) to extract sea ice motion data from SAR imagery (Kwok et al., 1990).At an initial time, a set of points forming a regular grid is located in the SAR data sets.Then, in the images resulting from subsequent satellite passes (approximately every 3-days), the original points are found again using area-based and feature-based tracking.This procedure provides displacements for each point.If the set of points in the original configuration is viewed as the vertices of square cells, then the motion of the points determines the deformation of the cells.With this interpretation, grid quantities such as the divergence, shear, and vorticity, can be calculated using the nodal displacements.
The RGPS deformation products are based on the assumption that the displacements and velocities are smooth functions of the spatial coordinates.However, if the dominant form of deformation of multiyear ice is in the opening, closing, and shearing of linear features or leads, then the displacements and velocities can be discontinuous.In Coon et al. (2007) we discuss the kinematics associated with strong discontinuities that describe possible jumps in displacement or velocity.Specifically, we determine a jump in displacement and the orientation of this jump that account best for the observed deformation of an RGPS cell.We use this treatment of the data for feature extraction and subsequent assimilation of the information.

Sea-ice model and simulations
The sea-ice simulations we perform to demonstrate our scheme use an elastic-decohesive constitutive model for the sea ice (Schreyer et al., 2006).The model was developed for predicting the initiation and opening of leads in the Arctic ice.Once the existence of leads is taken into account, the remaining motion of the ice has small deformations and is appropriately described as elastic.Several features were designed into the model.First, the model was constructed to transition from observed brittle failure under tension, to compressive brittle failure under moderate compression, and to a plastic-like faulting under large confinement (Schulson, 2004).The various modes of failure occur in the model, depending on the stress state in the material.Where the transitions occur in stress space depends on the material parameters and can be adjusted based on empirical data.Second, the model can handle multiple cracks at a point, and therefore can predict crack branching.Third, the numerical implementation of the model is accomplished similarly to standard plasticity models.Thus, in principle, modular codes that call a subroutine to implement the constitutive model can substitute the elastic-decohesion model if it proves worthwhile.A final aspect of the model is the ability to build in pre-existing planes of weakness that may be due to pre-existing, partially frozen leads, for example.It is this feature that we exploit in order to initialize simulations using observed data.More information about this model can be found in Appendix A.
The equation of motion for the ice is the balance of momentum equation that includes, in addition to the internal forces determined by the constitutive model, drag forces from the wind and ocean currents, and Coriolis forces.Six www.geosci-model-dev.net/3/669/2010/Geosci.Model Dev., 3, 669-677, 2010 hour wind fields from NCEP reanalysis are used to determine the wind drag and the ocean currents are updated daily using output from an ocean model (MITgcm, Marshall et al., 1997) run independently from the ice simulation.The momentum equation is solved using the materialpoint method (Sulsky et al., 2007).With the material-point method (MPM), a set of material points is identified in the body of fluid or solid that is tracked throughout the deformation process.Each material point has a mass, position, velocity and stress, as well as material parameters and internal variables as needed for constitutive models or thermodynamics.These material points provide a Lagrangian description of the material that is not subject to mesh tangling because no connectivity is assumed between the points.This Lagrangian frame naturally models convection and transport since the trajectory and history of each material point is followed.Each point carries material properties without error, and history variables can be integrated along the trajectory.However, computing gradients for solution of the momentum equation is complicated in this representation since the neighbors of a given point are not known a priori, and can change during a simulation.To keep the computational work linear in the number of material points, a second discretization is used for solving the momentum equations.This representation of the solution is often a regular, background mesh that covers the computational domain.
Run explicitly, the time step in MPM is governed by the CFL condition based on the background mesh size and the elastic wave speed.This step size is comparable to the step size used in sub-cycling the elastic-viscous-plastic model (Hunke and Lipscomb, 2004), making MPM with the elasticdecohesive model competitive in terms of computational efficiency with the best available algorithms for ice dynamics (Sulsky et al., 2007).

Assessment
A core idea in our validation scheme is the ability to assess model performance using fuzzy verification models.Consider, for example, the importance of resolving lower dimensional features in sea-ice for climate models.In the winter, the leads or discontinuities in sea ice are regions where the relatively warmer ocean is exposed to the cold atmosphere.New ice forms rapidly in these regions and as it forms, brine is ejected into the upper ocean.This denser water sinks and contributes to the conveyor-belt circulation.If forces change and leads close, the thinner ice in the leads can be forced up into ridges or down into keels, increasing the amount of ice that can be stored in a given area.Thus, leads and ridges (the linear kinematic features in the model and observations) are crucial for global climate models.However, their exact location is not important.Thus, our climate predictions would probably be fairly accurate if we predict leads of approximately the right size and in roughly the right location.Thus we do not want a metric that compares model and observation point-wise, rather we want a metric that tells us the simulation is satisfactory if we get the leads roughly right.That is what our assessment is meant to accomplish.Levy et al. (2008) define two metrics to evaluate model success in representing lower dimensional features.They treat features through a frequency distribution at predetermined spatial regions of the domain.We use one, the RMS index of agreement, in a general skill score (Murphy, 1988): Where A is a measure of accuracy and subscripts f, r, and p denote forecast, reference, and perfect, respectively.The reference state used here is that of persistence.The measure of accuracy we consider -the RMS Index of Levy et al. (2008, see Appendix B) -can take a value between 0 (no agreement) and 1 (perfect agreement, A p ).As we deal with lower dimensional features for which no climatology (and hence no error correlation information) exists, persistence serves as the reference.Any skill score value greater than zero indicates prediction skill over persistence, and can be directly converted to a percentage improvement in skill.

Experimental design
We conduct two simulations with the sea-ice model: a control run, and an experimental run.The control and experimental simulations are of ice behavior in a 831.600 km 2 region of the Beaufort Sea over a time interval of 16 days from day 54 (23 February, insert in Fig. 2) through day 70 (11 March) of 2004.The background mesh size in the MPM calculation is 10 km, with four material points per cell initially.Model parameters include the elastic moduli, ice density, Coriolis parameter, drag coefcients and the parameters in the decohesive model that set the failure strength of ice in tension τ nf and shear τ sf , the length scale over which decohesion occurs u 0 , and two parameters that describe the shape of the failure surface in stress space, κ and τ sm (Schreyer et al., 2006, and Appendix A).Values of the parameters used for the simulations are given in Table 1.During this same period in 2004, RADARSAT SAR observations processed through the RADARSAT Geophysical Processor System at 10-km resolution are available daily for validation (middle column in Fig. 3) and assimilation (insert in Fig. 2).
The simulations use the elastic-decohesive constitutive model described above.Six hour wind fields from NCEP reanalysis are used to determine the wind drag and the ocean currents and are updated daily.Both the control and the experimental runs use the same boundary conditions.Boundaries are either land boundaries, in which case the displacement of the ice is zero, or open boundaries in the Beaufort Sea.In the latter case, displacements of boundary points are specified according to their observed values obtained from the RGPS data.
Initial values of all model variables and parameters are also the same, except that the experimental simulation goes through a one-step assimilation during day 54, where information on the pattern of leads from RGPS at 10-km resolution is used to predispose the model.Initialization includes setting the initial velocity and stress to zero, and the initial ice thickness uniformly to 3 m.The experimental simulation is further initialized as determined through a kinematic analysis of the RGPS data over one day (Coon et al., 2007;Peterson and Sulsky, 2011).Specifically, we determine a jump in displacement and the orientation of this jump that account best for the observed deformation of an RGPS cell.Material points within that cell are initialized to have this initial jump at the observed orientation.The effect of this initialization is to reduce the strength of the ice anisotropically, and to predispose it to continue deforming with this oriented opening provided the forcing is consistent with this deformation.
We run accuracy assessment daily and consider the model response and the impact of the one-step assimilation on simulating LKFs.The assessment consists of determining the accuracies A f and A r against the daily observations for both the experimental and control runs using the RMS index of Levy et al. (2008, see Appendix B), and substituting them in the general skill score formulation, Eq. ( 1) above.We score the different simulations using the RMS Index against the LKFs interpreted from the RGPS observations of the Beaufort Sea.We consider the agreement in (1) the existence of 100 km×100 km grid cells containing LKFs in the domain; and (2) the existence of 100 km×100 km grid cells containing LKFs at the observed orientation using four cardinal orientations in the domain.We thus score M=5 features at N=1, the entire simulated/observed domain (the Beaufort Sea), in Eq. ( B1).Here we assume that all weights, w i and ω i equal one.Frequencies are defined in terms of cell counts, and cells with missing data are excluded from the statistics.The impact of the one-step assimilation on simulating LKFs is shown in terms of a skill score in Fig. 2, and visually in Fig. 3.

Results
Figure 2 shows the evolution of the skill score with respect to persistence for both the control and experimental simulations for the duration of the simulations from day 55 (24 February 2004), the first day after assimilation of RGPS data (in figure insert) in the experimental run, through day 70 (11 March 2004), the last day of the simulations.The impact of the one step assimilation is evident throughout the entire 16day simulation, although it is most dramatic during the first week.For visual assessment of the impact that the assimilation has on the LKF field throughout the domain, Fig. 3 provides snapshots of that field as simulated in the control and experimental runs, side by side with the observed field at key time steps of the 16-day simulation.
As the lower dimensional information assimilated in our tested assimilation scheme is used to nudge the model through preconditioning of the simulated field towards the observed field at a future time step, a positive skill score with respect to persistence is achieved only on 25 February (top panel in Fig. 3; day 56 in Fig. 2).This reflects the time it takes for the model physics to consistently respond to the www.geosci-model-dev.net/3/669/2010/Geosci.Model Dev., 3, 669-677, 2010 The RGPS data were processed assuming that all deformation should be accounted for by shearing, opening and closing of a discontinuity, which passes through the cell center (Coon et al., 2007).
preconditioning, as well as the relatively strong persistence and slow evolution of the lead (LKF) field.The skill score of the experimental run peaks at 0.6 (60% improvement over persistence; 100% improvement over the control run score) on day 58 (27 February), and remains positive for the rest of the simulation, and above 0.5 through day 66 (7 March).
In an adaptive assimilation cycle, one may consider another cycle on day 67.The true impact of the assimilation on the model performance is measured relative to the control run.
One day following the assimilation, on day 55 (24 February), the skill score of the experimental simulation is over 100% higher than that of the control run skill score, even prior to either simulation showing skill relative to persistence.At its peak on day 58, the experimental run skill score exhibits an improvement of 100% over the control run's skill.
The control run gains skill relative to persistence on day 59 (28 February; Fig. 3), when it is 80% lower than the skill of the experimental run.The control run skill relative to persistence remains positive for the rest of the simulation, peaks at 0.35 on days 63 (4 March in Fig. 3) and 64, but never exceeds the score of 0.5, and, as the skill of both simulations slowly decline at the end of the simulation (e.g., 10 March at the bottom of Fig. 3), it remains consistently at least 10% lower than that of the experimental simulation.

Conclusions
An algorithm for model validation and data assimilation that maintains the physical principles embodied in the model and can evaluate and assimilate lower dimensional features (e.g., discontinuities) contained within a bulk simulation is introduced and demonstrated with a sea-ice model and with remotely sensed observations of leads in a one step assimilation cycle.An underlying principle of the new algorithm is to use data as a guideline for the model by preconditioning model variables or material properties in a way that leads internal model physics to generate future model states that are in better agreement with the observed state.Similar to four dimensional data assimilation, this procedure produces results that are consistent with internal model physics and dynamics and thus avoids forcing unrealizable model states, largely resolving the problem of initialization.Furthermore, as tested here, this scheme works well for resolving and assimilating lower dimensional features, which do not match the model variables and are not directly measured by the observing system, without the need to transform these features into model variables.Computationally, this scheme is significantly more efficient than four-dimensional assimilation and could be regarded as a simplified 4-D assimilation.
We have tested the new scheme in a sixteen day simulation experiment.In this system, model skill (against persistence) is initiated three days earlier than the skill in the control run with no assimilation, and is consistently higher than the later throughout the simulation.In addition to this shortening of the model spin up time and the overall higher model skill, the assimilation improves model skill significantly during the first six days.Thus, the new scheme holds the potential of comparable improvements for the assimilation of lower dimensional features with similar persistence when infrequent observations exist.
In the sea-ice system tested, the properties needing physical adjustments are relatively clear and nudging is a natural implementation choice that is capable of a robust response to the adjustments.However, the same principle of physically based adjustments holds in the general case of other geophysical models and systems, the variables and lower dimensional features they resolve or represent, as well as in other modes of implementation.Meteorological examples include forecasting precipitation and tropical cyclone trajectories (e.g., Kurihara and Ross, 1993).Thus, this scheme could be extended to a generalized multivariate data assimilation and fuzzy verification system that would be physically based and capable of extracting lower dimensional information from observations and bulk simulations at different scales and for different geophysical (e.g., atmospheric, oceanic, coupled) models containing continuous and lower level features of significance (e.g., fronts, organized convection).

Appendix A
For the purposes of modeling sea ice, a 2-D, plane-stress description of failure has been formulated assuming cracks occur in the plane.The envelope of failure points in stress space is described by a failure function, F n (σ,n) where F n <0 implies no failure, F n =0 implies evolving failure and F n >0 is not allowed.This function is analogous to a plastic yield function in plasticity theories.The subscript n on F n indicates a separate failure function for each potential crack orientation, and F n depends on the stress σ , and the unit normal n to the crack surface.To consider all possible failure directions, a general failure function F is defined as F = max n F n .
Many classical failure criteria, such as the Rankine, Tresca and Mohr-Coulomb criteria, are expressed in terms of the traction on the failure surface (i.e., crack surface).The elastic-decohesion model extends these classic criteria by adding two new features: (1) a modification of the Rankine criterion for brittle failure to allow for the possibility that a compressive stress component may lower the resistance of the material to brittle failure, and (2) a transition from brittle to ductile failure within one criterion.If a local basis consisting of n, the unit normal to the crack, and t, a unit vector tangent to the crack, is introduced, then the traction on the failure surface has normal component Material parameters are τ nf , the tensile or normal failure stress and f c which denotes the failure stress in uniaxial compression.For the moment, take f n = 1.The new criterion for brittle failure is B n = 0.The McCauley bracket ( • ) is used to activate the normal component of stress σ tt only if it is negative.If the term involving σ tt were absent then failure would occur when the normal traction on the surface reaches the threshold τ nf , which is the Rankine criterion.With the σ tt term, this criterion is analogous to the Rankine criterion in that failure occurs in the direction of maximum principal stress, but the critical value of the normal traction component is potentially reduced when σ tt is compressive.The criterion allows for failure even if τ n is negative, and it is this aspect of the model that allows compressive brittle failure.Next, brittle and ductile aspects of failure are included by defining the failure function as The additional material parameter, τ sm , is the failure stress in shear when the material is under large compression (τ n → −∞).The parameter, κ, is derived from the condition that pure shear (B n = −1), the failure stress is τ sf , the failure stress under pure shear.If τ nf → ∞ and f c → ∞ then the criterion F = 0 reduces to the pure shear criterion of Tresca.Figure A shows a sketch of the decohesion failure envelope in stress space.The solid line represents the failure envelope F = 0. Along this solid line, the blue arrows indicate the direction of maximum principal stress and the red arrows indicate the normal to the crack surface.Under brittle failure the normal to the crack is in the direction of maximum principal stress.Under ductile and mixed-mode failure the normal to the crack is at an angle to the direction of maximum principal stress, with two orientations of the crack possible.Of the two, the orientation that preserves the sense of local rotation is chosen.The transition from brittle to ductile failure occurs at a point along the failure envelope determined by the ratio of τ nf to τ sf , and thus is a material property.
This failure envelope describes the model for lead initiation in the ice.Once the beginning of a crack has been identified, the evolution of the lead is required.The term decohesion or cohesive crack model refers to the reduction of the traction on the crack as the crack opens.This behavior is in contrast to Griffith's model (Griffith, 1921) where the traction is assumed to change discontinuously and instantaneously from a positive value to zero.Decohesion is included in the model by introducing a softening parameter, analogous to equivalent plastic strain in plasticity models, that drives the traction to zero as a crack continues to open.A dimensionless parameter, f n in Eq. (A1), starts with a value unity for undamaged material and reduces to zero as u n , the normal component of the jump in displacement, increases from zero.The crack is considered completely open when u n reaches the material-dependent value u 0 , at which point the traction on the crack surface has been reduced to zero and a free surface is thus formed.
The displacement discontinuity evolves according to a normal flow rule un = ω ∂F ∂τ n ut = ω ∂F ∂τ t . (A3) The displacement discontinuity is regularized into an effective decohesion strain, analogous to plastic strain, where L is a measure of the cell size in numerical simulations.(The value of L is chosen so that the physically correct energy is dissipated during fracture.)The stress is a function of the elastic strain e − e d .Thus, as a specimen of ice is loaded, we typically begin with F < 0; the stress is inside the failure envelope.We assume each loading step is elastic, giving a trial stress state.If the trial stress is outside the failure envelope (F > 0) then a jump in displacement is introduced to bring F back to zero.This procedure is identical to standard solution procedures for plasticity.The result is that as a crack opens we predict the amount of both the normal and tangential opening.Once a free surface has formed, the jump in displacement can continue to grow if the crack surfaces continue to separate, and the traction on the surface remains zero.
At each loading step we find the critical direction n for which F is largest.As a crack with a particular orientation begins to open, the softening makes it likely that this orientation will remain the critical direction.However, it is possible that a changing stress state will make another direction critical, in which case a second crack can form intersecting the first.In this manner, the model accommodates multiple cracks at a point.If weak areas are known to exist in the ice, the softening parameter f n can be initialized with a value less than one to account for this information.

Appendix B
The RMS index of agreement (Levy et al., 2008) we use as a measure of accuracy in the skill score, Eq. ( 1), treats features through a frequency distribution at predetermined spatial segments of the domain.It contains a term with the familiar format of standard error common in routine distance or root mean square error measures used for continuous variables: where w i and ω j are weights given to the features and the spatial segments, respectively, and D i,j is a normalized frequency difference function: where p i,j and o i,j are the predicted (simulated) and observed feature frequencies (cell count) of feature i in segment j .

Fig. 1 .
Fig. 1.Schematic of the validation and assimilation algorithm illustrated for the general case and (in parentheses) lower dimensional features and the sea ice model and data tested.The implementation and testing on lower dimensional features involve the assimilation of RGPS data processed to show regions of high deformation.The model state may be updated by changing material properties based on agreement of lower dimensional features and fuzzy metrics scores.For example, if the data indicates the presence/absence of LKFs in a region, a jump in displacement is determined to account for the observed RGPS deformation of the cell.

Fig. 2 .
Fig. 2. Skill score (against persistence of observed field shown in insert) evolution of the control and experimental simulations.A 1step DA algorithm was implemented on day 54 (23 February 2004) to the experimental run, whereby the ice model state was updated by changing material properties based on agreement of lower dimensional features deduced from RGPS data processed to show regions of high deformation (insert).Days are Julian days of 2004.

Fig. 3 .
Fig. 3. Comparisons of the experimental (left column) and control (right column) simulations with LKFs interpreted from the RADARSAT Geophysical Processor System (RGPS) data (center column) on (from top to bottom) 25 February 2004 (day 56, first day of positive skill score for experimental run), 28 February (day 59, first day of positive skill score for control run), 4 March (day 63) and 10 March (day 69).The RGPS data were processed assuming that all deformation should be accounted for by shearing, opening and closing of a discontinuity, which passes through the cell center(Coon et al., 2007).

Fig. A1 .
Fig. A1.Failure envelope in principal stress space for the elastic-decohesive model.