微观力学模型启发的颗粒材料各向异性弹塑性本构模型框架

2025/2026学年第一学期期末组会报告

汪海林

同济大学地下建筑与工程系

2026-01-31

大论文框架

flowchart LR
  micromechanics[各向异性微观本构模型] --> chap1[**第一章 绪论**]
  micromechanics --> elastic-stiffness[各向异性弹性刚度] --> chap2[**第二章 各向异性小应变刚度的<br>解析解及离散元验证**]
  micromechanics --> strength-anisotropy[各向异性强度] --> chap3[**第三章 基于微观力学模型的<br>各向异性强度破坏准则**]
  micromechanics --> critical-state[各向异性临界状态] --> chap4[**第四章 微观力学模型的<br>各向异性临界状态理论验证**]
  chap2 --> framework[各向异性弹塑性本构模型框架] --> chap5[**第五章 微观力学模型启发的<br>各向异性弹塑性本构模型框架**]
  chap3 --> framework
  chap4 --> framework
  chap5 --> chap6[**第六章 各向异性弹塑性本构模型<br>在地基承载力分析中的应用**]
  micromechanics --> chap7[**第七章 结论与展望**]

各向异性弹塑性本构模型

各向异性弹性刚度

静态假设假设下的应力应变关系： \[ \small \mathrm{d}\boldsymbol{\varepsilon}_{ij}\!=\!\mathbf{S}_{ijkl}\mathrm{d}\boldsymbol{\sigma}_{kl},\mathbf{S}_{ijkl}\!=\!\dfrac{1}{V}\sum_{c=1}^N{\left(\mathbf{k}_{jl}^c\right)^{-1} \mathbf{l}_m^c \mathbf{l}_n^c \mathbf{A}_{mi} \mathbf{A}_{nk}} \] Voigt 表示法下的柔度矩阵： \[ \small \mathbf{S}_{ij}^\prime=\begin{bmatrix} \dfrac{F_{11}^{[c_1]}}{E} & -\dfrac{\nu}{E}F_{33}^{\left[-\frac{2}{7}\right]} & -\dfrac{\nu}{E}F_{22}^{\left[-\frac{2}{7}\right]} & -\dfrac{6\nu}{7E}F_{12} & \dfrac{c_2}{E}F_{12} & -\dfrac{6\nu}{7E}F_{13} & \dfrac{c_2}{E}F_{13} & -\dfrac{2\nu}{7E}F_{23} & -\dfrac{2\nu}{7E}F_{23} \\ & \dfrac{F_{22}^{[c_1]}}{E} & -\dfrac{\nu}{E}F_{11}^{\left[-\frac{2}{7}\right]} & \dfrac{c_2}{E}F_{12} & -\dfrac{6\nu}{7E}F_{12} & -\dfrac{2\nu}{7E}F_{13} & -\dfrac{2\nu}{7E}F_{13} & -\dfrac{6\nu}{7E}F_{23} & \dfrac{c_2}{E}F_{23} \\ & & \dfrac{F_{33}^{[c_1]}}{E} & -\dfrac{2\nu}{7E}F_{12} & -\dfrac{2\nu}{7E}F_{12} & \dfrac{c_2}{E}F_{13} & -\dfrac{6\nu}{7E}F_{13} & \dfrac{c_2}{E}F_{23} & -\dfrac{6\nu}{7E}F_{23} \\ & & & \dfrac{F_{[1,2]}}{E} & -\dfrac{\nu}{E}F_{33}^{\left[-\frac{2}{7}\right]} & -\dfrac{2\nu}{7E}F_{23} & -\dfrac{2\nu}{7E}F_{23} & -\dfrac{2\nu}{7E}F_{13} & \dfrac{c_3}{E}F_{13} \\ & & & & \dfrac{F_{[2,1]}}{E} & -\dfrac{2\nu}{7E}F_{23} & \dfrac{c_3}{E}F_{23} & -\dfrac{2\nu}{7E}F_{13} & -\dfrac{2\nu}{7E}F_{13} \\ & & & & & \dfrac{F_{[1,3]}}{E} & -\dfrac{\nu}{E}F_{22}^{\left[-\frac{2}{7}\right]} & \dfrac{c_3}{E}F_{12} & -\dfrac{2\nu}{7E}F_{12} \\ & & & & & & \dfrac{F_{[3,1]}}{E} & -\dfrac{2\nu}{7E}F_{12} & -\dfrac{2\nu}{7E}F_{12} \\ & & & & & & & \dfrac{F_{[2,3]}}{E} & -\dfrac{\nu}{E}F_{11}^{\left[-\frac{2}{7}\right]} \\ \text{Sym.} & & & & & & & & \dfrac{F_{[3,2]}}{E} \\ \end{bmatrix} \]

