Non-proportional loading in sequentially linear solution procedures for quasi-brittle fracture_ A comparison and perspective on the mechanism of stress redistribution

Sequentially linear solution procedures provide a robust alternative to their traditional incremental-iterative counterparts for finite element simulation of quasi-brittle materials. Sequentially linear analysis (SLA), one such non-incremental (total) approach, has been extended to nonproportional loading situations in the past few years. Although the process of damage propagation and localisation is often dynamic in nature, the simulation being quasi-static poses a fundamental problem. This article gives an overview of the different approaches to address nonproportional loading in SLA and other sequentially linear methods, and their corresponding redistribution methodologies to address the dynamic phenomenon. Furthermore, the inherent differences between two such methods: SLA (total) and the Force-Release method (incremental), and their suitability to structural continuum models involving non-proportional loading, are illustrated using real-life concrete and masonry experimental benchmarks tested up to and beyond brittle collapse. In each illustration, SLA is shown to enforce equilibrium during dynamic failure by load reduction, using the intermittent proportional loading, while allowing for active damage propagation resulting in a relaxed failure mechanism which manifests as snap-back(s). Contrarily, the Force-Release method is shown to describe the collapse through states of disequilibrium.


Introduction
Finite element (FE) models, in use for simulation of materials characterised by brittle failure like concrete and masonry, often encounter problems related to snap back, bifurcation points, divergence or material softening leading to negative tangent stiffness and the consequent ill conditioning of the formulation. Convergence issues, resulting from multiple integration points being pushed simultaneously into strain softening, may lead to inaccurate results due to deviation of the simulation into alternative equilibrium paths. Use of path following procedures like the Arc-length control method help but requires extreme care from the user to steer the analysis (and possibly prior knowledge of crack locations). This gave rise to several alternate solution approaches which can generally be classed into three categories: purely total approaches wherein unloading and reloading is done non-proportionally, purely incremental approaches wherein the stress and loading history is explicitly tracked and finally a class of combined incremental-total approaches. The Sequentially Linear Analysis (SLA) is a feature, as a part or whole, of all three aforementioned categories.
known as the saw-tooth law. In principle, the material law is described as a series of successively reducing secant stiffnesses, starting from the initial elastic branch with the original Young's modulus of the material (E 0 ). Whenever there is breach of the stress limit in an integration point and it becomes the most critical in the FE model, the next secant relation with reduced strength and stiffness properties takes over from its previous secant branch. This process of reducing the stiffness upon attaining a stress limit in a single integration point is repeated until the stiffness of the structure has vanished, which corresponds to a state of complete damage. An overview of the approaches used thus far for saw-tooth approximations of typical tension and compression softening curves can be found in Reference [5]. Fig. 1(a) depicts the ripple bandwidth type saw-tooth law for a linear tension softening relation.
These uni-axial saw-tooth laws are used in the orthotropic fixed smeared crack framework, wherein as soon as the principal stress violates the allowable strength at an integration point, the isotropic stress strain relation = D transforms into an orthotropic relation as nst = D nst nst with the nst cracked coordinate system. The primary principal stress direction's Young's modulus and strength are damaged according to the uni-axial saw-tooth law. In the event that normal stresses in the tangential directions (secondary or tertiary) violate the corresponding allowable strengths, caused by stress rotations or redistribution of stresses or application of another load non-proportionally, damage is introduced in those directions similarly. So every integration point essentially requires two uni-axial saw-tooth laws each for tension and compression in the 2-D stress state and 3 in the case of a 3-D stress state. This aside, the shear behaviour in the fixed cracking model is represented using a variable stepwise shear retention function that takes into account the reduction of shear stiffness with increasing damage in the corresponding normal direction of the cracked plane [13]. Also, the Poisson's ratio is reduced at the same rate as the associated Young's modulus. For further information on the smeared cracking constitutive model used in SLA, the reader is referred to References [3,5,14,15]. Furthermore, SLA has been extended to other constitutive models as well like discrete cracking, bond slip and coulomb friction models for interface elements. For detailed information on these models, the reader is referred to Reference [5].
The studies in this work involve masonry and concrete examples. Masonry can be represented in FE models using micro, meso and macro modelling approaches depending on the required level of accuracy (local or global behavior of the structure) [16]. Macromodelling, which considers masonry as a homogeneous continuum, is used predominantly in this work because it is advantageous owing to reduced computational time and memory requirements as well as the user-friendly mesh generation. The material behaviour of the continuum is then described using continuum damage or plasticity-based constitutive laws. The masonry examples in this study are represented using a combination of the aforementioned fixed smeared crack/crush framework and discrete cracking interfaces, where necessary. A simple Rankine-type failure surface is used for both crack and crush initiation in the continuum. This is followed by linear tension softening and parabolic compression softening [17] relations for the orthogonal directions of fixed crack/crush system. The effects of lateral confinement or cracking on the compressive behaviour are currently under investigation for the SLA framework and therefore, discarded in this study. Concrete has also been represented using these formulations quite adequately thus far [18,19] and is therefore used in this study.

Work flow -Proportional loading
Firstly the saw-tooth laws are defined as shown in Fig. 1(a). Thereafter, the finite element model is loaded by a unit value of the total imposed load (L) which could either be forces or prescribed displacements (or a combination thereof). Subsequently, a linear analysis is performed and a load factor can be calculated for each integration point as the ratio of the allowable strength to that of the governing stress (considering the principal stress theory for a smeared fixed crack approach) as shown in the following.  where j refers to the analysis step and i either denotes an integration point number for an undamaged situation or alternatively upon damage, denotes events corresponding to tension/compression failure criteria along the 2 or 3 fixed damage directions, depending on the stress state.
i j gov, is the governing stress component for integration point i f , i j is the peak stress limit as defined by the current secant branch of the saw-tooth law and i j crit, the associated load multiplier. To ensure that only one integration point reaches its peak stress limit, the linear analysis is scaled with the minimum of all load factors, which is referred to as the critical load factor. The total factor and the applied load are then defined as The strength and stiffness of this critical integration point is then reduced based on the saw-tooth laws. Eventually, the results of the linear analysis; the stresses, strains and displacements, are scaled using the critical load multiplier.

