Significance of flow rule for the passive earth pressure problem

Determination of earth pressures is one of the fundamental tasks in geotechnical engineering. Although many different methods have been utilized to present passive earth pressure coefficients, the influence of non-associated plasticity on the passive earth pressure problem has not been discussed intensively. In this study, finite-element limit analysis and displacement finite-element analysis are applied for frictional materials. Results are compared with selected data from literature in terms of passive earth pressure coefficients, shape of failure mechanism and robustness of the numerical simulation. The results of this study show that passive earth pressure coefficients determined with an associated flow rule are comparable to the Sokolovski solution. However, comparison with a non-associated flow rule reveals that passive earth pressure coefficients are significantly over predicted when following an associated flow rule. Moreover, this study reveals that computational costs for determination of passive earth pressure are considerably larger following a non-associated flow rule. Additionally, the study shows that numerical instabilities arise and failure surfaces become non-unique. It is shown that this problem may be overcome by applying the approach suggested by Davis (Soil Mech 341–354, 1968).


Introduction
For the design of retaining walls (e.g. sheet pile, bored pile or diaphragm walls), gravity dams, cantilever walls or bridge abutments determination of earth pressures is a fundamental task. In contrast to the active earth pressure acting as driving force, the passive earth pressure contributes in general to the resisting forces in the system. As over prediction of the passive earth pressure may lead to an unsafe geotechnical design, estimation of the passive earth pressure is an important task for practical engineers.
Furthermore, since earth pressure coefficients are based on different theories and assumptions, it is highly recommended for practitioners to be familiar with the fundamental theories used to derive earth pressure coefficients.
In general, the passive earth pressure force E ph is calculated based on three terms which are related to soil weight c, effective cohesion c 0 and surcharge q. However, for a purely frictional material and no surcharge applied to the backfill, determination of E ph can be reduced to Eq. 1, with the height of the wall h and the passive earth pressure coefficient K ph .
The first theory to determine earth pressures was developed by Coulomb [8] using a translational wedge element to describe the active and passive zones. From the results earth pressure coefficients have been derived for variable soil and wall friction angles (d) as well as different inclinations of the soil-wall contact surface (a) and slopes of the backfill (b). Further milestones in estimation of earth pressures were achieved by Rankine [18], Caquot and Kérisel [3] and Sokolovski [22]. Since then, many contributions for a more realistic determination of earth pressure coefficients have been published for static and/or seismic conditions [4,5,9,[23][24][25]. Recently, Krabbenhoft [11] applied finite-element limit analysis (FELA) to derive earth pressure coefficients as mean value of rigorous upper and lower bound solutions. Krabbenhoft [11] reported that a maximum difference of less than 0.5% between upper and lower bound solutions was obtained, indicating an almost exact solution.
Although the determination of earth pressure coefficients using advanced methods and complex algorithms has become very accurate, in many studies the influence of the flow rule is not discussed sufficiently. In other words, the effect of a non-associated flow rule with a dilation angle w 0 smaller than the effective friction angle u 0 on the limit load is not taken into account, despite the fact that soils do not obey an associated flow rule. Depending on the degree of non-associativity (K ¼ u 0 À w 0 ) this may result in significantly smaller K ph values. Thus, assuming associated plasticity could yield results on the unsafe side. Particularly for large K values this effect is accompanied by numerical instabilities as well as non-unique failure mechanisms [12,14].
Due to the fact that in some circumstances the obtained failure load (or factor of safety) is significantly over estimated following an associated flow rule, Davis [6] introduced a reduction factor b to determine reduced shear parameters c Ã and u Ã which should account for the nonassociativity of the material. Applying the reduced shear parameters, as given in Eqs. 2-4, in combination with an associated flow rule should serve as an approximation of the response (e.g. factor of safety in slope stability problem) of a non-associated material.
To evaluate the influence of the flow rule for slope stability problems Tschuchnigg et al. [27] compared factors of safety (FoS) for different slopes using a frictional and a cohesive-frictional material following both an associated and a non-associated flow rule. In addition, calculations with an associated flow rule using reduced shear strength parameters according to Davis [6] approach were performed to analyze whether good agreement with the FoS using the non-associated material can be achieved. The comparison revealed that the Davis approach provides a conservative estimate to the FoS of an non-associated material, whereas an associated flow rule yields an over estimation of the FoS. Analyzing the differences between these three solutions, it was found that in many cases the gap between associated and non-associated flow rule was significantly larger than the gap between non-associated flow rule and Davis approach, putting emphasis on the use of the Davis approach. However, as Tschuchnigg et al. [27] only investigated the influence of the flow rule for the slope stability problem the question arises to what extent the flow rule influences, e.g. the passive earth pressure problem. Most importantly, K ph values have to be evaluated for different combinations of u 0 , w 0 and d to conclude whether (1) following an associated flow rule leads to significant over estimation and/or (2) applying the Davis approach yields approximations, which are too conservative for the use in design.
Therefore, the objective of this study is to analyze the influence of the flow rule on the passive earth pressure coefficient K ph for different combinations of u 0 , w 0 and d, the failure mechanism and the robustness of the numerical calculation for frictional materials. For this purpose, passive earth pressure coefficients derived from displacement finite-element analysis (FEA) using Plaxis 2D [2] will be compared to finite-element limit analysis (FELA) performed with OptumG2 [16] and selected data from literature.

Displacement finite-element analysis (FEA)
Within the last 20-30 years displacement finite-element analysis (FEA) has become very popular in geotechnical engineering. Due to its flexibility of modeling any geometrical problem, implementation of interface elements to approximate soil-structure interaction and advanced algorithms that enable fast calculations, FEA can be applied to a wide range of geotechnical problems. Combined with advanced constitutive models which enable to take into account hardening/softening, anisotropy, creep, non-associativity, barotropy, pyknotropy and other effects, the soil behavior can be approximated accurately. In this study displacement finite-element analysis is applied using the code Plaxis 2D [2].

Finite-element limit analysis (FELA)
Finite-element limit analysis (FELA) is based on the upper bound and lower bound theorems of plasticity introduced by Drucker et al. [7]. Finite-element lower-bound analysis was first applied by Lysmer [13] for plane-strain conditions. Bottero et al. [1] were among the first to apply finiteelement upper-bound analysis for geotechnical problems considering plane-strain conditions. Since then many contributions helped to continuously improve both methods.
For more details of the historical development as well as the theoretical background of both finite-element lowerand upper-bound analysis the reader is referred to Sloan [21]. In general, in FELA it is assumed that the material is rigid-perfectly plastic and obeys an associated flow rule. Furthermore, it only applies for small deformations. Using FELA in combination with the correct element types [21], rigorous upper and lower bound solutions can be obtained where the true collapse load is bracketed from above and below. In this study, finite-element limit analysis is applied using the code OptumG2 [16].

Sokolovski method
Lateral earth pressure coefficients according Sokolovski [22] are determined following the Rankine [18] method using the method of characteristics. Since a detailed explanation of this method can be found in [22], only a brief description is given in this section. Based on the static method for determining earth pressure, two differential equations [10] are solved for every infinitesimal element within the plastic zone. The soil mass within the plastic zone is assumed to be at failure and, thus, follows the Mohr-Coulomb failure criterion. Furthermore, plane-strain conditions are assumed. For a cohesionless soil the shear strength is governed by the effective friction angle which is fully mobilized within the plastic zone. Moreover, it should be mentioned that in contrast to limit equilibrium approaches, the method of characteristics does not require an assumption of the shape of the failure mechanism. However, Davis [6] noted that for solutions obtained using the method of characteristics the kinematic admissibility needs to be investigated additionally in order to examine whether all deformation boundary conditions of the specific problem are satisfied.
3 Passive earth pressure problem

Problem description and spatial discretizations adopted
To analyze the passive earth pressure problem, a numerical model is created, which is has the same dimensions and boundary conditions for FEA and FELA. Geometry and dimension are shown in Fig. 1a. The mechanical boundary conditions of the model are as follows: at the left and right boundaries deformations are restricted in normal direction, whereas at the bottom deformations are restricted in both normal and tangential direction. Deformations of the rigid block, which is attached at the top left part of the model, are restricted in vertical direction. No mechanical fixity is applied to the soil top surface of the model. The soil is assumed to be homogeneous and effects due to ground water are not considered. To simulate the passive earth pressure problem, the rigid block is pushed into the soil by applying either a force or displacement at the left vertical edge of the block. As vertical displacements of the rigid block are prohibited the block is translated in horizontal direction only. In displacement finite-element analysis (FEA) using Plaxis 2D the domain has been discretized using 15-node triangular elements with a shape function of 4th order. The spatial discretization adopted in FEA is presented in Fig. 1b. It can be seen from this figure that the overall element size in Plaxis 2D is small with further refinement of the elements in the area of the expected failure mechanism (2 h 9 4 h) and next to the interface elements between the rigid block and soil. This discretization has been found to be appropriate after performing a mesh study with different mesh element sizes.
To obtain rigorous lower and upper bound results with finite-element limit analysis 3-and 6-node triangular elements have been chosen, respectively. Due to the adaptive mesh refinement, the critical failure mechanism can be obtained with less elements and with less computational costs as has been shown, for instance by Oberhollenzer et al. [15] and Schmüdderich et al. [20] for the slope stability and bearing capacity problem, respectively. In FELA, the mesh has been refined in five adaptive steps using the shear dissipation as decisive parameter leading to a rather heterogeneous discretization with coarse elements towards the bottom and the right edge of the model and very fine elements close to the failure mechanism (Fig. 1c). In general, adaptive mesh refinement can also be implemented for displacement finite-element analysis in order to enable determination of precise results in terms of the approximation of the failure mechanism and the corresponding earth pressure force.

Material parameters
To analyze the earth pressure problem a linear elasticperfectly plastic material model using a Mohr-Coulomb (MC) failure criterion is applied. Apart from the soil weight c, the basic stiffness input parameters of the model are Young's modulus E and Poisson's ratio m. The strength of the soil is defined with the effective cohesion c 0 , the effective friction angle u 0 and the dilation angle w 0 . In this study, soil weight, Young's modulus, Poisson's ratio and effective cohesion are assumed to be constant (c ¼ 16 kN=m 3 , E ¼ 30 MPa, m ¼ 0:3 and c 0 ¼ 0 kN=m 2 ), while soil friction angle, dilation angle and soil-wall contact, expressed in terms of the interface friction angle d and modeled with an individual set of parameters, are varied. The importance of the flow rule (associated, non-associated) is investigated for the passive earth pressure problem for different soil friction angles u 0 and wall friction angles d. As reported among others by [26] for the slope stability problem, it should be noted that the degree of non-associativity K, defined as the difference between friction angle and dilation angle (K ¼ u 0 À w 0 ) is an important parameter which needs to be considered in the later evaluation. As simulations with the Davis approach are conducted, reduced shear strength parameters u Ã and c Ã are derived based on Eqs. 2-4.

Results and discussion
In this section, the results of the study are presented and general trends for different combinations of the friction angle u 0 , the dilation angle w 0 and the wall friction ratio d=u 0 are pointed out. Evaluation and discussion of the results are conducted in terms of passive earth pressure coefficients, failure mechanisms and performance of the numerical simulation within three subsequent steps. These steps account for (1) analyzing results obtained with an associated flow rule and (2) with a non-associated flow rule as well as (3) comparing the results presented in (1) and (2) with the estimation of K ph values obtained following the Davis [6] approach. For an easier distinction of the different calculations performed, it should be noted that letters A, NA and D refer to calculations based on an associated flow rule (A), non-associated flow rule (NA) and Davis approach (D), respectively.

Simulations with an associated flow rule
Results obtained using FEA, method of stress characteristics (following Sokolovski [22]) as well as lower (LB) and upper bound (UB) FELA are presented in Table 1 in terms of the passive earth pressure coefficient K ph for friction angles between u 0 ¼ 20 and u 0 ¼ 45 assuming an associated flow rule (w 0 ¼ u 0 ) and a wall friction ratio of d=u 0 ¼ 0:67. From results presented in this table it can be seen that, in general, there is a good agreement between all approaches. However, a more detailed analysis reveals that FELA provides smaller K ph values than FEA and even the UB solution is slightly smaller than FEA. An explanation for this observation can be given with respect to the adaptive mesh refinement which allows for more accurate approximation of the failure mechanism in FELA. Furthermore, a mesh fan with a fine subdivision of elements is utilized in FELA, which enables easier rotation of the failure surface at the bottom of the soil-wall contact.
Comparing the method of stress characteristics with the LB solution obtained with FELA, it is apparent that the Sokolovski [22] solution yields smaller K ph values at small effective friction angles (u 0 25 ) but larger K ph values at increasing effective friction angles, except for u 0 ¼ 45 , for which significantly smaller K ph values have been obtained with the first method. Although there is no clear trend to be discussed based on these observations, a possible explanation can be that the design equations [17] utilized to obtain K ph values in accordance with the Table 1 Comparison of passive earth pressure coefficient K ph obtained from FEA, method of stress characteristics and lower (LB) and upper bound (UB) FELA for friction angles between u 0 ¼ 20 and u 0 ¼ 45 , w 0 =u 0 ¼ 1:0 and d=u 0 ¼ 0:67 Equal to results by Krabbenhoft [11] Sokolovski [22] solution are less accurate for very large friction angles. However, the reader should keep in mind that the method of stress characteristics and the FELA yield K ph which compare well with the FEA solution assuming an associated flow rule. In addition to passive earth pressure coefficients, failure mechanisms obtained from both numerical methods are compared. As shown in Fig. 2, failure mechanisms obtained from (a) FEA (Plaxis 2D) and (b) FELA (Op-tumG2) are highly comparable in terms of the length (L x;A % 2:32 h) of the failure surface as well as the angle (m 1;A ¼ 45 À u 0 =2 ¼ 45 À w 0 =2 ¼ 30 ) at which the failure mechanism intersects the ground surface. Based on the performed studies it can also be concluded that an associated flow rule yields a unique failure mechanism which is (almost) independent of the numerical method applied. However, it should be noted that this statement is only valid if the discretization utilized is appropriate. Comparing the size of the failure mechanism with the dimensions of the numerical model, it is obvious that the model dimensions have been chosen large enough in order to avoid boundary effects. This observation also holds for larger friction angles.

Simulations with a non-associated flow rule
In order to investigate the significance of the flow rule for the passive earth pressure problem, FEA simulations have been performed for effective friction angles between u 0 ¼ 20 and u 0 ¼ 45 and ratios of non-associativity in the range of w 0 =u 0 ¼ 0:00, 0.25, 0.50 and 0.75 assuming a constant wall friction ratio of d=u 0 ¼ 0:67. For all simulations, the ratio of non-associativity, which is utilized in the soil domain, is also applied to the soil-wall contact (interface). The results are presented in terms of load-displacement curves expressed as u x =h (x-axis) versus K ph (yaxis) for u 0 ¼ 20 À 45 in Fig. 3a-f, respectively, where u x is the horizontal displacement of the rigid block. However, it should be noted that the solid line in Fig. 3f has been obtained for w 0 =u 0 ¼ 0:05 as FEA simulations with w 0 =u 0 ¼ 0:00 led to highly erratic results and terminated before sufficient displacements could be applied. Since only a small reduction in K ph is expected by decreasing from w 0 =u 0 ¼ 0:05 down to w 0 =u 0 ¼ 0:00, the results presented in this figure will be considered as the benchmark for u 0 ¼ 45 and w 0 =u 0 ¼ 0:00.
It is apparent from this figure that the K ph value is decreasing with increasing ratio of non-associativity (w 0 =u 0 ). This effect is more pronounced for large friction angles. In addition, it is clear from Fig. 3 that simulations with a non-associated flow rule show numerical oscillations, thus, do not converge to a single value after a certain applied displacement. This effect is more obvious with increasing degree of non-associativity K and increasing friction angle. These oscillations are a consequence of the non-associated flow rule.
To enable a fast evaluation of the influence of the flow rule for non-associated materials, results of Fig. 3 have been summarized in Table 2 with respect to the passive earth pressure coefficient K ph normalized by the value obtained with associated flow rule K ph ðw 0 ¼ u 0 Þ. In this table, results have been obtained as minimum and maximum values of the oscillating curves between 80% and 100% of the applied displacement. As the minimum values represent the conservative estimate and, thus, the more relevant results for practical applications, maximum values have been presented in brackets. As can be seen for u 0 ¼ 40 and w 0 =u 0 ¼ 0:00, the simulation with a non-associated flow rule yields a maximum reduction of the passive earth pressure coefficient of more than 40%. Still, for slightly more realistic cases such as u 0 ¼ 40 , w 0 =u 0 ¼ 0:25 (w 0 ¼ 10 ) and u 0 ¼ 30 , w 0 =u 0 ¼ 0:00, maximum reductions of approximately 26% and 15% are observed, respectively.
As the impact of the flow rule has been discussed in this study in terms of the passive earth pressure coefficient and the performance of the simulations, it is still not clear if there is any influence of w 0 \u 0 on the shape of the failure mechanism. Furthermore, it is still a point of discussions whether the failure mechanism varies/changes during the calculation. To investigate both effects, failure mechanisms obtained at four different calculation steps are presented in Fig. 4 for the case of u 0 ¼ 40 and w 0 =u 0 ¼ 0:00 (solid line in Fig. 3e). The calculation steps for this evaluation have been selected at different parts of the load-displacement curve (e.g. different calculation steps), thus, points located close to maximum or minimum values or intermediate points. From Fig. 4 it can be seen that the failure mechanism does not converge to a final form but varies throughout the calculation. The variation can be observed in terms of different geometrical factors, for instance length/extend of the failure mechanism (L x;NA ) or the intersection angles with the ground surface (m 1;NA ) and the bottom of the soil-wall contact surface (m 2;NA ). Furthermore, it can be seen in Fig. 4b and c that the failure surface is not continuous, but shear strains are identified close to the outer slip surface of the failure mechanism and inside the failing soil mass originating from the top left corner of the soil domain, which is the top point of the soil-wall contact.
Comparing the passive earth pressure coefficients of the four calculation steps highlighted in Fig. 3e it is clear from this figure that steps 855 and 1643 yield smaller K ph values than steps 899 and 1372. However, setting these results into context with the different failure mechanisms presented in Fig. 4 it appears that the failure mechanism with the largest lateral extend and the largest outer slip surface (step 1643) leads to the smallest K ph value. Furthermore, it can be recognized that the failure mechanism with a much Table 2 Influence of non-associativity ratio (w 0 =u 0 ) on normalized passive earth pressure coefficient K ph =K ph ðw 0 ¼ u 0 Þ for friction angles between u 0 ¼ 20 and u 0 ¼ 45 and wall friction ratios of d=u 0 ¼ 0:67

Simulations using the Davis approach
As has been shown in the previous section, K ph values are significantly smaller for non-associated materials compared to associated materials. Furthermore, FEA simulations show oscillations in the load-displacement behavior which lead to difficulties in determination of the passive earth pressure coefficient. Therefore, a major objective in this study is to evaluate whether the Davis [6] approach can provide a good estimate of the passive earth pressure coefficient in case of highly non-associated frictional materials. In order to discuss this issue, a parametric study has been conducted using FELA in combination with the Davis [6] approach for a wide range of friction angles, degrees of non-associativity and wall friction ratios. Although differences between FEA and FELA are small for simulations with associated flow rule, which is used in the Davis approach in combination with reduced shear strength parameters, the approximation of the failure mechanism is more accurate in FELA due to the adaptive mesh refinement.
The parametric study has been performed for friction angles in the range of u 0 ¼ 20 À 45 (in 6 steps), for degrees of non-associativity between K ¼ u 0 À w 0 ¼ 0 and K ¼ 45 (in 9 steps with respect to w 0 =u 0 ) and for wall friction ratios in the range of d=u 0 ¼ 0:0 À 1:0 (in 7 steps). Thus, a total number of 378 lower (LB) and upper bound (UB) FELA simulations have been performed and mean values from LB and UB simulations have been calculated. The results are presented in terms of the normalized passive earth pressure coefficient f ¼ K ph =K ph ðw 0 ¼ u 0 Þ (y-axis) with respect to the degree of non-associativity (x-axis) and the wall friction ratio (different curves) in Fig. 5a-f for friction angles of u 0 ¼ 20 À 45 , respectively. In general it can be observed from this figure that f is decreasing with increasing with friction angle, with increasing degree of non-associativity and with increasing wall friction ratio. Analyzing the curve trends and comparing these trends for different subplots, it is apparent that there is a strong influence of the degree of non-associativity on the normalized passive earth pressure coefficient. An influence of the friction angle and the wall friction ratio on f can also be observed, however, their influence is slightly smaller compared to the effect of K on f, especially for small degrees of non-associativity.
Moreover, it is observed for friction angles u 0 35 that all curves presented in Fig. 5a-d follow a concave trend, where the curvature of these concave curves is decreasing with increasing degree of non-associativity and decreasing wall friction ratio. However, for large friction angles and large wall friction ratios (e.g. u 0 ! 40 , d=u 0 ! 0:5) concave curve trends can only be observed up to approximately K % 25 . For larger values of K, the curvature changes from negative to positive and the curve trend is convex.
To evaluate whether the Davis approach provides a useful/realistic estimate for the passive earth pressure problem of a frictional material following a non-associated flow rule, direct comparison of all three solutions (A, NA, Davis) is necessary. Considering the minimum K ph values obtained for the non-associated material in FEA as the benchmark, the over-or underestimation due to an associated flow rule or the Davis approach is evaluated in terms of Ç according to Eq. 5. In this equation, K ph;A=D refers to the passive earth pressure coefficients obtained following an associated flow rule (A) or Davis approach (D) and K ph;NA refers to the passive earth pressure coefficient obtained using a non-associated flow rule (NA). The parameter Ç can then be considered as the relative deviation of the A/D results from the benchmark results (NA) with positive and negative values indicating over-and underestimation, respectively.
The results of the (A, NA and D) simulations presented in Figs. 3 and 5 are evaluated in terms of Ç in Fig. 6 for u 0 ¼ 20 À 45 and w 0 =u 0 ¼ 0:00, 0.25, 0.50 and 0.75 at a constant wall friction ratio of d=u 0 ¼ 0:67. It is apparent from this figure that jÇ j increases with increasing friction angle and increasing degree of non-associativity (decreasing w 0 =u 0 ). This trend holds for all results obtained with an associated flow rule and for most of the results obtained with the Davis approach. However, in case of u 0 ! 35 and w 0 =u 0 ¼ 0:00 the Davis approach shows less underestimation compared to w 0 =u 0 ¼ 0:25. Furthermore, it can be seen that the overestimation of K ph due to associated plasticity is significantly larger than the underestimation of K ph due to Davis approach, especially for large friction angles and small w 0 =u 0 ratios. For instance, comparing Ç values for u 0 ¼ 40 , w 0 =u 0 ¼ 0:25 (w 0 ¼ 10 ) and u 0 ¼ 30 , w 0 =u 0 ¼ 0:00 an overestimation of 35% and 18% is observed using an associated flow rule, whereas an underestimation of 6.5% and 10% is observed using the Davis approach, respectively. Furthermore, it should be noted that, as stated in the first paragraph of Sect. 4.2, the NA benchmark results utilized for u 0 ¼ 45 and w 0 =u 0 ¼ 0:00 have been obtained using w 0 =u 0 ¼ 0:05. Thus, the NA results for w 0 =u 0 ¼ 0:00 may be slightly smaller, which would also lead to an increase of the accuracy of the Davis approach estimate.  5 Influence of degree of non-associativity K on normalized K ph following Davis [6] approach for wall friction ratios between d=/ 0 ¼ 0:0 and d=/ 0 ¼ 1:0 and friction angles of a / 0 ¼ 20 , b / 0 ¼ 25 , c / 0 ¼ 30 , d / 0 ¼ 35 , e / 0 ¼ 40 and f / 0 ¼ 45 Based on the above mentioned observations it can be stated that the Davis approach provides a reasonable estimate for the passive earth pressure coefficient. Although slightly conservative compared to the NA benchmark solution, the K ph value determined using the Davis approach can be considered a more realistic estimate than a K ph value which has been derived using an associated flow rule. As the Davis approach yields good estimates of the passive earth pressure for a wide range of friction angles and degrees of non-associativity, the results presented in Fig. 5 may also serve as design charts to enable a fast determination of passive earth pressure coefficients. Thus, for practical use, K ph values can be taken from Table 1 or other earth pressure tables and corrections can be applied using the f plots given in Fig. 5. Still it should be noted that more FEA simulations with other wall friction ratios (d=u 0 6 ¼ 0:67) are necessary to further evaluate the applicability of the Davis approach for the passive earth pressure problem of non-associated materials.

Further discussion and additional remarks
Prior studies investigated the influence of the flow rule on the stability of slopes as well as the limit load for bearing capacity problems. However, less research was conducted with respect to the importance of this issue for the passive earth pressure problem. The results of this study show that also the determination of passive earth pressure coefficient is significantly influenced by the choice of the flow rule.
Taking another look on the results presented in Table 1, it is interesting to note that although the flow rule is not a direct input, results of Sokolovski [22] can clearly be categorized to an associated flow rule with w 0 ¼ u 0 . As expected, the passive earth pressure coefficients obtained using FELA and FEA for w 0 ¼ u 0 compare well with the approach from Sokolovski [17] when using the same friction angle. However, since dilation angle is not an input parameter, it needs to be questioned why Sokolovski results are similar to the ones obtained for associated flow rule. A possible explanation for this is that in the Sokolovski solution passive earth pressure coefficients are derived based on the method of characteristics considering a stress field which must nowhere violate the Mohr-Coulomb failure criterion. However, since it is generally accepted that the failure surface is better approximated by velocity characteristics rather than by stress characteristics, also shear strength parameters which correspond to the stress state on the velocity characteristics have to be utilized. For this reason, Davis [6] and Roscoe [19] argued that shear strength on a slip or failure surface is smaller than on stress characteristics and recommended to use reduced shear strength parameters. One possibility to account for reduced shear strength on velocity characteristics is to apply the Davis approach using u Ã and c Ã . To evaluate to what extent shear strength needs to be reduced, Roscoe [19] evaluated direct shear tests and confirmed the approach suggested by Davis [6]. Pregl [17] also came to this conclusion, but argued that the influence of not considering this shear strength reduction (strictly speaking the Davis approach) should be small compared to other influences (e.g. shape of the failure surface, distribution of normal forces, estimation of shear strength parameters, Fig. 6 Comparison of relative deviations of K ph values obtained using an associated flow rule and Davis [6] approach from the benchmark defined as K ph values obtained using a non-associated flow rule with (-) and (þ) values indicating under-and overestimation progressive evolution of failure mechanism etc.). Though, the findings of this study indicate that the impact of w 0 on the computed values of K ph may be larger than expected, as can be seen, for instance in Fig. 6.
Following an associated flow rule instead of non-associated flow rule, the passive earth pressure coefficient is significantly over predicted, especially for large friction angles and large degrees of non-associativity. However, it still needs to be investigated whether the impact of the flow rule is smaller for cohesive-frictional soils. For non-associated plasticity the observations match those described in earlier studies as numerical instabilities are likely to occur for materials with large friction angles and low values of w 0 . Consequences of this problem are oscillations in the load-displacement curves as well as multiple failure mechanisms which change during the calculation, making a rigorous assessment of K ph values difficult. Therefore, for simple cases as discussed in this contribution it can be recommended to follow the Davis approach which is in good agreement with FEA applying non-associated plasticity.

Summary and conclusions
In this contribution, the influence of the flow rule on the passive earth pressure problem has been studied. Calculations have been performed using displacement finite-element analysis (FEA) considering associated and nonassociated frictional materials as well as using finite-element limit analysis (FELA) with and without reduced shear strength parameters according to Davis [6]. Results have been presented in terms of passive earth pressure coefficients and compared with selected data from literature. Furthermore, the influence of the flow rule on the passive earth pressure coefficient, the shape of the failure mechanism and the robustness of the numerical simulation has been studied and a parametric study has been conducted with FELA in order to create design charts of a correction factor for the passive earth pressure coefficient in accordance with the Davis approach. Based on the analyses performed in the present study, the following general remarks can be drawn: 1. For an associated flow rule, passive earth pressure coefficients obtained with FEA and FELA are in good agreement with K ph values reported in literature, e.g. [11,22]. 2. Comparing results for an associated and a non-associated flow rule, it is apparent that K ph values are significant smaller for the latter case. This trend is significantly noticeable with increasing difference between internal friction and dilation angle (K ¼ u 0 À w 0 ). Furthermore, FEA simulations conducted for these materials may lead to highly erratic results as well as a variation of the failure mechanism at different calculation steps. 3. Following Davis [6] approach to take into account the non-associativity of real soils, passive earth pressure coefficients are slightly underestimated. Due to the combination of u Ã and c Ã with an associated flow rule, the robustness of the numerical simulations is enhanced compared to calculations following a nonassociated flow rule.
Funding Open Access funding enabled and organized by Projekt DEAL.

Declarations
Conflict of interest The authors declare that they have no conflict of interest.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons. org/licenses/by/4.0/.