解析解的验证

Figure 1: Comparison between analytical and numerical solutions of elastic properties: (a) contact normal density; (b)-(d) analytical solutions of \(E_{33}/E\), \(G_{31}/G\) and \(\nu_{31}/\nu\); (e)-(g) numerical solutions of \(E_{33}/E\), \(G_{31}/G\) and \(\nu_{31}/\nu\)

各向异性弹塑性模型

根据静态假设，提出以下应力张量，将各向异性材料转化为等效的各向同性材料： \[ \boldsymbol{\sigma}^\prime=\boldsymbol{\sigma}\mathbf{A}\left(\mathbf{A}^{\mathrm{iso}}\right)^{-1}=\boldsymbol{\sigma}\left(\mathbf{I}+\frac{2}{5}\mathbf{F}\right)^{-1}=\boldsymbol{\sigma}\mathbf{A}^{\mathrm{n}} \] 在传统弹塑性模型中使用上述应力张量，屈服函数和塑性势函数对应力的导数为： \[ \dfrac{\partial F}{\partial\sigma_{ij}}=\dfrac{\partial F}{\partial\sigma_{kl}^\prime}\dfrac{\partial\sigma_{kl}^\prime}{\partial\sigma_{ij}}, \dfrac{\partial G}{\partial\sigma_{ij}}=\dfrac{\partial G}{\partial\sigma_{kl}^\prime}\dfrac{\partial\sigma_{kl}^\prime}{\partial\sigma_{ij}} \] 弹塑性刚度矩阵为： \[ D_{ijkl}^{\mathrm{ep}}=D_{ijkl}-\dfrac{D_{ijkl}\dfrac{\partial G}{\partial\sigma_{mn}^\prime}\dfrac{\partial\sigma_{mn}^\prime}{\partial\sigma_{kl}}\dfrac{\partial F}{\partial\sigma_{kl}^\prime}\dfrac{\partial\sigma_{kl}^\prime}{\partial\sigma_{ij}}D_{ijkl}}{\dfrac{\partial F}{\partial\sigma_{kl}^\prime}\dfrac{\partial\sigma_{kl}^\prime}{\partial\sigma_{ij}}D_{ijkl}\dfrac{\partial G}{\partial\sigma_{mn}^\prime}\dfrac{\partial\sigma_{mn}^\prime}{\partial\sigma_{kl}}-\dfrac{\partial F}{\partial\kappa}\dfrac{\partial\kappa}{\partial\varepsilon_{ij}^{\mathrm{p}}}\dfrac{\partial G}{\partial\sigma_{mn}^\prime}\dfrac{\partial\sigma_{mn}^\prime}{\partial\sigma_{kl}}} \]

广义的强度破坏准则

摩尔库伦模型的临界应力比为： \[ M^\prime=\dfrac{6\sqrt{(b^\prime)^2-b^\prime+1}\sin\phi_c}{3+(2b^\prime-1)\sin\phi_c} \] 采用与Lode角相关的摩擦角，可得到广义强度破坏准则： \[ \small \begin{aligned} &\sin\phi_c^\prime=\omega_{\mathrm{MC}}\sin\phi_c+\omega_{\mathrm{Mises}}\sin\phi_{\mathrm{Mises}}+\omega_{\mathrm{Tresca}}\sin\phi_{\mathrm{Tresca}}+\omega_{\mathrm{SMP}}\sin\phi_{\mathrm{SMP}} \\ &\sin\phi_{\mathrm{Mises}}=\dfrac{3\sin\phi_c}{3\sqrt{\left(b^\prime\right)^2-b^\prime+1}-\left(2b^\prime+\sqrt{\left(b^\prime\right)^2-b^\prime+1}-1\right)\sin\phi_c} \\ &\sin\phi_{\mathrm{Tresca}}=\frac{3\sin\phi_c}{3-2b\prime\sin\phi_c} \\ &\sin\phi_{\mathrm{SMP}}=\dfrac{\tan\phi_{\mathrm{SMP}}}{\sqrt{1+\tan^2\phi_{\mathrm{SMP}}}} \text{ with } \tan\phi_{\mathrm{SMP}}=\sqrt{\dfrac{2}{1+\dfrac{(1-b^\prime)^2+\left(b^\prime\right)^2 R^\prime}{1-b^\prime+b^\prime R^\prime}}}\tan\phi_c \end{aligned} \]