Existing approaches to non-proportional loading in sequentially linear solution methods
The total framework in an SLA type of method works well under proportional loading conditions, i.e. when loads acting on a structure increase or decrease simultaneously. When there are non-proportional loads, some of which are constant over the structure unlike loads which vary over time, problems arise due to considerable stress rotations and consequently finding the critical load multiplier and integration point was not as straightforward anymore. The extension to non-proportional loading was intended to make the procedure more suitable to wider applications and the various strategies that have been proposed to incorporate this aspect are presented in this section. Additionally, all the incremental and combined incremental-total alternatives are reviewed in this section. The sequentially linear class of methods have alternatively been referred to as non-iterative methods in literature [21]. However, considering that SLA is a feature of all the methods, this class is referred to as sequentially linear methods hereon.

Total approaches
This class of sequentially linear methods relies on unloading and reloading a specimen, of all previously applied loads, while undergoing damage propagation. In line with this nature, they are also alternatively referred to as load-unload methods in literature [10,22]. Each of them trace the critical event differently but are all similar with respect to non-proportional loading.

DeJong's method
The method of DeJong et al. [7] was proposed for the smeared cracking approach, involving tensile failure only, under plane stress assumptions. For a system loaded by constant loads (L con ) and a unit variable load (L var ), the global stresses are expressed as the superposition of the stresses due to constant and the scaled variable loads ( con and var respectively) as shown in Eq. (3). These are then substituted in the expression for principal stress, to be limited by the tensile strength ( f ten ), as a function of the load multiplier for the variable load ( ).

= + con var
(3) The closed form solution for of Eq. (4) corresponds to crack or damage initiation for an integration point. For already damaged integration points, linear equations are solved along the orthogonal directions of the crack coordinate system, n and s, as shown in Eq. (5).
The selection criteria for the critical load multiplier are based on the stress caused by the variable load. In case the variable load causes tension, considering the crack opening effect, the associated t is considered a maximum load multiplier (in which the superscript j indicating the analysis step is dropped for readability). In case the variable load causes compression, considering the crack closing effect, the associated c is considered as a minimum load multiplier. The minimum of all t from the integration points is denoted min t and the maximum of all c from the integration points is denoted max c . In the scenario that > , min t max c min t is chosen as crit and the corresponding integration point is damaged based on the saw-tooth law ensuring that all other integration points adhere to the constitutive law. However, when < min t max c , no crit can be selected such that all integration points in the finite element model adhere to their respective constitutive relations. This situation hereon is referred to as limit point, wherein no constitutively admissible 3 crit can be found. This limit point is not to be confused with the traditional limit point encountered in an NLFEA response. DeJong's approach in such a situation proceeds forward by selecting max c as the critical load multiplier. Therefore, it allows for temporary violation of the constitutive law in one or more integration points, which may lead to rupture under invalid stress fields as has been pointed out in References [5,10].

Constrained maximisation analogy
The constrained maximisation approach proposed by Van de Graaf [5] follows the basics of the non-proportional loading concepts proposed by DeJong et al., with regards to the principal of superposition of stresses, the use of principal stresses for damage initiation M. Pari, et al. Engineering Fracture Mechanics 230 (2020) 106960 and calculation of the closed form solution for the critical load multiplier. The difference lies in the selection criteria for the critical load multiplier which is based on a constrained maximisation approach, wherein sets of admissible load multipliers are deduced. Damage initiation involves solving quadratic (2D) or cubic inequalities (3D) in load multiplier as the stress state may be [14], while damage propagation involves solving linear inequalities per direction of the orthogonal cracked system. The common global subset of these local sets would then reflect the set of constitutively admissible load multipliers for the entire FE model, the maximum of which is chosen as the critical load multiplier.
In the event that this common subset is empty in a certain analysis step j, the approach runs into the aforementioned limit point. In order to adhere to the constitutive law, Van de Graaf proposed the double load multiplier strategy which includes an intermittent proportional loading (IPL), wherein the last successful load combination L ipl of step ( j 1) (Eq. (6)) is scaled proportionally as shown in Eq. (7). It allows for redistribution by implicitly allowing for the reduction of the constant load and is deemed acceptable when the redistribution is followed by full recovery of the constant load. For further information on the work flow for the constrained maximisation analogy, the reader is referred to References [5,14]. Chenjie et al. [8] also proposed an algorithm for the non-proportional loading problem but confined their applications to problems wherein the determination of the critical load multiplier was not affected significantly by the constant loads. This class of problems was referred to as weakly nonlinear problems.

Improved total analysis
The inherent feature of unloading and reloading non-proportionally is the main cause for the problems associated with nonproportional loading and to supposedly overcome this, Alfaiate et al. proposed the Improved total analysis [9]. Considering a nonproportional loading case of loads L 1 and L 2 applied sequentially, the improved total analysis proposes an initial load close to the last successful applied load instead of reloading again sequentially as L 1 and L 2 . So at analysis step j the load is determined based on the preceding load combination as follows.
The term j could be negative or positive but the total second load L 2 applied remains positive. In principle, the approach claims to allow for partial unloading and thereby avoid reloading non-proportionally. However, it seems clear from a prestressed beam example in their work that the stress history is not kept in memory. The stress state corresponding to each load is indeed recomputed for the damaged state, the representation of using about 75% of the last successful load multiplier is only mathematical and the approach does not truly adhere to the physical loading history like in an incremental approach. The approach is therefore purely total and equivalent to the constrained maximisation analogy described above. Additionally, the approach does not indicate its take on the limit point situation where there is no possible load multiplier. However, simultaneously Alfaiate et al. [9] proposed an incremental version, allowing for partial unloading, called the secant-incremental analysis which is described in the following section.