广义的强度破坏准则

Figure 2: Effect of anisotropy degree and loading direction on strength anisotropy: (a)-(d) effect of anisotropy degree; (e)-(h) effect of loading direction

应用到一种简单的弹塑性模型

屈服函数： \[ F(\boldsymbol{\sigma}^\prime)=\dfrac{q^\prime}{p^\prime+c\cot\phi_c}-\dfrac{M_\mathrm{p}^\prime\varepsilon_{\mathrm{d}}^{\mathrm{p}}}{k_p+\varepsilon_{\mathrm{d}}^{\mathrm{p}}} \] 应力-剪胀关系： \[ \mathrm{D}^\prime=A_d\left(M_{\mathrm{pt}}^\prime-\dfrac{q^\prime}{p^\prime+c\cot\phi_c}\right), \dfrac{\partial G}{\partial p^\prime}=D^\prime, \dfrac{\partial G}{\partial q^\prime}=1 \] 非线性剪切模量： \[ G=G_0\dfrac{(2.97-e)^2}{1+e}\left(\dfrac{p^\prime+c\cot\phi_c}{p_{\mathrm{a}}}\right)^n \]

模拟结果

各向异性小应变刚度

Figure 3: Comparison between simulations and experimental results (Mitaritonna, Amorosi, and Cotecchia 2014) for anisotropic elastic properties of Lucera clay: (a) stress path; (b) shear stiffness ratio; (c) norm of fabric tensor; (d) Young’s modulus ratio

强度各向异性

Figure 4: Comparison between simulations and experimental results for true triaxial tests: (a) Ochiai and Lade (1983); (b) Kirkgard and Lade (1993); (c) Abelev and Lade (2004); (d) Haruyama (1981)

排水三轴实验

Figure 5: Comparison between simulations and experimental results (Ochiai and Lade 1983) for drained triaxial tests on Cambria sand: (a) deviatoric stress; (b) volumetric strain; (c) fabric anisotropy variable; (d) norm of fabric tensor

各向异性临界状态理论

新的应力张量

基于修正应力张量的塑性应变： \[ \mathrm{d}\boldsymbol{\varepsilon}^{\mathrm{p}}=\mathrm{d}\lambda\dfrac{\partial G}{\partial\sigma_{kl}^\prime}\dfrac{\partial\sigma_{kl}^\prime}{\partial\sigma_{ij}}=\mathrm{d}\lambda\left(\dfrac{\partial G}{\partial p^\prime}\dfrac{\partial p^\prime}{\partial\sigma_{kl}^\prime}\dfrac{\partial\sigma_{kl}^\prime}{\partial\sigma_{ij}}+\dfrac{\partial G}{\partial q^\prime}\dfrac{\partial q^\prime}{\partial\sigma_{kl}^\prime}\dfrac{\partial\sigma_{kl}^\prime}{\partial\sigma_{ij}}\right) \] 以此得到的应力-剪胀关系： \[ \dfrac{\mathrm{d}\varepsilon_{\mathrm{v}}^{\mathrm{p}}}{\mathrm{d}\varepsilon_{\mathrm{d}}^{\mathrm{p}}}=A_d\left(M_{\mathrm{pt}}^\prime-\dfrac{q^\prime}{p^\prime+c\cot\phi_c}\right)\dfrac{\mathrm{tr}\left(\mathbf{A}^{\mathrm{n}}\right)}{3}-{\color{red} {A^{\mathrm{n}}_q}} \] 上式中 \(A_q^n\) 的存在使得体应变会不断累积，即使达到临界状态。

因此，提出以下应力张量： \[ \boldsymbol{\sigma}^*=(1-A)\boldsymbol{\sigma}^\prime+A\boldsymbol{\sigma} \] 各向异性达到临界状态时（\(A=1\)），应力张量退化到原始应力张量。

模拟结果

Figure 6: Comparison between simulations of the anisotropic SIMSAND model with original and modified stress-dilatancy relationships for drained triaxial tests on Cambria sand (thicker lines: coaxial loading with \(\alpha=0^\circ\); thinner lines: non-coaxial loading with \(\alpha=90^\circ\)): (a)-(c) original stress-dilatancy relationship; (d)-(f) modified stress-dilatancy relationship

模拟结果

Figure 7: Comparison between simulations and experimental results (Yoshimine, Ishihara, and Vargas 1998) for torsional shear tests on Toyoura sand: (a) axial strain vs. deviatoric stress; (b) mean stress vs. deviatoric stress; (c) fabric anisotropy variable; (d) norm of fabric tensor

弹塑性模型的应力积分算法

切平面法CPA和最近点投影法CPPM

1. 弹性试算：\(\mathrm{d} \boldsymbol{\varepsilon}_{\mathrm{p}}=\mathbf{0}, \Delta \lambda=0, \boldsymbol{\sigma}=\boldsymbol{\sigma}+\mathbf{K}^{\mathrm{e}} \mathrm{d} \boldsymbol{\varepsilon}\)
2. 检查是否进入塑性区： \[ F=F(\boldsymbol{\sigma})<\operatorname{tol}_1,\|\mathbf{R}\|=\left\|-\mathrm{d} \boldsymbol{\varepsilon}_{\mathrm{p}}+\Delta \lambda \frac{\partial G}{\partial \boldsymbol{\sigma}}\right\|<\operatorname{tol}_2 \Rightarrow \text { 弹性区 } \]
3. 求解塑性乘子： \[ \mathbf{K}^{\mathrm{R}} =\left(\mathbf{I}+\Delta \lambda \mathbf{K}^{\mathrm{e}} \frac{\partial^2 G}{\partial \boldsymbol{\sigma}^2}\right)^{-1} \mathbf{K}^{\mathrm{e}}, \delta \lambda =\frac{F-\left[\frac{\partial F}{\partial \boldsymbol{\sigma}}\right]^T \mathbf{K}^{\mathrm{R}} \frac{\partial G}{\partial \boldsymbol{\sigma}}}{\left[\frac{\partial F}{\partial \boldsymbol{\sigma}}\right]^T \mathbf{K}^{\mathrm{R}} \frac{\partial G}{\partial \boldsymbol{\sigma}}-\sum_\kappa \frac{\partial F}{\partial \kappa} \frac{\partial \kappa}{\partial \boldsymbol{\varepsilon}_{\mathrm{p}}} \frac{\partial G}{\partial \boldsymbol{\sigma}}}, \Delta \lambda=\Delta \lambda+\delta \lambda \]
4. 更新应力应变： \[ \delta \boldsymbol{\sigma}=\mathbf{K}^{\mathrm{R}} \mathbf{R}+\delta \lambda \mathbf{K}^{\mathrm{R}} \frac{\partial G}{\partial \boldsymbol{\sigma}}, \boldsymbol{\sigma}=\boldsymbol{\sigma}-\delta \boldsymbol{\sigma}, \delta \boldsymbol{\varepsilon}_{\mathrm{p}}=\left(\mathbf{K}^{\mathrm{e}}\right)^{-1} \delta \boldsymbol{\sigma}, \mathrm{d} \boldsymbol{\varepsilon}_{\mathrm{p}}=\mathrm{d} \boldsymbol{\varepsilon}_{\mathrm{p}}+\delta \boldsymbol{\varepsilon}_{\mathrm{p}} \]

主要问题：

• 切平面法（目前使用的方法）
  • 无法得到一致切线刚度矩阵的解析解，传入ABAQUS的切线刚度矩阵为弹性刚度矩阵，导致收敛性较差
• 最近点投影法
  • 需要计算塑性势函数对应力的二阶导数，计算复杂度较高。且SIMSAND模型没有显式的塑性势函数表达式，难以实现

计算的主要问题

自定义Mohr-Coulomb模型：