Force-release method
Elias et al. [10] also pointed out that DeJong's approach essentially did not incorporate the correct loading history and additionally that an event may lead to a series of subsequent failures in the vicinity of a damaged element without increasing the prescribed displacement or forces. These regions normally appear as instabilities in a displacement controlled experiment which look smoothed out. However, in principle, these are snap backs which are referred to as avalanche of ruptures by Elias et al., driven by redistribution of the elastic energy released from damaged elements into their vicinity. Addressing this fundamental problem in simulating a dynamic phenomenon, like an avalanche of ruptures due to damage (cracking or crushing), using a quasi-static approach, Elias et al. [10] proposed the Force-Release method.
This method is a purely incremental approach wherein the non-proportional load path is discretised into a series of piece-wise proportional loading paths. They proposed to discretise each load into a series of n u load vectors with magnitudes L u , with = u n 1, 2, .. u , ensured to be non-decreasing, so that the proper loading/stress history is taken care of. For force loads, the magnitude is based on the euclidean norm of the prescribed forces. The variable t u , ranging from 0 to 1, represents how far the simulation has proceeded along a certain load magnitude L u and is essentially a fraction of the prescribed load. It is expressed as in Eq. (12), where L u are reference load increments analogous to the unit load in SLA. The stiffness matrix, displacements, reactions and elements stresses are denoted as K d R , , and . Linear analyses are performed either with the reference load increments, which would represent regions of static equilibrium in the response, or with disequilibrium forces representing the intermediate redistribution region. After each load increment or event, all quantities are updated with their incremental parts and stored.
The work flow is summarised as follows. To begin, linear analysis is performed with the reference load increment of the first load L 1 to compute the reactions R, displacements d and stresses . Thereafter the critical element, i, is deduced such that the following holds.
In Eq. (10), all quantities with are their corresponding incremental values caused by the load increment L 1 and f is the allowable strength. Once the critical element is evaluated, a decision between damage & redistribution or proceeding to the next load's increment L 2 is made based on the following (with appropriate values of the load sequence u).
If the inequalities of Eq. (11) are true, the critical element i is damaged. The stiffness matrix, reactions, displacements and stresses are updated with their corresponding scaled incremental quantities. The variable t u is also updated with its scaled increment as shown in Eq. (12). Subsequently, the system is loaded by disequilibrium forces, = S R Kd that arise due to the previous event. This region represents a state of disequilibrium and a new potential critical element is deduced based on Eq. (10). < 0 1 indicates that the unbalanced forces generated by the previous damage event i needs to be redistributed in its vicinity, by another rupture. The newly identified critical element is damaged, and the stiffness matrix, reactions, displacements and stresses are again updated. Thereafter, the system is loaded by the newly computed disequilibrium forces to determine if there needs to be yet another rupture due to the redistribution. Thereby, the avalanche of ruptures or the so-called redistribution loop continues, as long as < In other words, the system is in equilibrium and further external load is added. The internal state parameters are updated and the simulation moves to the next load + L u 1 similar to the continuation at the end of the avalanche. The Force-Release method, in summary, uses the constitutive model as in SLA performing linear analyses with reference load increments which are fractions of the total prescribed load. Each reference load increment corresponds to an analysis which may or may not lead to damage. Upon damage, the stress from a damaged element is released gradually through a sequentially linear redistribution loop wherein other elements may be damaged. The approach is suitable for non-proportional loading paths but cannot handle snap backs because of its incremental nature. Additionally, although it is claimed that these imposed loads could either be forces or prescribed displacements, for simulations involving global softening, the Force-Release procedure would work only under prescribed displacements since imposed forces can't reduce in magnitude.

General method
This approach is shown to be a generalisation, of which the load-unload and the Force-Release methods are extreme cases [12]. It also overcomes the problems of non-observance of snap back in the Force-Release method and the lack of avalanche of ruptures due to a single damage event in the load-unload method. The general method is based on the redistribution of stresses from damaged elements as is done in the Force-Release method and the simultaneous scaling of external load characteristic of the load-unload method.
The formulation involves defining two reference stress variables L and S corresponding to the external load increment and the disequilibrium forces. The method is suited to either force or displacement loading or both (depending upon appropriate boundary conditions). The total reference stress state is expressed as a sum of the two reference stress variables, enhanced by a slight modification of the external load by a parameter , as shown in Eq. (13). determines the ratio between the effects of the external load and the disequilibrium forces to the reference stress. Substituting this in the governing equation of the Force-Release method, Eq. (10), the expression for the general method is deduced (Eq. (14)).
The method now is general because the choice of the parameter allows the process to be steered as desired. By setting = 0, the external load is kept constant which would be equivalent to the Force-Release method. It can also be led to alternate equilibrium states by choosing values for arbitrarily. If the is set such that the redistribution finishes exactly when the external load diminishes to zero, the path would lead back towards the origin and this corresponds to the load-unload method. The choice of the parameter to control the external load is direct while indirect controls, like the crack mouth opening displacements (CMOD) which are often used to steer experiments, are also possible in this set-up of the general method. The approach in shown be suited to both the proportional and non-proportional loads [12].
Simultaneously, Liu and Sayed [22] proposed another general approach which also treats the load-unload and Force-Release methods as special cases. However, this is based on a simple linear interpolation of the trial displacement fields of both the methods as against the scaled combination of the reference stress states due to the external load and residual forces as in Elias' general approach [12]. Elias claims that the latter seems a more reasonable formulation considering the underlying physics of the problem since the displacement field interpolation lacks a strong physical motivation [12].

Secant-incremental analysis
In addition to the previously described Improved total analysis, a purely incremental procedure referred to as the secant incremental solution was proposed by Alfaiate et al. [9]. Inspired from the former, the latter approach allows for partial unloading with the secant stiffness K j 1 of the previous step thereby avoiding sudden load drops as in SLA. The load is kept close to the previously applied load, a reloading secant-incremental stiffness K is determined based on a heuristically evaluated local constitutive D matrix, and the minimum load factor is determined to reach the envelope using the defined stiffness relation. The new secant stiffness K j is then determined a posteriori. Additionally, the approach is shown to allow for non-secant unloading, without allowing the damage to increase, by discretising the unloading path using intermediate multi-linear branches. Unlike the improved total analysis, this approach is incremental in principle and therefore overcomes the non-proportional loading problem naturally. The approach however is computationally intensive owing to the need of evaluating the global stiffness matrix twice per analysis and there is also no proof that there always exists an admissible solution for .