1. 屈服函数： \[ \small \begin{aligned} & F=R_{m c} q-p \tan \phi-c=0 \\ & R_{m c}(\Theta, \phi)=\frac{1}{\sqrt{3} \cos \phi} \sin \left(\Theta+\frac{\pi}{3}\right)+\frac{1}{3} \cos \left(\Theta+\frac{\pi}{3}\right) \tan \phi \end{aligned} \]
2. 塑性势函数： \[ \small \begin{aligned} & G=\sqrt{\left(\left.\epsilon c\right|_0 \tan \psi\right)^2+\left(R_{m w} q\right)^2}-p \tan \psi \\ & R_{m w}(\Theta, e)=\frac{4\left(1-e^2\right) \cos ^2 \Theta+(2 e-1)^2}{2\left(1-e^2\right) \cos \Theta+(2 e-1) \sqrt{4\left(1-e^2\right) \cos ^2 \Theta+5 e^2-4 e}} R_{m c}\left(\frac{\pi}{3}, \phi\right) \\ & R_{m c}\left(\frac{\pi}{3}, \phi\right)=\frac{3-\sin \phi}{6 \cos \phi}, e=\frac{3-\sin \phi}{3+\sin \phi} \end{aligned} \]

• 基于CPA的自定义Mohr-Coulomb模型的收敛性比Abaqus内置Mohr-Coulomb模型差，需要10倍以上的迭代次数才能收敛，且在较大变形时无法收敛。
• 各向异性越大，模型收敛性越差。
• 地面超载较小（γ=16kN/m³，d=0.5m，p=8kPa，试验无超载），模型收敛性较差，应力容易减小到0到拉应力。

Abaqus内置Mohr-Coulomb模型（w=0.5b）

自定义Mohr-Coulomb模型（w=0.5b）

各向异性 SIMSAND模型 （w=0.18b）

基于Runge-Kutta算法的积分算法

1. 弹性试算：\(\boldsymbol{\sigma}^{\text{elastic}}=\boldsymbol{\sigma}+{\color{red}\alpha}\mathbf{K}^{\mathrm{e}} \mathrm{d} \boldsymbol{\varepsilon}\)
2. 塑性积分：\(\boldsymbol{\sigma}=\boldsymbol{\sigma}^{\mathrm{elastic}}+(1-{\color{red}\alpha})\mathbf{K}^{\mathrm{ep}} \mathrm{d} \boldsymbol{\varepsilon}\)

Runge-Kutta算法： \[ \begin{aligned} &\boldsymbol{\sigma}_{n+1}=\boldsymbol{\sigma}_n+\sum_{i=1}^s{b_i\mathbf{K}_i^{\mathrm{ep}}}\mathrm{d} \boldsymbol{\varepsilon} \\ &\mathbf{K}_i^{\mathrm{ep}}=\mathbf{K}^{\mathrm{ep}}\left(\boldsymbol{\varepsilon}_n+c_i\mathrm{d}\boldsymbol{\varepsilon},\boldsymbol{\sigma}_n+\sum_{j=1}^s{a_{ij}\mathbf{K}_j^{\mathrm{ep}}}\mathrm{d} \boldsymbol{\varepsilon}\right)\mathrm{d} \boldsymbol{\varepsilon} \end{aligned} \]

自动误差控制的Runge-Kutta算法¹：

1. 使用\(s\)阶和\(s+1\)阶RK算法：\(\boldsymbol{\sigma}_{n+1}^{(s)}\)和\(\boldsymbol{\sigma}_{n+1}^{(s+1)}\)。
2. 计算误差：\(\text{E}={\left\Vert\boldsymbol{\sigma}_{n+1}^{(s+1)}-\boldsymbol{\sigma}_{n+1}^{(s)}\right\Vert}\Big/{\left\Vert\boldsymbol{\sigma}_{n+1}^{(s+1)}\right\Vert}\)。
3. 自动调整步长使误差\(E\)满足预设的容许误差。

1. Liao, D., & Yang, Z. X. (2021). Effect of fabric anisotropy on bearing capacity and failure mode of strip footing on sand: An anisotropic model perspective. Computers and Geotechnics, 138, 104330. https://doi.org/10.1016/j.compgeo.2021.104330

工作总结与计划

地基承载力的模拟

地基承载力的模拟计划

本学期工作总结及接下来的工作计划

• 工作总结
  • 编写微观力学模型的ABAQUS子程序模拟地基极限承载力（2025年9月-10月）
  • 撰写了微观力学模型启发的颗粒材料各向异性弹塑性本构模型框架的英文文章，目前仍在修改（2025年11月-2026年1月）
• 工作计划
  • 使用显式Runge-Kutta算法结合子步增量法编写ABAQUS子程序，传入加权平均切线刚度到ABAQUS（2026年2月）
  • 完成地基承载力模拟并撰写英文小论文（2026年3月）
  • 撰写大论文（2026年4月）
  • 完成大论文答辩（2026年5-8月）

学术成果

• 汪海林, 刘航宇, 顾晓强, & 宋许根. (2021). 基于多元概率分布模型的珠海黏土多参数预测. 岩土工程学报, 43(S2), 193–196. https://doi.org/10.11779/CJGE2021S2046
• Wang, H.-L., Jin, Y.-F., Yin, Z.-Y., & Gu, X.-Q. (2024). Uncertainty Quantification of Soil–Structure Interface Properties with an Enhanced Hypoplastic Interface Model. International Journal of Geomechanics, 24(6), 04024103. https://doi.org/10.1061/IJGNAI.GMENG-9248
• Wang, H.-L., Yin, Z.-Y., Jin, Y.-F., & Gu, X.-Q. (2024). Evaluation of Soil–Structure Interface Models. Journal of Geotechnical and Geoenvironmental Engineering, 150(7), 04024050. https://doi.org/10.1061/JGGEFK.GTENG-11486
• Wang, H.-L., Yin, Z.-Y., Gu, X.-Q., & Jin, Y.-F. (2024). Evaluation of Soil–Structure Interface Models Considering Cyclic Loading Effect. International Journal for Numerical and Analytical Methods in Geomechanics, 48(17), 4257–4273. https://doi.org/10.1002/nag.3831
• Wang, H.-L., Gu, X., Yin, Z.-Y., & Zhao, C.-F. (2025). Integration of critical state theory into a micromechanical model for granular materials accounting for fabric evolution. Computers and Geotechnics, 186, 107379. https://doi.org/10.1016/j.compgeo.2025.107379
• Wang, H.-L., Gu, X., & Yin, Z.-Y. (2026). A generalized anisotropic failure criterion for granular materials based on a multi-scale micromechanical approach. International Journal of Solids and Structures, 325, 113696. https://doi.org/10.1016/j.ijsolstr.2025.113696
• Wang, H.-L., Gu, X., Yin, Z.-Y., & Zhao, C.-F. (2025). An Anisotropic Micromechanical Model for Granular Materials: Development and Validation of a Novel Contact Law. IACMAG2025.
• 汪海林, 顾晓强, & 尹振宇. (2021). 基于微观力学模型的砂土非共轴特性研究. 力学学报（返修）
• Wang, H.-L., Gu, X., & Yin, Z.-Y. Micromechanical Modelling of Anisotropic Elastic Properties for Granular Materials with Two Enhanced Homogenization Methods, Géotechnique (Under Review)
• Wang, H.-L., Gu, X., & Yin, Z.-Y. A Micromechanics-Inspired Anisotropic Elastoplastic Framework for Granular Materials (Pending Submission)

参考文献

Abelev, Andrei V., and Poul V. Lade. 2004. “Characterization of Failure in Cross-Anisotropic Soils.” Journal of Engineering Mechanics 130 (5): 599–606. https://doi.org/10.1061/(ASCE)0733-9399(2004)130:5(599).

Haruyama, Motohisa. 1981. “Anisotropic Deformation-Strength Characteristics of an Assembly of Spherical Particles Under Three Dimensional Stresses.” Soils and Foundations 21 (4): 41–55. https://doi.org/10.3208/sandf1972.21.4_41.

Kirkgard, M. M., and P. V. Lade. 1993. “Anisotropic Three-Dimensional Behavior of a Normally Consolidated Clay.” Canadian Geotechnical Journal 30 (5): 848–58. https://doi.org/10.1139/t93-075.

Mitaritonna, G., A. Amorosi, and F. Cotecchia. 2014. “Experimental Investigation of the Evolution of Elastic Stiffness Anisotropy in a Clayey Soil.” Géotechnique 64 (6): 463–75. https://doi.org/10.1680/geot.13.P.191.

Ochiai, Hidetoshi, and Poul V. Lade. 1983. “Three-Dimensional Behavior of Sand with Anisotropic Fabric.” Journal of Geotechnical Engineering 109 (10): 1313–28. https://doi.org/10.1061/(ASCE)0733-9410(1983)109:10(1313).

Yoshimine, Mitsutoshi, Kenji Ishihara, and William Vargas. 1998. “Effects of Principal Stress Direction and Intermediate Principal Stress on Undrained Shear Behavior of Sand.” Soils and Foundations 38 (3): 179–88. https://doi.org/10.3208/sandf.38.3_179.