Combined incremental-total approaches
Two combined incremental-total approaches were first presented by Graca-e-Costa et al. and were called the Automatic Method and the NIEM (Non-Iterative Energy based method) [21]. Both approaches, in contrast to SLA, use multi-linear material laws to ensure that all non-linearities are due to loading or unloading. The idea of both these methods it to employ a non-iterative incremental solution and then to shift to the total approach intermediately when bifurcation points arise, to guarantee a unique admissible path which may not be traced otherwise. The incremental part is controlled by increment of loads leading to largest energy release at the structural level in the FE model, in line with the Energy release control proposed by Gutierrez et al. [23]. Locally, at the integration point level, for a certain load increment L, the evolution of the constitutive law is forced to follow the valid path, either positive (loading) or negative (unloading), that ensures the release of the largest amount of energy. To begin with, the critical load factor crit is evaluated in a trial step, which corresponds to reaching the nearest point in the multi-linear law for an integration point. Subsequently, the true step ensuring = L L true crit trial is performed and the material stiffness is updated. This carries on as long as bifurcation points do not arise. When they do, the total approach takes over. The stiffness update in each total step could be done in two ways and accordingly two approaches were deduced.
In the first approach, referred to as the Automatic Method, similar to SLA, the secant stiffness is reduced by a predefined factor as in the saw-tooth law. Subsequently, only one of the points becomes critical and reaches the original envelope while all others are on their respective current secant paths. And after intermediately traversing through this region, the method shifts back to the incremental approach using tangent stiffness matrix D. In the second approach, referred to as the NIEM (Non-Iterative Energy based method), the a priori definition of the secant stiffness reductions is avoided. When an integration point encounters a bifurcation, a true incremental step is performed, although it is incorrect. This is done to determine the new updated stiffness for all non critical points that may be on the envelope or unloading, and compute their secant stiffnesses as: while keeping the secant stiffness fixed for the critical bifurcation points. This step in any case is deemed null and the total approach is adopted with the computed secant stiffnesses. Eventually, the approach reverts back to the incremental part in a similar way to the Automatic Method. Although the material loading history is taken care of in the incremental part, the total part of both approaches still required a non-proportional loading approach which in turn was based on DeJong's approach of principle of superposition of stresses followed by determination of load factors. However, the authors mention that the adopted procedure overcomes the development of new nonlinearities during the initial loading (non-proportional loading) due to the preferential usage of the incremental solution. Both methods have been validated against several benchmarks concluding that the NIEM leads to more consistent results in terms of energy dissipated, peak capacity and the consistency condition. Both methods were able to properly track the material loading history as against a total method like SLA due to the inherent incremental nature. Additionally, no material parameters were to be regularised, as in SLA to ensure mesh objectivity, due to the use of multi-linear constitutive relations.

Incremental Sequentially linear analysis
The Incremental Sequential Linear Analysis (ISLA) was proposed recently [24], in which an incremental procedure is used in combination with internal iterations intermittently. In contrast to the traditional incremental-iterative methods, this approach uses the secant stiffness relation globally and therefore the solution path will have jumps similar to the combined incremental-total approaches. The imposed load is either force or displacement controlled. In general, for a non-proportional loading case with two loads L 1 and L 2 , for each load step the loading and stress history in stored because of the incremental nature. Three variations namely: load control, damage control and load & damage control were proposed depending on the method of stiffness reductions and how the corresponding loading factors , 1 2 are treated. In the load control, the 1 corresponding to L 1 is fixed for all load steps (once fully applied) and the 2 corresponding to L 2 is fixed in a certain load step j 1. A critical element is identified based on the ratio of the strength to the total governing stress accumulated thus far. For the subsequent step j, the second load factor 2 is increased and if the stress in the critical point exceeds the allowable strength, damage is increased based on the saw-teeth definition and the stiffness matrix is updated. The displacements and other system variables of the previous step j 1 are restored and the internal forces are calculated to check for equilibrium. If not, damage is again induced to the same structure-state as step j 1, for the same or another critical element, until the algorithm reaches a temporary equilibrium states using Newton-Raphson iterations. In the damage control method, once a critical element is damaged, the damage is kept constant and the load multiplier 2 corresponding to L 2 is repeatedly reduced based on a constant reduction value or a function, based on N-R iterations to arrive at an equilibrium state. In the load & damage control method, once a critical element is damaged, the damage is also varied (in the same or different critical elements) along with the load multiplier 2 corresponding to L 2 , which is also iteratively reduced, to arrive at an equilibrium state.
The load control method in ISLA can be seen as an equivalent of the Force-Release method which also traverses through disequilibrium states to arrive at equilibrium. Additionally, they are similar in the sense that post peak softening behaviour can be realised only if the imposed loads are prescribed displacements.

Discussion on the methods
The sequentially linear approaches thus differ on two main accounts. Firstly, the application of load and the tracking of its history is different under a non-proportional loading situation. The class of total approaches unload and reload non-proportionally. However, as long as the current damage situation allows for all previously applied loads (or the so called constant loads) to be applied without new failures, the constant loads are truly kept constant and this can be interpreted as unloading only until the beginning of the variable load. When the limit point situation arises, previously described in Section 3.1, the constant load could potentially cause premature failure problems which brings about the need for an explicit redistribution. Of all the aforementioned total approaches, only the constrained maximisation analogy [5] provides an alternative in such a situation by allowing to unload and reload from the beginning of the constant load. Or, in other words, an intermittent proportional loading is performed to allow for redistribution by scaling the last successful load combination proportionally and potentially reach a damage state which allows for full recovery of the constant loads. In case of the class of incremental methods, the inherent incremental nature keeps track of the loading history thereby avoiding the aforementioned problems with regards to non-proportional loading.
Secondly, the avalanche of ruptures or multiple failure events are addressed differently by both classes of approaches. In case of the total approaches like SLA, these manifest as sharp snap backs or sudden decay of loads. SLA brings about an instantaneous stiffness reduction in the damaged element thereby relieving it of the stress. However, the neighbouring elements that are almost loaded up to their respective strengths are unable to take over the released stress. In such a scenario, the SLA set-up allows for a reduction in load factor as the only option to redistribute which brings about a snap back in the response. This happens multiple times inducing failures in the surrounding elements before allowing for a full recovery of the load factor close to the true equilibrium path. While in the case of the incremental approaches like the Force-Release method, the load is kept constant and the avalanche of ruptures are achieved in states of disequilibrium using a sequentially linear redistribution loop. Herein, the unbalance forces that may arise due to a damage event are gradually reduced. Thus, both the approaches address the time scales of redistribution differently as has been previously pointed out by Liu and Sayed [26]. It is understood that when the internal forces of the damaged elements are released instantaneously, the total methods are preferable. Contrarily, when the internal forces are released gradually, the elemental failure sequence becomes crucial, and the Force-Release method becomes suitable. Liu and Sayed [26] conclude from their illustrations with lattice models that damage evolution in the Force-Release method ensures compatibility of controlled displacement loadings throughout the analysis and that the area under the force displacement curve is the exact amount of work done by the external loads and the released internal forces (Table 1).
With respect to constitutive relations, incremental formulations allow for discretised laws with elastic stiffnesses for evolving stress points. This means unloading under arbitrary stiffness is possible, which allows for plastic strains unlike the total approaches which rely purely on damage mechanics. In other words, total approaches unload to zero-stress-strain state. Furthermore, such modified constitutive models would make it possible for incremental approaches to be applied to cyclic loading problems.
In the current study, it is sought to understand the differences between both classes of approaches at a structural level using simulations involving continuum elements. The validity of unloading and reloading from the beginning of the constant loads, upon reaching the limit point, as in the constrained maximisation approach in the SLA is looked into further in this study. This mechanism of redistribution is compared to the gradual static stress redistribution brought about in a Force-Release simulation and conclusions Table 1 Constitutive modelling and applications of the discussed approaches.

Method Constitutive models & illustrations
Dejong's method [7] Continuum models -Smeared fixed crack model and discrete cracking interfaces using uniaxial sawteeth relations Force-Release (F-R) method [10] Lattice models -elastic brittle relations, Continuum models -Smeared rotating crack model [25] General method [12] Lattice models with meso-level internal definition -Mohr Coulomb surface with tension cut-off Constrained maximisation analogy [5] In addition to the illustrations by DeJong [7], bond-slip relations for reinforcement, coulomb friction model for interfaces Automatic method [21] Similar to the illustrations by DeJong [7] except using uniaxial multi-linear relations NIEM method [21] Similar to the illustrations by DeJong [7] except using uniaxial multi-linear relations Improved total analysis [9] Mixed mode traction based model -discrete cracking Secant-Incremental analysis [9] Mixed mode traction based model -discrete cracking Incremental sequentially linear analysis (ISLA) [9] Continuum & discrete models -with softening, plasticity etc. on suitability of either methods are drawn. For simplicity, the constrained maximisation variant of SLA is simply referred to as SLA hereon.

Comparative study: Sequentially Linear Analysis vs Force-Release method for non-proportional loading
The mechanism of stress redistribution in a total approach and its incremental counterpart is compared in this section using structural continuum finite element models. All FE simulations in this study are performed in the commercial FEA package DIANA and under 2D under plane stress assumptions. The constitutive model is that of the smeared fixed crack/crush model using a simple Rankine-type failure surface for crack/crush initiation. The shear retention is damage based [13] and the uniaxial saw-tooth relations are of the ripple band type illustrated as in Fig. 1(a) with p the saw-teeth discretisation factor, unless stated otherwise. In the presence of interfaces, similar ripple band saw-teeth relations are used to simulate discrete cracking and for reinforcements, an elasto-plastic saw-teeth relation is used to simulate yielding. The latter is kept damage based as in Reference [5] for suitability to SLA. Alternatively, for the F-R, the discretised elasto-plastic constitutive law can be made with elastic stiffnesses to attain true plastic unloading but this is not included herein. Pure sliding failures, which are possible in typical masonry mechanism, may not be relevant for masonry walls which fail in rocking like in Section 4.2.1. But the facade study in Section 4.3 which involves the diagonal failure mechanism in piers requires it in principle. This is possible using a meso-model in combination with coulomb friction interface formulation for SLA [5] but, considering the scale of the facade, the study is limited to the orthogonal fixed cracking approach under macro modelling as previously done in Reference [15].
The loading condition for all the chosen examples is non-proportional and is composed of two stages: firstly a precompression or prestress is added as force loads, followed by a second stage of imposed/prescribed displacements. The force loads and prescribed displacements are depicted in orange and green/yellow colours respectively in the FE model descriptions for consistency.

Experiment
Several three-point bending tests were performed on plain concrete beams in an experimental campaign by Zhao et al. [27]. These were previously analysed using SLA with plane stress assumptions [5]. Of these beams, SG2-B1, of span 600 mm, 150 mm depth, and 120 mm thickness is chosen as reference for this study. The beam is tested in three point bending with a notch of depth 60 mm and 2 mm width and a schematic representation of the experiment is shown in Fig. 2(a). The middle of the beam is in pure bending and if such a beam were to be subject to axial loads (purely for an academic purpose) like a prestress force, the loading would become nonproportional. The effect of compressive axial loads will increase the capacity of a beam and has been a proven benchmark study for testing non-proportional loading strategies in SLA thus far [7]. This simple case study is used to compare the performance of SLA and the Force-Release method.

Finite element model
The simply supported beam is modelled with geometry as shown in Fig. 2(a) except that the notch width is assumed to spread over an entire element length. The FE model is made using linear 4-noded iso-parametric plane stress elements, with 2 translational degrees of freedom per node and a regular 2x2 Gaussian integration scheme. All elements are square and approximately 10 mm in size. Concrete has been modelled as a linear-elastic material everywhere except the notched column of gray coloured elements, see Fig. 2(b), where all the physical non-linearity is lumped (only tensile). A linear tensile softening saw-tooth law with material parameters mentioned in Table 2 is used. Two cases are considered in this section. First, the three point bending test of SG2-B1 beam as in the experiment is simulated. Since the notch extends over one whole linear element, the point load is assumed to be equivalent to two prescribed displacements applied either side of the central element at the top of the notch as shown in Fig. 2(b). Second, the beam is firstly assumed to be subject to a uniform compressive prestress, of magnitude 0.42 MPa, at the ends of the beam before proceeding to the three point bending.

Results & discussion.
The force displacement curves from SLA and the Force-Release method, for the first case without axial loads, match up well with those from the experiment. Furthermore, both simulations show the propagation of the tensile crack through the height of the notch resulting in the reduction of effective cross section of the beam. This behaviour results in loss of flexural capacity and eventually is observed as the global softening in the force displacement curves shown in Fig. 2(c). And similar to the case without axial loads, the SLA and Force-Release response are quite similar for the non-proportional case of axial loads followed by the bending load. The presence of axial loads like the compressive prestress tends to delay the cracking and consequently results in a higher capacity. However, there are small variations between SLA and Force-Release in both sets of responses with and without the prestress force. The SLA shows small snap backs in the response which is indicative of the load-unload nature of the approach and the way multiple failure events are handled. Since one event occurs at a time, a potential event in a neighboring integration point occurs only by allowing for a reduction in the load in the next linear analysis. Conversely, in case of the Force-Release method, vertical drops are observed in the force displacement response rather than snap-backs. Since the previously applied load cannot be modified, for a constant prescribed displacement, an additional event or more occurs due to the gradually released unbalanced forces due to the previous event. This is the previously described avalanche of ruptures denoting a region of disequilibrium. Additionally, since the constant load i.e. the compressive axial load isn't too high, SLA doesn't run into the intermittent proportional loading and the associated redistribution. That is, SLA could keep the prestress load (constant) throughout the entire simulation, while fully satisfying the constitutive laws. This case study shows the inherent difference in load application between the two approaches. The small differences in element failure sequence is not very apparent due to the problem definition, wherein the overall final crack is confined to a single band of the elements, and is therefore not deliberated upon.

Slender masonry wall 4.2.1.1. Experiment.
TUDCOMP-20 was one of the several quasi-static in-plane cyclic tests performed recently at the TU Delft laboratories. This is part of an extensive experimental campaign to characterise the response of unreinforced Dutch masonry buildings to induced seismicity [28]. The specimen was a calcium-silicate brick masonry wall, roughly 2.7 m high and 1.1 m long, with a thickness of 102 mm. The wall was tested under cantilever boundary conditions, i.e. it was clamped at the bottom and free to rotate at the top. Firstly, an overburden stress of 0.63 MPa is applied on the wall through a combination of actuators and a steel beam as shown in Fig. 3(a). This load is maintained constant throughout the experiment while allowing for the rotation of the top of the wall. Subsequently, a cyclic lateral displacement is applied to the top of the wall using a horizontal actuator connected to the top steel beam. This wall exhibits pure rocking movement and for large displacements, splitting cracks were observed at the bottom-left and bottom-right corner of the wall, see Fig. 5(a). This was followed by the phenomenon of toe crushing when the base shear force dropped substantially for the same imposed displacement, starting at the red dot in the experimental response, see Fig. 3(b). Considering the scope of this article, the explosive failure in this non-proportionally loaded experiment makes for an interesting case to be investigated using the SLA and Force-Release methods. The explosive failure led to instability of the wall and its eventual collapse (See Fig. 5). For further information about the experiment, the reader is referred to Reference [29].  Fig. 4(b). The finite element model of the wall, as shown in Fig. 4(a), is made with 4-noded iso-parametric plane stress elements, roughly 55x55 mm in size, with linear interpolation shape functions and a 2x2 Gaussian integration scheme. The compressive non-linearities are in these continuum elements and are described by a smeared fixed crush framework, while behaving linearly in tension. The uniaxial saw-tooth compressive softening law is of a parabolic type [17] with material parameters as described in Table 3. Additionally, a zero-thickness interface (shown in green in Fig. 4(a)) is provided at the bottom of the wall to simulate the discrete cracking leading to the rocking behaviour. This is made using 2+2 noded interface elements, roughly 55 mm in length each, with linear interpolation shape functions and a 2 point Newton-Cotes integration scheme. The tensile non-linearities in these interface elements are described using a linear tension softening type saw-teeth law with material parameters mentioned in Table 3. The FE model also includes another row of continuum elements at the top of the wall (not shown herein) with stiffer properties behaving linear elastically. This is to ensure the cantilever boundary conditions. The bottom of the wall is fixed as shown in Fig. 4(a). The overburden load is applied first as force loads at the top edge of the wall and subsequently, a prescribed displacement is applied at the left top corner. The assumption in the FE model to separate tensile and compressive non-linearities into the discrete cracking interface and the continuum plane stress elements is done to simplify the problem while accurately capturing the failure mode.

Results & discussion.
Firstly, both the SLA and Force-Release simulations exhibit the rocking behaviour, which is localised in the discrete cracking interface, followed by compressive softening in the toe region of the wall. This leads to dissipation of energy  resulting in a rather ductile response as shown in Fig. 4(b). Both SLA and the Force-Release simulations are rather identical but showcase minor differences as previously observed in the prestressed beam case. However, the responses show a good qualitative agreement to the experimental backbone curve until the collapse mechanism begins. The explosive failure is captured adequately by both methods, but the mechanism of redistribution differs. The point of onset of the mechanism is denoted by a green mark in the Fig. 4(b) and it occurs around the same imposed displacement for both simulations. The damage patterns are also identical as seen in Fig. 6(a) where DmSS is the amount of damage accumulated in the discrete cracking interface and the crushing continuum. DmSS ranges from 0 to 1, which correspond to undamaged and fully damaged conditions. It is clear that the two elements at the bottom right corner of the wall are fully crushed and that all interface elements to the left of this region are completely cracked, leaving a very tiny portion which effectively supports the wall. The ensuing mechanism is described by both approaches differently. On the one hand, SLA runs into the limit point situation described previously in Section 3.1.1, where there is no constitutively admissible critical load multiplier. The intermittent proportional loading (IPL) commences and the last successful load combination is scaled proportionally. Firstly, the IPL occurs a little before the onset of collapse as well but recovers back to the conventional non-proportional loading. However, once the collapse begins the IPL never recovers which is evident from the amount of precompression, the first load applied (constant), that remains on the structure in the rest of the simulation, see Fig. 4(c). The IPL implicitly reduces the constant load, thereby describing the entire dynamic brittle collapse mechanism while maintaining equilibrium. On the other hand, the Force-Release method runs into an avalanche of ruptures while going through disequilibrium states. Since the previously applied load cannot be altered, for the same imposed displacement and the full value of precompression, the Force-Release method attempts to allow for redistribution due to successive failure events by gradually releasing the stresses. The ongoing failure is therefore captured differently by both approaches as seen in Fig. 6(b).

Fig. 5.
Failure pattern of masonry wall in the TUDCOMP20 experiment before collapse (left) and after collapse (right) [29].

Table 3
Modelling parameters -TUDCOMP20: identifed by sensitivity studies based on material properties obtained from the experiments [28,29]. The difference between the approaches in describing collapse may be interpreted as two extremes of the time scales for a dynamic redistribution process [12]. SLA in this situation is essentially assuming that the loading equipment is fast enough to react to the collapse mechanism, alter the load and consequently release the stresses quickly to avoid further failures. This is clear as the eventual failure pattern additionally only involves the crushing of the tiny effective portion supporting the wall before the onset of collapse.
The crush zone appears to be more realistic, in a quasi-static sense, wherein SLA gives room for damage propagation while quasistatically releasing the loads. The Force-Release method on the other hand stays true to the displacement controlled experiment, and realises the full collapse by gradually releasing the stresses in a sequentially linear redistribution loop. Although the process is dynamic, the Force-Release method effectively neglects all inertial forces. Since vertical equilibrium is not possible anymore the simulation could be interpreted to have been completed and the wider crush-zone is indicative of the instability. In summary, both approaches adequately describe a real non-proportionally loaded experiment involving true brittle collapse, in terms of the failure patterns and the eventual mechanism, differing only in their respective approaches to the latter.

Squat reinforced concrete wall 4.2.2.1. Experiment.
The third benchmark is a reinforced concrete shear wall, named SW13, from the experiments of Lefas et al. [30].   5 kN). At 40% of the horizontal load (130 kN) the first inclined crack appeared and simultaneously the flexural cracks had already spread at a slight inclination within the wall web. This was followed by growth in the flexural and inclined cracks, almost reaching the compressive edge. The first yielding of the tensile reinforcement occurred at 75% of the horizontal load (250 kN) which was followed by the failure. Fig. 8(b) shows the failure that occurred by way of near-vertical splitting of the compressive zone at a peak load of 330 kN. For further information on the experiment, the reader is referred to Reference [30]. Fig. 8(a) shows the discretised model of the shear wall used in the SLA and the Force-Release simulations. 6 noded iso-parametric triangular continuum elements of 75 mm, with quadratic interpolation shape functions and a 3point Gaussian integration scheme, were used. Reinforcements are modelled as embedded units and are shown in yellow in Fig. 8(a). The constant vertical load of 355 kN was applied as two uniformly distributed loads over a distance of 125 mm on the top edge of the wall, to replicate the load application to the spreader beam as in the experiment. The horizontal load was applied to the right side of the header beam at mid-height as a prescribed displacement. The foot of the wall was not modelled because its influence was found to be negligible [32]. Instead, the base of the wall was fixed with translational restraints in both directions to model the clamped wall base. Table 4 shows the material parameters used in the simulations for concrete, according to the fib Model Code 2010 [31], and reinforcement properties obtained from the experiment.   Fig. 9, exhibit the decisive failure as a drop of load. Although the response from both the approaches is stiffer than the experiment, have a delayed onset of cracking and exhibit a lower peak load, the qualitative behaviour is quite consistent with the experiment. The aforementioned problems isn't a feature of the SLA or the Force-Release methods but also arose in previous NLFEA simulations and is more a feature of the constitutive modelling used [33,32].The material parameters in Table 4 are based on previously done NLFEA studies [33,32] to obtain the best fit. The flexural crack progression is followed by the inclined cracks until the peak load (Fig. 10) and eventually the toe crushing combined with yielding of reinforcements is captured well (Figs. 11 and 12). The toe region of this wall undergoes crushing of the concrete (see Fig. 13), compressive yielding of the reinforcements (see Fig. 14) and additionally exhibits splitting cracks at the onset of collapse. The splitting cracks in this region are due to the lateral compressive softening, in addition to the presence of the yielding vertical reinforcements that are close to rupture. The collapse mechanism begins with the crushing of concrete at the toe and the yielding of the reinforcement as seen in Figs. 13& 14 which further propagates through the simulation. Both approaches capture this adequately. However, the drop of load corresponding to the eventual instability is described by the SLA and the Force-Release methods in diametrically opposite ways with regards to the time scales for the redistribution. This is in line with the differences observed between the approaches to the explosive failure in the previous case study and is clear from how the loading is modified in case of SLA ( Fig. 9(b)) during collapse. SLA again describes the collapse mechanism while maintaining equilibrium by reducing the constant load while the Force-Release method addresses it using the avalanche of ruptures in disequilibrium. As the complex stress state in this zone evolves and undergoes redistribution, SLA experiences stress reversal problems. The evolving splitting cracks (due to lateral compression softening) in this complex zone causes unloading in surrounding elements resulting in crack closure in some elements, with incorrect damaged stiffness. This is a characteristic to the approach, wherein damaged stiffness is carried onto tension or compression upon stress  reversal due to redistribution. This therefore hinders the splitting crack formation at the peak and close to the onset of collapse in the case of SLA, see Figs. 10 and 11. On the other hand, at a similar stage in the simulation, Force-Release captures the splitting cracks better as observed in Fig. 11. Both the approaches experience stress reversal issues but the problem is much less accentuated in case of the Force-Release method because of the inherent incremental nature of the approach in keeping track of the stress history.   Nevertheless, a dedicated stress reversal event and a corresponding load multiplier should be incorporated in line with previous proposals for SLA in continuum models [34] and Force-Release method in lattice models [35].

Experiment
A full scale masonry building tested by Magenes et al. [36] is a well known benchmark study to investigate the seismic response of masonry buildings [15]. It consists of four walls: walls A + B + C which are connected by an interlocking brick pattern around the corner and wall D which is disconnected from perpendicular walls A and C ( Fig. 15(a)). This detached wall can be simulated using 2D plane stress elements and thus fits the scope of the current study. Although the experiment was a cyclic test, this benchmark is simulated with monotonically increasing lateral loads and the backbone curve of the cyclic response is used as reference for comparison. The geometry of the facade is shown in Fig. 15(b) and has a thickness of 0.25 m. The facade is loaded in two stages: First, the floors are loaded by a distributed load of 10 kN/m 2 , resulting in total vertical loads of 248.4 kN and 236.8 kN at the first and second floor respectively. The floors rest on beams, which are connected to walls B and D, such that half of the floor load is carried by each facade. Second, the building is loaded by a monotonically increasing lateral load at both floors to simulate the static equivalence of seismic action. The damage pattern corresponding to the cyclic test at an ultimate displacement of 23 mm is as shown in Fig. 15(c).

Finite element model
The FE model of wall D is made of 8 noded iso-parametric quadrilateral plane stress elements with quadratic interpolation shape functions and 2x2 Gaussian integration scheme. All elements in white, representing masonry, are provided tensile and compressive non-linearities with parameters mentioned in Table 5. The ones of the lintel and those along the floor levels in gray are kept linear elastic since cracks are not expected here [15]. All elements are of size 230 mm. Loading of the first stage contains vertical loads applied as line loads along the height of the two stories and additionally gravity loads. The second stage involves applying prescribed displacements to an auxiliary frame, as shown in Fig. 16, such that equal monotonically increasing lateral forces are applied at both floors (using the master-slave connection) to simulate the static equivalence of seismic action. This is in line with a mass-proportional  pushover analysis (assuming the floors have the lumped mass), contrary to a modal pushover analysis. The facade is supported at the bottom in both horizontal and vertical direction to simulate the fixed boundary condition.

Results & discussion
The results of the SLA and Force-Release simulations of the pushover test are shown in Fig. 17. The results show good qualitative agreement against the experimental backbone curve. Since the experiment was cyclic in nature, the damage accumulation over the test i.e. the progressive degradation in stiffness affects the peak load achieved as against a monotonic counterpart. Accordingly, the simulations show a slightly higher peak load. However, the first stage of the rocking of the piers occurs, followed by diagonal splitting type cracks in the center and right piers is in accordance to the experimentally observed damage crack patterns. Around peak load, the first drop of load corresponds to the development of the aforementioned diagonal cracks in the central and right piers. This is shown in Fig. 18. SLA and Force-Release methods differ very slightly in the crack patterns and this corresponds to minor differences in the elemental sequence failure that is characteristic of both approaches. However, around this first drop of load, SLA requires the intermittent proportional loading for redistribution due to the dynamic crack propagation in the central and right piers. These are also seen as drops of the constant load between 2000th and 4000th step (linear analysis) in Fig. 17(b). Although such a redistribution driven by the intermittent proportional loading is premature, i.e. much ahead of the actual structural collapse observed in previous case studies, the redistribution allows for recovery of the full magnitude of the constant loads and is therefore acceptable. Force- Release, on the contrary, goes through three large avalanche of ruptures for corresponding constant imposed displacements through this region. This again illustrates the inherent differences in the approaches (Figs. 19 and 20).
At an ultimate displacement of around 18 mm, SLA runs again into the intermittent proportional loading. However this time the actual collapse begins with the right pier completely failing by diagonal splitting failure and also lateral splitting near the toe region due to the compression softening. Force-Release, owing to the inherent incremental approach and the manner of redistribution, proceeds further to an ultimate displacement of 24 mm leading to the collapse while going through disequilibrium states. The collapse mechanism is also captured similarly by both approaches varying again only in the manner of load modification. SLA enforces a relaxing mechanism in equilibrium while Force-Release describes the explosive failure as an instability.

Conclusions
Finite element simulation of quasi-brittle materials using incremental-iterative methods often suffers from robustness issues. Consequently, several sequentially linear methods were devised to address this problem. This article, firstly, provides an extensive overview on all such alternatives available in recent literature [5,21,9,10,12,24]. Furthermore, the mechanism of redistribution due to dynamic failure propagation in a quasi-static set-up, as addressed by the two classes of sequentially linear methods: namely the  total and incremental approaches, is qualitatively discussed based on previous perspectives to this topic from the field of lattice modelling.
Although several studies have been made in the lattice modelling community to elucidate the inherent differences between these classes of methods [10,22], their suitability to structural applications has been touched upon only sparsely. Of all the available total approaches in literature, only the constrained maximisation analogy [5] addresses the limit point situation wherein no constitutively admissible critical load multiplier can be found. This method is therefore chosen as the reference total approach to be compared against an incremental counterpart, the Force-Release method. The inherent differences between these methods is illustrated using continuum FE models of structural examples involving non-proportional loading in this article. To this end, an academic example of a prestressed beam and three experimental benchmarks: a slender masonry wall, a squat reinforced concrete wall and a full-scale 2D masonry facade, all subject to different magnitudes of precompression followed by pushover are studied. In general, it is observed in all 4 case studies that the non-proportional loading strategy in a total approach like SLA is almost equivalent to the incremental solution obtained using a Force-Release method, in terms of force-displacement curves and damage patterns, until the limit point situation is reached. Prior to this, although the unloading and reloading are both done non-proportionally in SLA, it is acceptable to interpret that unloading happens only until the beginning of the variable load since the all previously applied (constant) loads can be reloaded fully without causing new damage. It is also interesting to note that the differences in the sequences of events, which can significantly affect the eventual failure mechanisms as in the case of lattice modelling examples [22], does not greatly affect the results at continuum level of all 4 presented case studies. However, upon reaching the limit point, which denotes an oncoming dynamic failure propagation region in the true response, the process needs to be driven by a sudden release of energy which inherently requires multiple failure events or the so called avalanche of ruptures. The SLA and Force-Release approaches differ in the redistribution mechanism at this juncture. On the one hand, SLA allows for such multiple failures sequentially, only by reduction of all previously applied loads which brings about snap-backs. Since every damaged element's stress is released instantaneously, the neighbouring elements whose stresses are close to their respective  allowable strengths subsequently become critical in SLA at a considerably lower load. This is done by scaling the last successful load proportionally which in turn, unlike total unloading, allows for damage propagation by reduction of forces thereby giving the structure the time to relax. In other words, SLA lets the damage progress quasi-statically by load reduction during a dynamic event. This is acceptable under proportional loading situations and also under non-proportional loading conditions if the experiments are controlled truly quasi-statically [37]. When the damage process zone is controlled for quasi-static evolution of failure in an experiment [37], SLA is preferable. This can also be simulated using arc-length control procedures in an implicit NLFEA set-up. On the other hand, the Force-Release method keeps track of the load history incrementally and allows for multiple failure events by gradual redistribution of stresses through disequilibrium states. The Force-Release method is suitable for typical displacement controlled experiments which actually exhibit instabilities but is not suitable for physical processes which exhibit snap backs or for truly quasistatic experiments. Thus, both approaches are essentially extremes with regards to time scales for redistribution, which is more apparent only in major dynamic failure processes in the case of continuum models. These dynamic failure processes can either be intermediate in the entire process of damage evolution or the final collapse mechanism in itself. In the former case of an intermediate failure process, with regards to SLA, if the simulation allows for a full recovery of all previously applied loads, the results are deemed acceptable. Such an instance was observed around first drop of peak load in the pushover case study of the 2D masonry facade. However, during these intermediate failure processes, if all previously applied loads have to be gradually reduced to extremely low values to simulate the failure process, then the results of SLA correspond to alternate equilibrium paths of damage propagation that do not culminate in the actual expected collapse mechanism. Or, in other words, premature loss of previously applied load much ahead of the actual mechanism is not equivalent to simulating the experiment truly. In the latter case where the dynamic failure is the actual collapse mechanism as observed in the pushover case studies of the masonry and reinforced concrete walls, SLA tries to enforce equilibrium by load reduction and allows for a relaxed failure mechanism. While the Force-Release method describes both the intermediate and the explosive ultimate failures through avalanches of events that are states of disequilibrium. The Force release method stays true to the instability in the experiment but the excessive damage propagation at collapse as observed in the wider-crush zone of the slender masonry wall case study (Section 4.2.1) needs further deliberation including a localisation analysis. In summary, complex dynamic failure processes which cause robustness issues in incremental iterative methods also affect sequentially linear methods under nonproportional loading and are addressed differently by the two sub-classes: the total and incremental approaches, with their suitability depending on the experiment being simulated.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.