# General-to-Detailed GAN for
Infrequent Class Medical Images

Tatsuki Koga
<EMAIL_ADDRESS>Naoki Nonaka
<EMAIL_ADDRESS>Jun Sakuma
<EMAIL_ADDRESS>Jun Seita
<EMAIL_ADDRESS>Medical Sciences Innovation Hub Program, RIKEN and Center
for Advanced Intelligence Project, RIKENMedical Sciences Innovation Hub
Program, RIKENCenter for Advanced Intelligence Project, RIKEN and Department
of Computer Science, University of TsukubaMedical Sciences Innovation Hub
Program, RIKEN and Center for Advanced Intelligence Project, RIKEN

###### Abstract

Deep learning has significant potential for medical imaging. However, since
the incident rate of each disease varies widely, the frequency of classes in a
medical image dataset is imbalanced, leading to poor accuracy for such
infrequent classes. One possible solution is data augmentation of infrequent
classes using synthesized images created by Generative Adversarial Networks
(GANs), but conventional GANs also require certain amount of images to learn.
To overcome this limitation, here we propose General-to-detailed GAN (GDGAN),
serially connected two GANs, one for general labels and the other for detailed
labels. GDGAN produced diverse medical images, and the network trained with an
augmented dataset outperformed other networks using existing methods with
respect to Area-Under-Curve (AUC) of Receiver Operating Characteristic (ROC)
curve.

##  1 Introduction

Deep convolutional networks have made a breakthrough in the field of image
classification. In the medical imaging field, many applications have been
developed utilizing open datasets such as Chest X-ray dataset from National
Institutes of Health (NIH) [1].
The typical model for medical imaging diagnosis takes images and auxiliary
data as inputs, and outputs whether the data indicates diseases. While the
task itself does not so much differ from the task in other fields, the model
may fail to learn features well due to medical image specific issues. One main
issue causing less accurate models is imbalance of data availability. The
number of data without diseases is much larger than those with diseases. In
addition, while some diseases commonly appear in the dataset, others appear
very rarely. In general, it is harder for models to learn well on the
relatively infrequent classes due to the amount of the information given to
it. To overcome the imbalance, oversampling and undersampling are typical
techniques [2], and derivative algorithms have been proposed [3]. The most
straightforward way to oversample is to randomly duplicate some data with
minor classes and add this to a train dataset. The way to undersample is to
remove some data with major classes from the dataset. In general, oversampling
leads to overfitting and undersampling leads to loss of the information.
Generative Adversarial Networks (GAN), originally proposed by Goodfellow et al
[4], generate realistic synthetic samples from high-dimensional distributions.
Some GANs, such as auxiliary classifier GAN (ACGAN) [5], can generate images
with a given class. However, conventional GANs also require certain amount of
train dataset.
In order to deal with the imbalance of the medical image, we propose in this
paper a General-to-detailed GAN (GDGAN) to generate medical images with
infrequent classes. GDGAN contains two stages of GANs. The first GAN generates
medical images with general classes, such as gender, and the second GAN adds
disease specific detailed classes to the images generated by the first GAN.

##  2 Related Works

Dataset augmentation for non-medical datasets utilizing GANs has been reported
in 2017 [6]. Although the datasets used did not have class imbalance, the
augmentation generally improved the performance of deep convolutional
networks.
The synthesis has been done for relatively small medical datasets as well [7]
[8]. They aimed to increase the total number of data in the dataset. Xin, Ekta
and Paul tackled the imbalance of medical images by using GAN for a future
extraction [9]. However, they tested the effect of GAN only by support vector
machine.

##  3 Method

Figure 1: Flowchart of GDGAN.

In general, medical image datasets contain both labels for general features,
e.g., “male”, “female”, and labels with disease specific detailed features. To
use this nature effectively, we designed a GAN named General-to-detailed GAN
(GDGAN). GDGAN consists of two serially connected GANs. The first GAN takes a
uniform noise distribution Z and general labels, e.g., gender, as inputs, then
generates diverse images. The second GAN takes images generated by the first
GAN, and detailed labels such as disease classes as inputs, and generates the
final images (Figure 1).
The loss function of the first GAN is partly taken from Wasserstein GAN (WGAN)
with gradient penalty [10] [11], and also has the log-likelihood of the
correct basic classes. The second GAN also has the loss of WGAN with gradient
penalty and, in addition to that, it has log-likelihoods of the correct
disease classes and basic classes (using the discriminator of the first GAN),
and mean squared error loss between inputs and outputs.
The Chest X-ray dataset from NIH was used for the experiment. It contains
112,120 images in total and each image has one or more labels: age, gender,
view position and 14 disease-specific phenotypes. We split the dataset into a
train dataset (about 70%), a validation dataset (about 10%) and a test dataset
(about 20%). Two “general” labels, gender and view position, were fed to the
first GAN, and “detailed” labels, such as Cardiomegaly (2,772 images, out 2.5%
of total images), were fed to the second GAN. Images were resized to a size of
$64\times 64$. Each GAN in GDGAN was trained separately. The first GAN was
trained initially, then training of the second GAN followed. For each
training, we used Adam [12] as an optimizer with a learning rate of 0.0002. To
represent single stage GAN, ACGAN was trained, too.
The diversity of images synthesized by ACGAN and GDGAN were compared to the
diversity of input images using the inception score [13], which uses an
Inception v3 Network [14] pre-trained on ImageNet [15].
Following the original proposal [11], 50,000 images were randomly selected
from Chest X-ray dataset, or generated by ACGAN and GDGAN, then split into 10
batches (5,000 images each), and the inception score was computed for each
batch. Welch’s t-test was performed to determine whether the scores were
significantly different.
Next, the VGGNet-19 [16] model was built for a disease classifier, with a
little modification to deal with grayscale input images. The network was
trained with Chest X-ray dataset augmented in five ways: no augmentation
($N$=78,468$\pm$0, $N_{Cardiomegaly}$=1,937$\pm$13 for training), undersampled
($N$=64,276$\pm$243, $N_{Cardiomegaly}$=1,937$\pm$13), oversampled
($N$=100,716$\pm$119, $N_{Cardiomegaly}$=14,103$\pm$62), supplied by ACGAN
($N$=100,716$\pm$119, $N_{Cardiomegaly}$=14,103$\pm$62) and supplied by GDGAN
($N$=100,716$\pm$119, $N_{Cardiomegaly}$=14,103$\pm$62). For each method, we
calculated an Area-Under-Curve (AUC) of ROC. We performed this experiment
three times with different splits of the dataset.

##  4 Results

Figure 2: Randomly selected Chest X-ray dataset image, resized to a size of
$64\times 64$ (left) and synthetic image by ACGAN (center) and GDGAN (right).

###  4.1 Inception score

Table 1: Mean inception scores and their standard deviation (SD) of Chest X-ray dataset and synthetic images by ACGAN and GDGAN (N=10). | Chest X-ray dataset | Images by ACGAN | Images by GDGAN
---|---|---|---
Mean score | 2.1258 | 2.7716 | 2.3220
SD | 0.0103 | 0.0172 | 0.0148

Table 1 shows inception scores using 50,000 images of Chest X-ray dataset and
synthetic datasets, by ACGAN and by GDGAN (examples shown in Figure 2). As the
images were grayscale and only composed of X-ray images, both scores were not
so high as other datasets such as CIFAR-10 [17]. However, the score of
synthetic images by GDGAN was higher than that of Chest X-ray dataset
(p=4.3748e-16). The score of synthetic images by ACGAN was the highest
(p=4.6421e-22 for Chest X-ray dataset and p=9.4031e-22 for synthetic images by
GDGAN).

###  4.2 AUCs of ROC

Figure 3: VGGNet-19 classification ROC curves of Cardiomegaly trained with Chest X-ray dataset augmented in different ways. Each curve represents the result of one of the experiments (blue: 1st, orange: 2nd, green: 3rd). Table 2: Comparison of total number of images for training ($N_{total}$), total number of images with Cardiomegaly ($N_{Cardiomegaly}$), VGGNet-19 classification mean AUC of ROC of Cardiomegaly and their standard deviation (SD) among different augmentation methods (N=3). | Original | Augmentation
---|---|---
| Chest X-ray | Undersample | Oversample | ACGAN | GDGAN
$N_{total}$ | 78,468$\pm$0 | 64,276$\pm$243 | 100,716$\pm$119 | 100,716$\pm$119 | 100,716$\pm$119
$N_{Cardiomegaly}$ | 1,937$\pm$13 | 1,937$\pm$13 | 14,103$\pm$62 | 14,103$\pm$62 | 14,103$\pm$62
Mean AUC | 0.8074 | 0.7919 | 0.8084 | 0.7963 | 0.8207
SD | 0.0177 | 0.0095 | 0.0140 | 0.0212 | 0.0066

Figure 3 demonstrates VGGNet-19 classification ROC curves of Cardiomegaly with
different augmentation methods (no augmentation, undersampling, oversampling,
supplying synthetic images by ACGAN and supplying synthetic images by GDGAN)
to the Chest X-ray dataset. Table 2 shows the total number of images used in
training, the total number of images with Cardiomegaly and the AUC of ROC
(mean and standard deviation) corresponding to Figure 3. In comparison with
the model trained without data augmentations, the models trained with
undersampled Chest X-ray dataset and the dataset with synthesis images by
ACGAN did not outperform. In contrast, the models trained with oversampled
Chest X-ray and the dataset with synthesis images by GDGAN outperformed. And
the model trained with Chest X-ray dataset with synthesis images by GDGAN
scored the best AUC for Cardiomegaly.

##  5 Discussion

The comparison of inception scores revealed that GDGAN successfully generated
infrequent class images as diverse as the original Chest X-ray dataset.
Furthermore, supplying infrequent class images by GDGAN improved the
performance of VGGNet-19 deep convolutional network quantified by AUCs of ROC.
One note is that we only tested relatively low resolution of images ($64\times
64$) and one classification network. When generating and/or supplying higher
resolution images by GDGAN and/or using other deep convolutional networks for
evaluation, performance may alter. Taken together, GDGAN, combination of a GAN
for general classes and a GAN for detailed classes may open up a new strategy
for deep learning with infrequent classes.

#### Acknowledgments

We gratefully acknowledge the support of NVIDIA Corporation for the donation
of the GPU used for this research.

## References

  * [1] Xiaosong Wang, Yifan Peng, Le Lu, Zhiyong Lu, Mohammadhadi Bagheri, and Ronald M Summers.  Chestx-ray8: Hospital-scale chest x-ray database and benchmarks on weakly-supervised classification and localization of common thorax diseases.  In Computer Vision and Pattern Recognition (CVPR), 2017 IEEE Conference on, pages 3462–3471. IEEE, 2017.
  * [2] Vaishali Ganganwar.  An overview of classification algorithms for imbalanced datasets.  International Journal of Emerging Technology and Advanced Engineering, 2(4):42–47, 2012.
  * [3] Hui Han, Wen-Yuan Wang, and Bing-Huan Mao.  Borderline-smote: a new over-sampling method in imbalanced data sets learning.  In International Conference on Intelligent Computing, pages 878–887. Springer, 2005.
  * [4] Ian Goodfellow, Jean Pouget-Abadie, Mehdi Mirza, Bing Xu, David Warde-Farley, Sherjil Ozair, Aaron Courville, and Yoshua Bengio.  Generative adversarial nets.  In Advances in neural information processing systems, pages 2672–2680, 2014.
  * [5] Augustus Odena, Christopher Olah, and Jonathon Shlens.  Conditional image synthesis with auxiliary classifier gans.  arXiv preprint arXiv:1610.09585, 2016.
  * [6] Antreas Antoniou, Amos Storkey, and Harrison Edwards.  Data augmentation generative adversarial networks.  arXiv preprint arXiv:1711.04340, 2017.
  * [7] John T Guibas, Tejpal S Virdi, and Peter S Li.  Synthetic medical images from dual generative adversarial networks.  arXiv preprint arXiv:1709.01872, 2017.
  * [8] Maayan Frid-Adar, Idit Diamant, Eyal Klang, Michal Amitai, Jacob Goldberger, and Hayit Greenspan.  Gan-based synthetic medical image augmentation for increased cnn performance in liver lesion classification.  arXiv preprint arXiv:1803.01229, 2018.
  * [9] Xin Yi, Ekta Walia, and Paul Babyn.  Unsupervised and semi-supervised learning with categorical generative adversarial networks assisted by wasserstein distance for dermoscopy image classification.  arXiv preprint arXiv:1804.03700, 2018.
  * [10] Martin Arjovsky, Soumith Chintala, and Léon Bottou.  Wasserstein gan.  arXiv preprint arXiv:1701.07875, 2017.
  * [11] Ishaan Gulrajani, Faruk Ahmed, Martin Arjovsky, Vincent Dumoulin, and Aaron C Courville.  Improved training of wasserstein gans.  In Advances in Neural Information Processing Systems, pages 5767–5777, 2017.
  * [12] Diederik P Kingma and Jimmy Ba.  Adam: A method for stochastic optimization.  arXiv preprint arXiv:1412.6980, 2014.
  * [13] Tim Salimans, Ian Goodfellow, Wojciech Zaremba, Vicki Cheung, Alec Radford, and Xi Chen.  Improved techniques for training gans.  In Advances in Neural Information Processing Systems, pages 2234–2242, 2016.
  * [14] Christian Szegedy, Vincent Vanhoucke, Sergey Ioffe, Jon Shlens, and Zbigniew Wojna.  Rethinking the inception architecture for computer vision.  In Proceedings of the IEEE conference on computer vision and pattern recognition, pages 2818–2826, 2016.
  * [15] Jia Deng, Wei Dong, Richard Socher, Li-Jia Li, Kai Li, and Li Fei-Fei.  Imagenet: A large-scale hierarchical image database.  In Computer Vision and Pattern Recognition, 2009. CVPR 2009. IEEE Conference on, pages 248–255. Ieee, 2009.
  * [16] Karen Simonyan and Andrew Zisserman.  Very deep convolutional networks for large-scale image recognition.  arXiv preprint arXiv:1409.1556, 2014.
  * [17] Shane Barratt and Rishi Sharma.  A note on the inception score.  arXiv preprint arXiv:1801.01973, 2018.


$\displaystyle\forall\,\tilde{\alpha}\in\mathrm{H}^{1}(\Omega;\mathbb{R}).$

In order to avoid significant deviations from the original problem,
$\eta_{\mathrm{vd}}$ must be chosen as small as possible. In the simulations
presented in section 4, $\Delta t$ is defined as the loading increment, while
$\eta_{\mathrm{vd}}=1\times 10^{-7}$ MPa$\,\cdot\,$s. If required for
convergence in certain time steps, $\eta_{\mathrm{vd}}$ is increased up to a
maximum $\eta_{\mathrm{vd}}=1\times 10^{-5}$ MPa$\,\cdot\,$s. We note that
another possibility to deal with numerical difficulties associated with brutal
crack propagation is the use of a dissipation-based path-following constraint
[123], a topic worth considering in future research.

###  3.5 Spatial discretization

The linearized forms (84) and (105) (or the viscous version (107)) are
suitable for spatial discretization. To this end, we consider standard finite
elements, where the non-local primary fields $\bm{u}$ and $\alpha$ are
interpolated using bilinear shape functions. On the other hand, the local
primary field $\bm{\varepsilon}^{\mathrm{p}}$ is evaluated at Gauss
integration points. This procedure is straightforward and is therefore not
presented for the sake of brevity.

It is worth mentioning that the chosen finite elements cannot describe the
discontinuities embedded in the function spaces (63) and (64). In this
context, the use of discontinuous finite element techniques appears as a more
suitable choice to be considered for future works. Nevertheless, it has been
shown that for ductile phase-field models with perfect plasticity, a strong
concentration of plastic strains is attained upon damage localization,
representing a regularized version of the discontinuous response, with the
element size playing the role of a convergence parameter [78]. As will be
shown in the following section, the same behavior is observed for the present
model. Note that this remark is only relevant for the compressive/shear
regime, which approaches a purely frictional sliding stage with perfect
plasticity as $\alpha\to 1$.

##  4 Numerical simulations

This section presents numerical simulations that highlight the main features
of the model described in section 2.3.3. Specifically, 2D finite element
simulations are performed, aiming to capture different failure modes including
tensile, shear, and mixed-mode fracture. Depending on the example, the results
are compared with numerical and experimental observations, as well as
analytical results from fracture mechanics.

In all examples, plane-strain conditions are assumed. Moreover, a low-level
initial damage $\alpha_{0}=1\times 10^{-5}$ is uniformly distributed in the
domain to allow for plastic-damage evolution in the compressive/shear regime.

###  4.1 Biaxial compression tests

The first example highlights the ability of the model to describe shear
banding and shear fracture. To this end, we consider plane-strain specimens
subjected to biaxial compression. These tests are of particular interest for
the analysis of failure in geomaterials [124, 125, 126]; plane-strain
conditions are encountered often in geotechnical engineering problems, such as
the analysis of underground excavations in rock [125].

[scale=0.95]ps/bicomp_scheme

Figure 6: Schematic representation of the biaxial compression test: (a)
initially homogeneous specimen and (b) specimen with a central hole of $3.4$
mm diameter. The confining pressure is fixed at $p_{0}=5$ MPa during the
displacement loading stage.

For the present study, we consider the two scenarios depicted in figure 6: an
initially homogeneous specimen and a specimen with a central hole,
representative of experiments on sandstone reported in Labuz et al. [125]. In
the simulations, uniform meshes of $\sim$43000 quadrilateral elements with a
characteristic element size $h_{\mathrm{c}}=0.3$ mm are employed. The test is
divided into two loading stages. In the first stage, confining pressure is
gradually applied until $p_{0}=5$ MPa, keeping the vertical displacements at
the bottom edge fixed, where only the center node is also fixed horizontally.
For the second loading stage, the lateral pressure $p_{0}$ is fixed while
vertical displacements $\bar{u}$ are imposed downwards in increments of 0.01
mm.

Labuz et al. [125, 127] report the following experimental material parameters
for the sandstone specimens: Young’s modulus $E=14$ GPa; Poisson’s ratio
$\nu=0.31$; friction angle $\varphi=44^{\circ}$; dilation angle
$\theta=30^{\circ}$, showing non-associative behavior; mode I critical stress
intensity factor $K_{\mathrm{cI}}=0.4$ MPa$\,\cdot\,$m1/2; and mode II
critical stress intensity factor $K_{\mathrm{cII}}=3$ MPa$\,\cdot\,$m1/2. The
same parameters are chosen for the simulations, with the friction and dilation
coefficients computed from $\varphi$ and $\theta$ according to the shear
approximation of the Mohr-Coulomb failure criterion:

$A_{\varphi}=\frac{\sqrt{3}\sin\varphi}{3},\quad
A_{\theta}=\frac{\sqrt{3}\sin\theta}{3}.$ (108)

Moreover, for in-plane self-similar crack growth, $G_{\mathrm{Ic}}$ and
$G_{\mathrm{cII}}$ are computed from $K_{\mathrm{cI}}$ and $K_{\mathrm{cII}}$:

$G_{\mathrm{cI}}=\frac{K_{\mathrm{cI}}^{2}}{E^{\prime}},\quad
G_{\mathrm{cII}}=\frac{K_{\mathrm{cII}}^{2}}{E^{\prime}},$ (109)

with $E^{\prime}=E/(1-\nu^{2})$ for plane strain. Finally, a degradation
parameter $b=25$ and a length scale $\ell=0.97$ mm were chosen. As such,
$\ell/h_{\mathrm{c}}\approx 3.23$, which showed mesh-converging results in
preliminary simulations not reproduced here for the sake of brevity.

Figure 7: Force-displacement curves for the biaxial compression tests: (a)
initially homogeneous specimen and (b) perforated specimen including a
comparison with the experimental curve reported in Labuz et al. [125]. Here,
the specimens are assumed to have an out-of-plane thickness of 100 mm. Figure
8: Shear fracture process in the biaxial compression test for the initially
homogeneous specimen: (a) pre-failure hydrostatic generalized stress
$\operatorname{tr}\bm{s}^{\mathrm{p}}$ [MPa], (b) pre- and post-failure
equivalent plastic strain, and (c) pre- and post-failure damage profile.
Figure 9: Kinematics of the shear fracture process in the biaxial compression
test for the initially homogeneous specimen, showing the pre- and post-failure
values for (a) the horizontal displacement component and (b) the vertical
displacement component.

Figure 7(a) shows the force-displacement curve for the initially homogeneous
specimen, exhibiting failure for $\bar{u}>0.96$ mm. At $\bar{u}=0.96$ mm,
figure 8 shows a homogeneous hydrostatic generalized stress
$\operatorname{tr}\bm{s}^{\mathrm{p}}<0$, as well as homogeneous profiles for
the equivalent plastic strain field $\kappa$ and the damage field $\alpha$.
Further, figure 9 shows smooth and symmetric displacements. At the next time
step, the load abruptly drops, while shear fracture brutally takes place in a
V-shaped pattern. Note from figure 8(b2) that, as expected, the plastic
strains concentrate in a very narrow band. This is also evident from the (mode
II) kinematics shown in figure 9, where the displacement field closely
approximates a jump across the sliding failure planes. At this point, the
material points with $\alpha\approx 1$ approach a perfectly plastic stage that
corresponds to a frictional sliding mechanism between the (approximately)
separate specimen blocks. As such, a residual strength is observed in figure
7(a), modulated by the friction coefficient $A_{\varphi}$ and the confining
pressure.

For the perforated specimen, a lower peak load than that of the initially
homogeneous specimen is observed in figure 7(b), with failure occurring at
$\bar{u}>0.60$ mm. Therein, the force-displacement curve is compared with the
experimental results reported in [125], showing a very close agreement. We
note that, having fixed all material parameters and the fracture length scale,
the curve was tuned by the single degradation parameter $b=25$. At
$\bar{u}=0.60$ mm, figure 10 shows a hydrostatic generalized stress
$\operatorname{tr}\bm{s}^{\mathrm{p}}$ where the tensile region
($\operatorname{tr}\bm{s}^{\mathrm{p}}=0$, i.e., open microcracks) can be
clearly distinguished from the compressive/shear region
($\operatorname{tr}\bm{s}^{\mathrm{p}}<0$, i.e., closed microcracks with
frictional sliding). The corresponding equivalent plastic strain and damage
profiles already hint shear banding, even prior to the peak load, with crack
nucleation at the hole edges at $\bar{u}=0.60$ mm. Nevertheless, the response
remains symmetric and weakly localized. The same observation can be made from
the displacement profiles in figure 10. At the next time step, the load
abruptly drops, while two shear fractures brutally propagate from the hole
along the same failure plane. As in the previous case, figure 10(b2) shows a
strong concentration of plastic strains in a very narrow band along the
failure plane. Accordingly, the displacement profiles in figure 11 clearly
show the kinematics of mode II failure.

Figure 10: Shear fracture process in the biaxial compression test for the
perforated specimen: (a) pre-failure hydrostatic generalized stress
$\operatorname{tr}\bm{s}^{\mathrm{p}}$ [MPa], (b) pre- and post-failure
equivalent plastic strain, and (c) pre- and post-failure damage profile.
Figure 11: Kinematics of the shear fracture process in the biaxial compression
test for the perforated specimen, showing the pre- and post-failure values for
(a) the horizontal displacement component and (b) the vertical displacement
component.

In both the initially homogeneous and the perforated specimens, the failure
plane is found to have an inclination of $62^{\circ}$ with respect to the
horizontal axis, in very close agreement with
${45^{\circ}+(\varphi+\theta)/4=63.5^{\circ}}$, as obtained using bifurcation
analyses of Mohr-Coulomb–based frictional plasticity [14] as well as in
experiments on granular materials [93, 94]. This result further motivates the
need for a non-associative law, highlighting the influence of the dilation
angle in the inclination of the failure plane. Moreover, it is remarkable that
the experimental results of Labuz et al. [125, 127] for these specimens also
report an inclination of about $62^{\circ}$. Nevertheless, it is worth noting
that these experiments showed a kinked fracture with a second, steeper plane
of about $77^{\circ}$. As discussed in Labuz et al. [125], Mohr-Coulomb-based
theories are unable to account for this steep portion, requiring an
unrealistically high friction angle.

###  4.2 Brazilian splitting tests

Having shown the ability of the proposed model to describe shear fracture, we
now turn our attention to the occurrence of tensile fracture under mixed-mode
loading conditions. For this purpose, we consider centrally cracked Brazilian
disk specimens subjected to diametral compression, representing benchmark
experiments in rock-like materials [128, 129, 130, 131, 132, 133, 134].
Particular emphasis is put on the crack propagation angle with respect to the
initial flaw, as typically studied in problems of this type.

[scale=1]ps/brazil

Figure 12: Schematic representation of the Brazilian tests: (a) geometry and
boundary conditions and (b) crack initiation angle. The initial flaw has a
thickness of 0.75 mm, and displacements are imposed on a small region spanning
a circular arc of $8.5^{\circ}$.

Figure 12 shows the geometry and boundary conditions of the pre-cracked disk
specimens. We consider ten different inclinations of the initial flaw, with
$\beta_{0}\in\\{0^{\circ},6.75^{\circ},13.5^{\circ},20.25^{\circ},27^{\circ},40^{\circ},52^{\circ},64^{\circ},76^{\circ},90^{\circ}\\}$.
In the simulations, a mesh of quadrilateral elements with refinement in a
central region was employed, with a minimum characteristic element size
$h_{\mathrm{min}}=0.1$ mm. This element size led to, e.g., $\sim$43000
elements for $\beta_{0}=0^{\circ}$ and $\sim$80000 elements for
$\beta_{0}=27^{\circ}$. In order to avoid numerical difficulties, the
displacements are imposed at the top and bottom caps of the boundary over a
small region, chosen for simplicity as a circular arc of $8.5^{\circ}$.
Vertical displacements $\bar{u}$ are imposed at the top with increments of
$2\times 10^{-4}$ mm, while the bottom cap is fixed vertically, with only the
center node also fixed horizontally.

The material parameters for sandstone specimens are chosen as follows [134]:
Young’s modulus $E=20$ GPa, Poisson’s ratio $\nu=0.25$, mode I fracture
toughness $G_{\mathrm{cI}}=15.68$ N/m, and length scale $\ell=0.5$ mm. We
further assume the degradation parameter $b=1$, the mode II fracture toughness
$G_{\mathrm{cII}}=10\,G_{\mathrm{cI}}$, and the same friction and dilation
parameters from the previous example, i.e., $\varphi=44^{\circ}$ and
$\theta=30^{\circ}$.

The loading mode is characterized by the stress intensity factors
${K_{\mathrm{I}}}$ and ${K_{\mathrm{II}}}$ through the mixity parameter

$M_{\mathrm{e}}=\frac{2}{\pi}\tan^{-1}\bigg{(}\frac{K_{\mathrm{I}}}{K_{\mathrm{II}}}\bigg{)}=\frac{2}{\pi}\tan^{-1}\bigg{(}\frac{Y_{\mathrm{I}}}{Y_{\mathrm{II}}}\bigg{)},$
(110)

where ${Y_{\mathrm{I}}}$ and ${Y_{\mathrm{II}}}$ are dimensionless geometry
parameters depending on $\beta_{0}$ and on the ratio between the diameter of
the specimen and the length of the initial flaw. For the present case,
${Y_{\mathrm{I}}}$ and ${Y_{\mathrm{II}}}$ were extracted from Ayatollahi and
Aliha [129] for all considered $\beta_{0}$ values, with ${Y_{\mathrm{II}}}=0$
for $\beta_{0}=0^{\circ}$ and ${Y_{\mathrm{I}}}=0$ for $\beta_{0}\approx
27^{\circ}$. Thus, pure mode I loading ($M_{\mathrm{e}}=1$) and pure mode II
loading ($M_{\mathrm{e}}=0$) correspond to $\beta_{0}=0^{\circ}$ and
$\beta_{0}=27^{\circ}$, respectively, while intermediate angles
$0^{\circ}<\beta_{0}<27^{\circ}$ imply a transition from mode I to mode II.
Beyond $\beta_{0}=27^{\circ}$, mixed-mode loading is again achieved,
recovering pure mode I loading at $\beta_{0}=90^{\circ}$.

Figure 13 shows the resulting fracture patterns for varying $\beta_{0}$ and
the corresponding initial crack propagation angle $\theta_{0}$. Focusing first
on the transition from mode I to mode II loading (top row), we observe that
cracks nucleate at the tips of the initial flaw, with $\theta_{0}$ increasing
from $0^{\circ}$ at $\beta_{0}=0^{\circ}$ to $61^{\circ}$ at
$\beta_{0}=27^{\circ}$. Thus, pure mode I loading results in self-similar
crack propagation and tensile splitting, while mixed-mode and mode II loading
result in wing-shaped cracks. For higher $\beta_{0}$
($40^{\circ}\leq\beta_{0}\leq 90^{\circ}$; bottom row), crack nucleation
shifts from the tips to the edges of the flaw, with an increasing distance
from the tip as $\beta_{0}$ increases. Moreover, the type of failure shifts
from wing-shaped cracks to tensile splitting, with $\theta_{0}=90^{\circ}$ at
$\beta_{0}=90^{\circ}$.

Figure 13: Fractured specimens showing the final damage profiles for different
load mixities obtained by varying the inclination of the initial flaw. In all
cases, tensile fracture is achieved. Figure 14: Tensile fracture process at
different time steps for the specimen with $\beta_{0}=0^{\circ}$, showing (a)
$\operatorname{tr}\bm{s}^{\mathrm{p}}$ in MPa and (b) the corresponding damage
profile. Figure 15: Tensile fracture process at different time steps for the
specimen with $\beta_{0}=27^{\circ}$, showing (a)
$\operatorname{tr}\bm{s}^{\mathrm{p}}$ in MPa and (b) the corresponding damage
profile. Figure 16: Tensile fracture process at different time steps for the
specimen with $\beta_{0}=90^{\circ}$, showing (a)
$\operatorname{tr}\bm{s}^{\mathrm{p}}$ in MPa and (b) the corresponding damage
profile.

A key observation is that, in all cases shown in figure 13, the resulting
fractures are tensile. Thus, mode I fracture is predicted by the model for all
loading modes, with fracture patterns in agreement with both numerical and
experimental observations reported in the literature [131, 132, 133, 134].
Figures 14, 15, and 16 show the fracture process for $\beta_{0}=0^{\circ}$,
$\beta_{0}=27^{\circ}$, and $\beta_{0}=90^{\circ}$, respectively. In all
cases, at a given time step, a marked tensile region physically linked to
opening microcracks ($\operatorname{tr}\bm{s}^{\mathrm{p}}=0$) can be observed
ahead of the crack tip, paving the way for tensile crack propagation at the
obtained orientations $\theta_{0}$. This region is clearly distinguished from
compressive/shear regions ($\operatorname{tr}\bm{s}^{\mathrm{p}}<0$) in other
parts of the domain, particularly at the caps where displacements are imposed,
which in turn show a Boussinesq-like generalized stress distribution.
Moreover, for $\beta_{0}=27^{\circ}$ and $\beta_{0}=90^{\circ}$, a marked
compressive/shear region is also present at the crack tips, but in directions
different from $\theta_{0}$. These results highlight an important feature of
the model: the ability to preclude unrealistic fracture in compressive regions
without resorting to the usual heuristic energy splits.

Finally, we proceed to assess the initial crack propagation angles in relation
to the orientations predicted by classical fracture criteria, for which we
consider the maximum tangential stress (MTS) criterion [135] and the
generalized maximum tangential stress (GMTS) criterion [136]. In both cases,
mode I crack propagation is predicted in the direction of maximum tangential
stress $\sigma_{\theta\theta}$ near the crack tip. However, the GMTS criterion
employs an approximation of $\sigma_{\theta\theta}$ that retains the non-
singular term, and thus considers the maximum $\sigma_{\theta\theta}$ at a
critical radial distance from the crack tip. As a consequence, the GMTS
criterion takes into account the effect of both the $T$-stress and the
fracture process zone. For pre-cracked Brazilian splitting tests, it has been
shown experimentally [130, 131, 133] and numerically [134] that the MTS
criterion strongly overestimates the initial crack propagation angle, while
the GMTS criterion provides a much better agreement. Therefore, it is
interesting to verify if the present model agrees with these observations.

According to the MTS criterion, the initial crack propagation angle
$\theta_{0}$ is computed from

${Y_{\mathrm{I}}}\sin\theta_{0}+{Y_{\mathrm{II}}}(3\cos\theta_{0}-1)=0.$ (111)

On the other hand, the GMTS can be expressed in terms of dimensionless
parameters as [134]

${Y_{\mathrm{I}}}\sin\theta_{0}+{Y_{\mathrm{II}}}(3\cos\theta_{0}-1)-\frac{64}{3}T^{*}\sqrt{\frac{r_{\mathrm{c}}}{R}}\cos\theta_{0}\sin\frac{\theta_{0}}{2}=0.$
(112)

In this expression, $R=25$ mm denotes the radius of the specimen (figure 12),
$T^{*}$ denotes the dimensionless form of the $T$-stress, and $r_{\mathrm{c}}$
is a critical fracture zone distance, which can be estimated as

$r_{\mathrm{c}}=\frac{1}{2\pi}\bigg{(}\frac{{K}_{\mathrm{cI}}}{\sigma_{\mathrm{t}}}\bigg{)}^{\\!2}.$
(113)

Here, ${K}_{\mathrm{cI}}$ is the mode I fracture toughness, computed from
equation (109). Moreover, $\sigma_{\mathrm{t}}$ is the uniaxial tensile
strength, here computed as in standard phase-field models, i.e.,
$\sigma_{\mathrm{t}}=\sqrt{27\,G_{\mathrm{cI}}E^{\prime}/(256\,\ell)}$ [60,
62].

To compute $\theta_{0}$ from the MTS criterion (111) and the GMTS criterion
(112), the values of $Y_{\mathrm{I}}$, $Y_{\mathrm{II}}$, and $T^{*}$ for
varying $\beta_{0}$ were extracted from previous works [129, 134]. Therein,
$Y_{\mathrm{I}}$, $Y_{\mathrm{II}}$, and $T^{*}$ were computed for the same
geometry (figure 12) using stress intensity factors obtained from finite
element analyses.

The mode mixity for $0^{\circ}\leq\beta_{0}\leq 27^{\circ}$ was computed from
equation (110) and plotted versus the initial crack propagation angle obtained
from the simulations (measured from figure 13), where $M_{\mathrm{e}}=1$ for
$\beta_{0}=0^{\circ}$ (pure mode I), $M_{\mathrm{e}}=0.19$ for
$\beta_{0}=6.75^{\circ}$, $M_{\mathrm{e}}=0.42$ for $\beta_{0}=13.5^{\circ}$,
$M_{\mathrm{e}}=0.68$ for $\beta_{0}=20.25^{\circ}$, and $M_{\mathrm{e}}=0$
for $\beta_{0}=27^{\circ}$ (pure mode II). The results are compared with the
values predicted by both the MTS and the GMTS criteria, showing close
agreement with the latter. As shown in figure 17, this result is consistent
with the simulations performed in Wang and Zhou [134]. Therein, a standard
phase-field model with a spectral energy split [58] was employed, as typically
done in the phase-field literature for Brazilian tests [137]. Consequently, we
conclude that, by virtue of the micromechanics-based free energy (24), the
proposed model can recover the behavior expected from classical mixed-mode
fracture mechanics, as well as the predictions of standard phase-field models
that are enhanced with a typical energy split. Moreover, compared to Wang and
Zhou [134], the profiles shown in figure 13 show much thinner cracks, although
these differences may be attributed to differences in the numerical treatment
of the governing equations.

Figure 17: Crack initiation angles for the different load mixities and
comparison with the MTS and GMTS criteria, as well numerical results reported
in Wang and Zhou [134] for a standard phase-field model with a spectral energy
split.

###  4.3 Brittle-to-ductile transition in a dog-bone–shaped specimen

The last example addresses the response of dog-bone–shaped specimens under
tensile loading at different confining pressure levels. Ramsey and Chester
[138] studied this problem experimentally, revealing that tensile fracture
occurs under little confining pressure, and that a continuous transition to
shear fracture is observed as the confining pressure increases. The problem
was studied numerically under plane-strain conditions in Choo and Sun [3]
using a non-variational phase-field model at finite strains. Therein, a
spectral energy split [58] was employed for the elastic energy density, while
a split based on the Jacobian of the plastic deformation gradient was proposed
for the plastic energy density.

Figure 18 shows the geometry and boundary conditions considered in the present
study. A mesh of $\sim$11500 quadrilateral elements with refinement in the
central region was employed, with a minimum characteristic element size
$h_{\mathrm{min}}=0.15$ mm. The test is divided into two loading stages. In
the first stage, confining pressure is gradually applied up to $p_{0}$,
keeping the vertical displacements at the bottom edge fixed, where only the
center node is also fixed horizontally. For the second loading stage, the
lateral pressure $p_{0}$ is fixed while vertical displacements $\bar{u}$ are
imposed upwards on the top edge. The loading increments are selected depending
on the level of confining pressure, with $5\times 10^{-4}$ mm for $p_{0}\leq
5$ MPa and $5\times 10^{-3}$ mm otherwise.

The material parameters are chosen as follows: Young’s modulus $E=5.96$ GPa,
Poisson’s ratio $\nu=0.15$, mode I fracture toughness $G_{\mathrm{cI}}=16$
N/m, mode II fracture toughness $G_{\mathrm{cII}}=250$ N/m, fracture length
scale $\ell=1$ mm, degradation parameter $b=1$, friction angle
$\varphi=30^{\circ}$, and dilation angle $\theta=25^{\circ}$.

Figure 19 shows the fracture process for $p_{0}=5$ MPa. The hydrostatic
generalized stress $\operatorname{tr}\bm{s}^{\mathrm{p}}$ in figure 19(a)
shows tensile regions ($\operatorname{tr}\bm{s}^{\mathrm{p}}=0$, i.e., open
microcracks) that can be clearly distinguished from compressive/shear regions
($\operatorname{tr}\bm{s}^{\mathrm{p}}<0$, i.e., closed microcracks with
frictional sliding). Accordingly, the equivalent plastic strain in figure
19(b) vanishes in the tensile regions and strongly localizes in the critical
compressive/shear regions, respectively corresponding to tensile fracture and
shear fracture in figure 19(c). At $\bar{u}=0.085$ mm,
$\operatorname{tr}\bm{s}^{\mathrm{p}}$ shows a clear tensile region at the
center of the specimen, while $\kappa$ remains null and $\alpha$ remains low
and homogeneous. In the next time step, a central tensile crack propagates
horizontally in figure 19(c2). At this point, the
$\operatorname{tr}\bm{s}^{\mathrm{p}}$ profile in figure 19(a2) shows
distinctive tensile and compressive/shear regions ahead of the crack tip. Due
to the influence of confining pressure, secondary shear fractures propagate
from the crack tips in the compressive/shear regions. We interpret this
response as a mixed or hybrid fracture pattern, exhibiting tensile (mode I)
and shear (mode II) cracks in different parts of the domain. Thus, we expect a
transition to purely tensile fracture at lower confining pressure, and a
transition to purely shear fracture at higher confining pressure.

[scale=1.25]ps/b2d_scheme

Figure 18: Schematic representation of the tensile test in the dog-bone–shaped
specimen under confining pressure. The cross-section decreases radially from
$20$ mm to $13.6$ mm at the center of the specimen. Figure 19: Fracture
process at different time steps for the dog-bone–shaped specimen with
$p_{0}=5$ MPa, showing (a) $\operatorname{tr}\bm{s}^{\mathrm{p}}$ in MPa, (b)
the corresponding equivalent plastic strains, and (c) the corresponding damage
profiles. Figure 20: Force-displacement diagram for varying confining
pressures, exhibiting brittle-to-ductile transition.

To address this hypothesis, we subject the specimen to varying confining
pressures $0\leq p_{0}\leq 80$ MPa. Figure 20 shows the resulting force-
displacement curves, while figure 21 shows the corresponding profiles for the
hydrostatic generalized stress $\operatorname{tr}\bm{s}^{\mathrm{p}}$, the
equivalent plastic strain field $\kappa$, and the damage field $\alpha$. Note
that the main topological changes in the fracture pattern are observed between
$p_{0}=0$ MPa, $p_{0}=5$ MPa, and $p_{0}=10$ MPa, as expected. At $p_{0}=0$
MPa, figure 20 shows a typical brittle response, while figure 21 shows a mode
I fracture, which propagates brutally in the horizontal direction in the
absence of plastic strains. In contrast, for $p_{0}>5$ MPa, the response
becomes increasingly ductile, capturing the brittle-to-ductile transition.
Thus, figure 20 shows force-displacement curves with residual strength, while
figure 21 shows inclined shear fractures accompanied by localized plastic
strains. As $p_{0}$ increases, the peak load, the corresponding failure
displacement, and the residual strength increase, ending with a very limited
softening stage for $p_{0}=80$ MPa. Accordingly, the plastic strains become
increasingly diffuse (but still localized), while the fracture orientation
angle with respect to the horizontal axis increases. This progressive increase
in the fracture orientation angle can be viewed as a continuous transition
from tensile fracture at low confining pressure to shear fracture at high
confining pressure, in agreement with experimental observations [138].

Figure 21: Fractured specimens showing the post-failure
$\operatorname{tr}\bm{s}^{\mathrm{p}}$ [MPa] profiles (top), the corresponding
equivalent plastic strain field (middle), and the corresponding damage field
(bottom) for different confining pressures (see also figure 20).

Moreover, these results are in agreement with numerical results reported in
the literature for a similar problem [3]. However, in contrast with Choo and
Sun [3], the proposed model is able to predict strongly localized plastic
strains and much more delineated crack profiles, allowing us to clearly
distinguish tensile fracture from shear fracture. Moreover, the brittle-to-
ductile transition includes a clear transitional mode at $p_{0}=5$ MPa. We
further note that most previous works dealing with similar problems adopt a
heuristic split for the strain energy and/or the crack driving force, while in
the proposed model, the different failure modes are the consequence of the
micromechanics-based formulation.

##  5 Conclusions

We have presented a micromechanics-based gradient-damage/phase-field model for
fracture in quasi-brittle geomaterials. Two distinctive behaviors were
embedded in the formulation: a brittle tensile regime, corresponding to the
growth of opening microcracks and mode I fracture, and a ductile
compressive/shear regime, corresponding to the growth of closed microcracks
under frictional sliding and mode II fracture. By virtue of the
micromechanical arguments, the model was constructed with a limited number of
parameters and field variables, all of which can be linked to physical lower-
scale mechanisms. A direct relation was thus established between the gradient-
damage/phase-field variable and a microcrack density parameter, as well as
between plastic strains and the frictional sliding of closed microcracks.
Moreover, the constitutive hardening/softening functions and parameters in the
free energy density were defined as functions of the elastic material
properties and a single degradation function. As a key feature, a non-
associative plasticity law was considered, including the effect of the
dilation angle. This feature was crucial to ensure a non-vanishing energy
dissipation due to frictional sliding. Moreover, the model was constructed in
variational form using the energetic formulation for rate-independent systems,
where a careful treatment of the non-associative law was employed by adopting
a generalized principle of maximum dissipation. Finally, the numerical
implementation procedure was described in detail.

Several numerical simulations were conducted, highlighting the ability of the
model to capture the expected failure modes under different loading
conditions. Biaxial compression tests were performed in specimens with and
without imperfections. The results from these tests showed (mode II) shear
fracture orientations in agreement with experimental observations and
bifurcation analyses of Mohr-Coulomb–based models. Moreover, mixed-mode
loading tests were performed on pre-cracked Brazilian specimens subjected to
diametral compression for a wide range of initial flaw inclinations. The
results showed a first transition from self-similar crack growth at pure mode
I loading to wing-shaped cracks at pure mode II loading, and a second
transition to tensile splitting for higher initial flaw inclinations with
respect to the loading axis. The (mode I) tensile cracks obtained in the pre-
cracked Brazilian tests were in agreement with classical mixed-mode fracture
mechanics, as expected from experimental observations. Finally, dog-
bone–shaped specimens were subjected to tension under different confining
pressure levels. In qualitative agreement with experimental evidence, the
results showed a brittle-to-ductile transition with increasing confining
pressure, from tensile fracture to shear fracture, including a hybrid
transitional mode. We highlight that the results of the different tests were
obtained as a consequence of the micromechanics-based formulation, providing a
physically meaningful alternative to the wide variety of heuristic
modifications that have been proposed for the strain energy density and/or the
crack driving force in the phase-field literature.

The present study paves the way for several topics of future research. As
presented in appendix A.1, a straightforward extension to damage-dependent
friction and dilation allows for decreasing the frictional resistance and
reaching a constant volume state in the post-critical stage. Further, an
enhancement of the formulation to model the response under cyclic loading is
discussed in appendix A.2. Another topic of interest consists of relaxing the
assumption of isotropic behavior, which may be achieved by taking the more
general micromechanical framework with multiple crack families [139] as a
point of departure, and extending the phase-field model to multiple damage
variables (cf. [110]) and multiple plastic strains. Other topics include
performing 3D simulations, e.g., to capture the brittle-to-ductile transition
in triaxial cylinders, the inclusion of a compression cap mechanism to model
plastic compaction, and the extension to fluid-driven fracture in porous
media. Finally, a thorough experimental verification study is worth
considering to quantitatively assess the performance of the model in relation
to standard test procedures in both plane-strain and triaxial conditions.

## Appendix

##  Appendix A On possible extensions

This appendix presents extensions of the proposed model that may prove useful
for future work.

###  A.1 Variable friction and dilation

So far, the friction and dilation coefficients have been considered constant.
A straightforward extension is to consider functions $\alpha\mapsto
A_{\varphi}(\alpha)$ and $\alpha\mapsto A_{\theta}(\alpha)$, with
$A_{\Box}(0)=A_{\Box}^{\mathrm{peak}}$ and
$A_{\Box}(1)=A_{\Box}^{\mathrm{res}}$. A simple example is

$A_{\varphi}(\alpha)\coloneqq
A_{\varphi}^{\mathrm{res}}-(1-\alpha)^{2}(A_{\varphi}^{\mathrm{res}}-A_{\varphi}^{\mathrm{peak}}),\quad
A_{\theta}(\alpha)\coloneqq
A_{\theta}^{\mathrm{res}}-(1-\alpha)^{2}(A_{\theta}^{\mathrm{res}}-A_{\theta}^{\mathrm{peak}}).$
(114)

The homogeneous response for such a model is shown in figure 22 with fixed
$A_{\varphi}^{\mathrm{peak}}=0.15$ and $A_{\theta}^{\mathrm{peak}}=0.0375$,
while $A_{\varphi}^{\mathrm{res}}$ and $A_{\theta}^{\mathrm{res}}$ vary. The
main observation is that as $A_{\varphi}^{\mathrm{res}}\to 0$, the residual
strength vanishes and the stress drops to the imposed confining pressure. We
interpret this response as the degradation of the surface asperities in
frictional contact. On the other hand, as $A_{\theta}^{\mathrm{res}}\to 0$,
the volumetric strain approaches a constant value resembling a _critical
state_ , a crucial concept in soil mechanics [140].

Figure 22: Homogeneous response for uniaxial compression with different
residual parameters $A_{\varphi}^{\mathrm{res}}$ and
$A_{\theta}^{\mathrm{res}}$: (a) axial strain vs. axial stress curves and (b)
the corresponding total volumetric strains. The responses are shown for
$A_{\Box}^{\mathrm{res}}\in\\{A_{\Box}^{\mathrm{peak}},0.5\,A_{\Box}^{\mathrm{peak}},0\\}$.

###  A.2 Opening/closure transition for cyclic loading

As presented in this paper, the proposed model relies on the opening/closure
condition (36) to distinguish between the tensile and compressive/shear
regimes. Thereby, $\operatorname{tr}\bm{s}^{\mathrm{p}}$ is known from the
solution of the plasticity evolution equations (table 2). Under monotonic
loading, the solution of these equations properly predicts
$\operatorname{tr}\bm{s}^{\mathrm{p}}=0$ and
$\operatorname{tr}\bm{s}^{\mathrm{p}}<0$ in the tensile and compressive/shear
regimes, respectively. The same observation can be made about unloading from a
compressive/shear state and further loading in the tensile regime.

Consider now the case in which the material is first loaded in tension and
then unloaded, as shown in figure 23 (points A, B, and C). Because the tensile
regime is brittle, we expect the unloading branch to return to the origin in
figure 23(a). However, the plasticity evolution equations predict
$\operatorname{tr}\bm{s}^{\mathrm{p}}<0$ as soon as unloading takes place. As
a consequence, the unloading branch follows the dotted green line and yields
residual strains. The reason for this unexpected behavior is that, according
to the Drucker-Prager model, $\operatorname{tr}\bm{\varepsilon}^{\mathrm{p}}$
is irreversible. We thus conclude that, for cyclic loading, the KKT plasticity
conditions should not be verified in the tensile regime. Of course, this
renders the opening/closure condition (36) insufficient to characterize the
response. We propose to remedy this issue by complementing condition (36) with
a new history parameter $\bar{\rho}$ intended to track the microcrack opening
state. To this end, during the tensile regime, we prescribe the microcrack
opening measure $\operatorname{tr}\bm{\varepsilon}^{\mathrm{p}}$ from equation
(38) and define

$\bar{\rho}(\bm{x},t)\coloneqq\begin{dcases}\int_{t_{*}}^{t}\operatorname{tr}\dot{\bm{\varepsilon}}^{\mathrm{p}}(\bm{x},s)\,\mathrm{d}s,\quad\operatorname{tr}\bm{\varepsilon}^{\mathrm{p}}=\big{(}1-g_{K}(\alpha)\big{)}\operatorname{tr}\bm{\varepsilon}&\quad\text{if
\sf{open}},\\\ 0&\quad\text{if \sf{closed}}.\end{dcases}$ (115)

Here, $t_{*}$ corresponds to the last time step in which a transition from the
compressive/shear to the tensile regime took place. Therefore, in a full
loading-unloading process in tension, the integral on top vanishes. Thus, we
prevent the verification of the KKT plasticity conditions when $\bar{\rho}>0$.
By incorporating this condition in the solution of the governing equations,
unloading from the tensile regime returns to the origin (point C) in figure
23(a), as expected. For illustrative purposes, the material is further loaded
in compression up to point E, and then unloaded elastically to point F. Beyond
this point, the stress path approaches the apex at point G, signaling the
transition to the tensile regime. Then, the material is fully damaged in
tension.

For simplicity of presentation, and in view of the monotonic loading
conditions considered in the numerical examples, condition (115) was not
considered in the present work, but it seems to be crucial for computing the
response under cyclic loading. This topic will be considered in future
studies. Moreover, the effect of residual stress [141] in the response emerges
as an interesting topic for further research.

Figure 23: Homogeneous response for uniaxial cyclic loading: (a) axial stress
vs. axial strain curves and (b) the corresponding stress path plotted in
$\big{(}\sqrt{3}\operatorname{tr}\bm{s}^{\mathrm{p}}/3,\mathrm{sign}(s^{\mathrm{p}}_{\mathrm{dev}\,zz})\|\bm{s}^{\mathrm{p}}_{\mathrm{dev}}\|\big{)}$
space. The dotted green line corresponds to unloading from the tensile stage
if condition (115) is not considered.

## References

  * Andrieux et al. [1986] S. Andrieux, Y. Bamberger, and J.-J. Marigo.  Un modèle de matériau microfissuré pour les bétons et les roches.  _Journal of Theoretical and Applied Mechanics_ , 5(3):471–513, 1986.
  * Borja and Aydin [2004] R.I. Borja and A. Aydin.  Computational modeling of deformation bands in granular media. I. Geological and mathematical framework.  _Computer Methods in Applied Mechanics and Engineering_ , 193(27-29):2667–2698, 2004.
  * Choo and Sun [2018] J. Choo and W. Sun.  Coupled phase-field and plasticity modeling of geological materials: From brittle fracture to ductile flow.  _Computer Methods in Applied Mechanics and Engineering_ , 330:1–32, 2018.
  * Zhu et al. [2011] Q.-Z. Zhu, J.-F. Shao, and D. Kondo.  A micromechanics-based thermodynamic formulation of isotropic damage with unilateral and friction effects.  _European Journal of Mechanics-A/Solids_ , 30(3):316–325, 2011.
  * Hajiabdolmajid et al. [2002] V. Hajiabdolmajid, P.K. Kaiser, and C.D. Martin.  Modelling brittle failure of rock.  _International Journal of Rock Mechanics and Mining Sciences_ , 39(6):731–741, 2002.
  * Renani and Martin [2018] H.R. Renani and C.D. Martin.  Cohesion degradation and friction mobilization in brittle failure of rocks.  _International Journal of Rock Mechanics and Mining Sciences_ , 106:1–13, 2018.
  * Peng and Cai [2019] J. Peng and M. Cai.  A cohesion loss model for determining residual strength of intact rocks.  _International Journal of Rock Mechanics and Mining Sciences_ , 119:131–139, 2019.
  * Zhao et al. [2018a] L.-Y. Zhao, Q.-Z. Zhu, and J.-F. Shao.  A micro-mechanics based plastic damage model for quasi-brittle materials under a large range of compressive stress.  _International Journal of Plasticity_ , 100:156–176, 2018a.
  * Pensée et al. [2002] V. Pensée, D. Kondo, and L. Dormieux.  Micromechanical analysis of anisotropic damage in brittle materials.  _Journal of Engineering Mechanics_ , 128(8):889–897, 2002.
  * Zhu et al. [2008] Q.-Z. Zhu, D. Kondo, and J.-F. Shao.  Micromechanical analysis of coupling between anisotropic damage and friction in quasi brittle materials: role of the homogenization scheme.  _International Journal of Solids and Structures_ , 45(5):1385–1405, 2008.
  * Zhu et al. [2016] Q.-Z. Zhu, L.-Y. Zhao, and J.-F. Shao.  Analytical and numerical analysis of frictional damage in quasi brittle materials.  _Journal of the Mechanics and Physics of Solids_ , 92:137–163, 2016.
  * Xie et al. [2011] N. Xie, Q.-Z. Zhu, L.-H. Xu, and J.-F. Shao.  A micromechanics-based elastoplastic damage model for quasi-brittle rocks.  _Computers and Geotechnics_ , 38(8):970–977, 2011\.
  * Jia et al. [2020] C. Jia, Q. Zhang, and S. Wang.  Experimental investigation and micromechanical modeling of elastoplastic damage behavior of sandstone.  _Materials_ , 13(15):3414, 2020.
  * Vermeer and de Borst [1984] P.A. Vermeer and R. de Borst.  Non-associated plasticity for soils, concrete and rock.  _HERON, 29 (3), 1984_ , 1984.
  * Lubliner et al. [1989] J. Lubliner, J. Oliver, S. Oller, and E. Oñate.  A plastic-damage model for concrete.  _International Journal of solids and Structures_ , 25(3):299–326, 1989.
  * Chazallon and Hicher [1998] C. Chazallon and P.Y. Hicher.  A constitutive model coupling elastoplasticity and damage for cohesive-frictional materials.  _Mechanics of Cohesive-frictional Materials_ , 3(1):41–63, 1998.
  * Chiarelli et al. [2003] A.-S. Chiarelli, J.-F. Shao, and N. Hoteit.  Modeling of elastoplastic damage behavior of a claystone.  _International Journal of Plasticity_ , 19(1):23–45, 2003.
  * Salari et al. [2004] M.R. Salari, S. Saeb, K.J. Willam, S.J. Patchet, and R.C. Carrasco.  A coupled elastoplastic damage model for geomaterials.  _Computer Methods in Applied Mechanics and Engineering_ , 193(27-29):2625–2643, 2004.
  * Wu et al. [2006] J.-Y. Wu, J. Li, and R. Faria.  An energy release rate-based plastic-damage model for concrete.  _International Journal of Solids and Structures_ , 43(3-4):583–612, 2006.
  * Parisio et al. [2015] F. Parisio, S. Samat, and L. Laloui.  Constitutive analysis of shale: a coupled damage plasticity approach.  _International Journal of Solids and Structures_ , 75:88–98, 2015.
  * Lanoye et al. [2013] E. Lanoye, F. Cormery, D. Kondo, and J.-F. Shao.  An isotropic unilateral damage model coupled with frictional sliding for quasi-brittle materials.  _Mechanics Research Communications_ , 53:31–35, 2013.
  * Marigo and Kazymyrenko [2019] J.-J. Marigo and K. Kazymyrenko.  A micromechanical inspired model for the coupled to damage elasto-plastic behavior of geomaterials under compression.  _Mechanics & Industry_, 20(1):105, 2019.
  * Rudnicki and Rice [1975] J.W. Rudnicki and J.R. Rice.  Conditions for the localization of deformation in pressure-sensitive dilatant materials.  _Journal of the Mechanics and Physics of Solids_ , 23(6):371–394, 1975.
  * Needleman [1979] A. Needleman.  Non-normality and bifurcation in plane strain tension and compression.  _Journal of the Mechanics and Physics of Solids_ , 27(3):231–254, 1979.
  * Sabet and de Borst [2019] S.A. Sabet and R. de Borst.  Structural softening, mesh dependence, and regularisation in non-associated plastic flow.  _International Journal for Numerical and Analytical Methods in Geomechanics_ , 43(13):2170–2183, 2019.
  * Jirásek [2000] M. Jirásek.  Comparative study on finite elements with embedded discontinuities.  _Computer Methods in Applied Mechanics and Engineering_ , 188(1-3):307–330, 2000.
  * Ortiz et al. [1987] M. Ortiz, Y. Leroy, and A. Needleman.  A finite element method for localized failure analysis.  _Computer Methods in Applied Mechanics and Engineering_ , 61(2):189–214, 1987.
  * Belytschko et al. [1988] T. Belytschko, Jacob F., and Bruce E. E.  A finite element with embedded localization zones.  _Computer methods in applied mechanics and engineering_ , 70(1):59–89, 1988.
  * Oliver et al. [1999] J. Oliver, M. Cervera, and O. Manzoli.  Strong discontinuities and continuum plasticity models: the strong discontinuity approach.  _International Journal of Plasticity_ , 15(3):319–351, 1999.
  * Oliver et al. [2004] J. Oliver, A.E. Huespe, E. Samaniego, and E.W.V. Chaves.  Continuum approach to the numerical simulation of material failure in concrete.  _International Journal for Numerical and Analytical Methods in Geomechanics_ , 28(7-8):609–632, 2004.
  * Regueiro and Borja [2001] R.A. Regueiro and R.I. Borja.  Plane strain finite element analysis of pressure sensitive plasticity with strong discontinuity.  _International Journal of Solids and Structures_ , 38(21):3647–3672, 2001.
  * Foster et al. [2007] C.D. Foster, R.I. Borja, and R.A. Regueiro.  Embedded strong discontinuity finite elements for fractured geomaterials with variable friction.  _International Journal for Numerical Methods in Engineering_ , 72(5):549–581, 2007.
  * Chen et al. [2011] Q. Chen, J.E. Andrade, and E. Samaniego.  AES for multiscale localization modeling in granular media.  _Computer Methods in Applied Mechanics and Engineering_ , 200(33):2473–2482, 2011.
  * Wu and Cervera [2016] J.-Y. Wu and M. Cervera.  A thermodynamically consistent plastic-damage framework for localized failure in quasi-brittle solids: Material model and strain localization analysis.  _International Journal of Solids and Structures_ , 88:227–247, 2016.
  * Zhao et al. [2018b] L.-Y. Zhao, J.-F. Shao, and Q.-Z. Zhu.  Analysis of localized cracking in quasi-brittle materials with a micro-mechanics based friction-damage approach.  _Journal of the Mechanics and Physics of Solids_ , 119:163–187, 2018b.
  * Oka et al. [1995] F. Oka, T. Adachi, and A. Yashima.  A strain localization analysis using a viscoplastic softening model for clay.  _International Journal of Plasticity_ , 11(5):523–545, 1995.
  * Cervera et al. [1996] M. Cervera, J. Oliver, and O. Manzoli.  A rate-dependent isotropic damage model for the seismic analysis of concrete dams.  _Earthquake engineering & structural dynamics_, 25(9):987–1010, 1996.
  * de Borst and Duretz [2020] R. de Borst and T. Duretz.  On viscoplastic regularisation of strain-softening rocks and soils.  _International Journal for Numerical and Analytical Methods in Geomechanics_ , 44(6):890–903, 2020.
  * Mühlhaus and Vardoulakis [1987] H.-B. Mühlhaus and I. Vardoulakis.  The thickness of shear bands in granular materials.  _Geotechnique_ , 37(3):271–283, 1987.
  * de Borst and Sluys [1991] R. de Borst and L.J. Sluys.  Localisation in a cosserat continuum under static and dynamic loading conditions.  _Computer Methods in Applied Mechanics and Engineering_ , 90(1-3):805–827, 1991.
  * Bažant and Pijaudier-Cabot [1988] Z.P. Bažant and G. Pijaudier-Cabot.  Nonlocal continuum damage, localization instability and convergence.  _Journal of Applied Mechanics_ , 55(2):287–293, 1988.  ISSN 0021-8936.
  * Bažant and Jirásek [2002] Z.P. Bažant and M. Jirásek.  Nonlocal integral formulations of plasticity and damage: survey of progress.  _Journal of Engineering Mechanics_ , 128(11):1119–1149, 2002.
  * Jirásek [2004] M. Jirásek.  Non-local damage mechanics with application to concrete.  _Revue française de Génie Civil_ , 8(5-6):683–707, 2004.
  * Grassl and Jirásek [2006] P. Grassl and M. Jirásek.  Plastic model with non-local damage applied to concrete.  _International Journal for Numerical and Analytical Methods in Geomechanics_ , 30(1):71–90, 2006.
  * Yoshioka et al. [2019] K. Yoshioka, F. Parisio, D. Naumov, R. Lu, O. Kolditz, and T. Nagel.  Comparative verification of discrete and smeared numerical approaches for the simulation of hydraulic fracturing.  _GEM-International Journal on Geomathematics_ , 10(1):13, 2019.
  * Mühlhaus and Aifantis [1991] H.-B. Mühlhaus and E.C. Aifantis.  A variational principle for gradient plasticity.  _International Journal of Solids and Structures_ , 28(7):845–857, 1991.
  * De Borst and Pamin [1996] R. De Borst and J. Pamin.  Gradient plasticity in numerical simulation of concrete cracking.  _European Journal of Mechanics-A/Solids_ , 15:295–320, 1996\.
  * Peerlings et al. [1998] R.H.J. Peerlings, R. de Borst, W.A.M. Brekelmans, and M.G.D. Geers.  Gradient-enhanced damage modelling of concrete fracture.  _Mechanics of Cohesive-frictional Materials_ , 3(4):323–342, 1998.
  * Comi [1999] C. Comi.  Computational modelling of gradient-enhanced damage in quasi-brittle materials.  _Mechanics of Cohesive-frictional Materials_ , 4(1):17–36, 1999.
  * Pamin et al. [2003] J. Pamin, H. Askes, and R. de Borst.  Two gradient plasticity theories discretized with the element-free galerkin method.  _Computer Methods in Applied Mechanics and Engineering_ , 192(20-21):2377–2403, 2003.
  * Zreid and Kaliske [2018] I. Zreid and M. Kaliske.  A gradient enhanced plasticity–damage microplane model for concrete.  _Computational Mechanics_ , 62(5):1239–1257, 2018\.
  * Hageman et al. [2021] T. Hageman, S.A. Sabet, and R. de Borst.  Convergence in non-associated plasticity and fracture propagation for standard, rate-dependent, and cosserat continua.  _International Journal for Numerical Methods in Engineering_ , 122(3):777–795, 2021.
  * Ulloa et al. [2021a] J. Ulloa, R. Alessi, J. Wambacq, G. Degrande, and S. François.  On the variational modeling of non-associative plasticity.  _International Journal of Solids and Structures_ , 217:272–296, 2021a.
  * Bourdin et al. [2000] B. Bourdin, G.A. Francfort, and J.-J. Marigo.  Numerical experiments in revisited brittle fracture.  _Journal of the Mechanics and Physics of Solids_ , 48(4):797–826, 2000.
  * Francfort and Marigo [1998] G.A. Francfort and J.-J. Marigo.  Revisiting brittle fracture as an energy minimization problem.  _Journal of the Mechanics and Physics of Solids_ , 46(8):1319–1342, 1998.
  * Bourdin et al. [2008] B. Bourdin, G.A. Francfort, and J.-J. Marigo.  The variational approach to fracture.  _Journal of Elasticity_ , 91(1-3):5–148, 2008\.
  * Amor et al. [2009] H. Amor, J.-J. Marigo, and C. Maurini.  Regularized formulation of the variational brittle fracture with unilateral contact: Numerical experiments.  _Journal of the Mechanics and Physics of Solids_ , 57(8):1209–1229, 2009.
  * Miehe et al. [2010a] C. Miehe, M. Hofacker, and F. Welschinger.  A phase field model for rate-independent crack propagation: Robust algorithmic implementation based on operator splits.  _Computer Methods in Applied Mechanics and Engineering_ , 199(45):2765–2778, 2010a.
  * Pham et al. [2011] K. Pham, H. Amor, J.-J. Marigo, and C. Maurini.  Gradient damage models and their use to approximate brittle fracture.  _International Journal of Damage Mechanics_ , 20(4):618–652, 2011.
  * Marigo et al. [2016] J.-J. Marigo, C. Maurini, and K. Pham.  An overview of the modelling of fracture by gradient damage models.  _Meccanica_ , 51(12):3107–3128, 2016.
  * Kristensen et al. [2021] P.K. Kristensen, C.F. Niordson, and E. Martínez-Pañeda.  An assessment of phase field fracture: crack initiation and growth.  _Philosophical Transactions of the Royal Society A_ , 379(2203):20210021, 2021.
  * Tanné et al. [2018] E. Tanné, T. Li, B. Bourdin, J.-J. Marigo, and C. Maurini.  Crack nucleation in variational phase-field models of brittle fracture.  _Journal of the Mechanics and Physics of Solids_ , 110:80–99, 2018.
  * De Lorenzis and Maurini [2021] L. De Lorenzis and C. Maurini.  Nucleation under multi-axial loading in variational phase-field models of brittle fracture.  _International Journal of Fracture_ , May 2021.
  * Miehe et al. [2010b] C. Miehe, F. Welschinger, and M. Hofacker.  Thermodynamically consistent phase-field models of fracture: Variational principles and multi-field FE implementations.  _International Journal for Numerical Methods in Engineering_ , 83(10):1273–1311, 2010b.
  * van Dijk et al. [2020] N.P. van Dijk, J.J. Espadas-Escalante, and P. Isaksson.  Strain energy density decompositions in phase-field fracture theories for orthotropy and anisotropy.  _International Journal of Solids and Structures_ , 196:140–153, 2020.
  * Lancioni and Royer-Carfagni [2009] G. Lancioni and G. Royer-Carfagni.  The variational approach to fracture mechanics. a practical application to the french panthéon in paris.  _Journal of Elasticity_ , 95(1-2):1–30, 2009\.
  * Freddi and Royer-Carfagni [2010] F. Freddi and G. Royer-Carfagni.  Regularized variational theories of fracture: a unified approach.  _Journal of the Mechanics and Physics of Solids_ , 58(8):1154–1174, 2010.
  * Freddi and Royer-Carfagni [2016] F. Freddi and G. Royer-Carfagni.  Phase-field slip-line theory of plasticity.  _Journal of the Mechanics and Physics of Solids_ , 94:257–272, 2016.
  * Ambati et al. [2015a] M. Ambati, T. Gerasimov, and L. De Lorenzis.  A review on phase-field models of brittle fracture and a new fast hybrid formulation.  _Computational Mechanics_ , 55(2):383–405, 2015a.
  * Wu and Nguyen [2018] J.-Y. Wu and V.P. Nguyen.  A length scale insensitive phase-field damage model for brittle fracture.  _Journal of the Mechanics and Physics of Solids_ , 119:20–42, 2018.
  * Kumar et al. [2020] A. Kumar, B. Bourdin, G.A. Francfort, and O. Lopez-Pamies.  Revisiting nucleation in the phase-field approach to brittle fracture.  _Journal of the Mechanics and Physics of Solids_ , 142:104027, 2020.
  * Shen and Stephansson [1994] B. Shen and O. Stephansson.  Modification of the G-criterion for crack propagation subjected to compression.  _Engineering Fracture Mechanics_ , 47(2):177–189, 1994.
  * Zhang et al. [2017] X. Zhang, S.W. Sloan, C. Vignes, and D. Sheng.  A modification of the phase-field model for mixed mode crack propagation in rock-like materials.  _Computer Methods in Applied Mechanics and Engineering_ , 322:123–136, 2017.
  * Spetz et al. [2020] A. Spetz, R. Denzer, E. Tudisco, and O. Dahlblom.  Phase-field fracture modelling of crack nucleation and propagation in porous rock.  _International Journal of Fracture_ , 224(1):31–46, 2020.
  * Bryant and Sun [2018] E.C. Bryant and W. Sun.  A mixed-mode phase field fracture model in anisotropic rocks with consistent kinematics.  _Computer Methods in Applied Mechanics and Engineering_ , 342:561–584, 2018.
  * Fei and Choo [2021] F. Fei and J. Choo.  Double-phase-field formulation for mixed-mode fracture in rocks.  _Computer Methods in Applied Mechanics and Engineering_ , 376:113655, 2021.
  * Fei and Choo [2020] F. Fei and J. Choo.  A phase-field model of frictional shear fracture in geologic materials.  _Computer Methods in Applied Mechanics and Engineering_ , 369:113265, 2020.
  * Alessi et al. [2015] R. Alessi, J.-J. Marigo, and S. Vidoli.  Gradient damage models coupled with plasticity: variational formulation and main properties.  _Mechanics of Materials_ , 80:351–367, 2015.
  * Kuhn et al. [2016] C. Kuhn, T. Noll, and R. Müller.  On phase field modeling of ductile fracture.  _GAMM-Mitteilungen_ , 39(1):35–54, 2016.
  * Miehe et al. [2016] C. Miehe, S. Teichtmeister, and F. Aldakheel.  Phase-field modelling of ductile fracture: a variational gradient-extended plasticity-damage theory and its micromorphic regularization.  _Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Science_ , 374(2066):20150170, 2016.
  * Rodríguez et al. [2018] P. Rodríguez, J. Ulloa, C. Samaniego, and E. Samaniego.  A variational approach to the phase field modeling of brittle and ductile fracture.  _International Journal of Mechanical Sciences_ , 144:502–517, 2018.
  * Ambati et al. [2015b] M. Ambati, T. Gerasimov, and L. De Lorenzis.  Phase-field modeling of ductile fracture.  _Computational Mechanics_ , 55(5):1017–1040, 2015b.
  * Miehe et al. [2015] C. Miehe, M. Hofacker, L.-M. Schänzel, and F. Aldakheel.  Phase field modeling of fracture in multi-physics problems. Part II. Coupled brittle-to-ductile failure criteria and crack propagation in thermo-elastic–plastic solids.  _Computer Methods in Applied Mechanics and Engineering_ , 294:486–522, 2015.
  * Alessi et al. [2018a] R. Alessi, M. Ambati, T. Gerasimov, S. Vidoli, and L. De Lorenzis.  Comparison of phase-field models of fracture coupled with plasticity.  In _Advances in Computational Plasticity_ , pages 1–21. Springer, 2018a.
  * Noii et al. [2021] N. Noii, A. Khodadadian, J. Ulloa, F. Aldakheel, T. Wick, S. Francois, and P. Wriggers.  Bayesian inversion for unified ductile phase-field fracture.  _Computational Mechanics_ , 68:943–980, 2021.
  * You et al. [2020] T. You, Q.-Z. Zhu, P.-F. Li, and J.-F. Shao.  Incorporation of tension-compression asymmetry into plastic damage phase-field modeling of quasi brittle geomaterials.  _International Journal of Plasticity_ , 124:71–95, 2020\.
  * You et al. [2021] T. You, H. Waisman, and Q.-Z. Zhu.  Brittle-ductile failure transition in geomaterials modeled by a modified phase-field method with a varying damage-driving energy coefficient.  _International Journal of Plasticity_ , 136:102836, 2021\.
  * Kienle et al. [2019] D. Kienle, F. Aldakheel, and M.-A. Keip.  A finite-strain phase-field approach to ductile failure of frictional materials.  _International Journal of Solids and Structures_ , 172:147–162, 2019.
  * Späth et al. [2021] M. Späth, C. Herrmann, N. Prajapati, D. Schneider, F. Schwab, M. Selzer, and B. Nestler.  Multiphase-field modelling of crack propagation in geological materials and porous media with drucker-prager plasticity.  _Computational Geosciences_ , 25(1):325–343, 2021\.
  * Aldakheel [2020] F. Aldakheel.  A microscale model for concrete failure in poro-elasto-plastic media.  _Theoretical and Applied Fracture Mechanics_ , 107:102517, 2020.
  * Kienle and Keip [2021] D. Kienle and M.-A. Keip.  A variational minimization formulation for hydraulically induced fracturing in elastic-plastic solids.  _International Journal of Fracture_ , pages 1–25, 2021.
  * Bryant and Sun [2021] E.C. Bryant and W. Sun.  Phase field modeling of frictional slip with slip weakening/strengthening under non-isothermal conditions.  _Computer Methods in Applied Mechanics and Engineering_ , 375:113557, 2021.
  * Arthur et al. [1977] J.R.F. Arthur, T. Dunstan, Q.A.J.L. Al-Ani, and A. Assadi.  Plastic deformation and failure in granular media.  _Geotechnique_ , 27(1):53–74, 1977.
  * Vardoulakis [1980] I. Vardoulakis.  Shear band inclination and shear modulus of sand in biaxial tests.  _International Journal for Numerical and Analytical Methods in Geomechanics_ , 4(2):103–119, 1980.
  * Mielke [2006] A. Mielke.  A mathematical framework for generalized standard materials in the rate-independent case.  _Multifield Problems in Solid and Fluid Mechanics_ , 28:399–428, 2006.
  * Mielke and Roubícek [2015] A. Mielke and T. Roubícek.  _Rate-Independent systems. Theory and application_.  Springer, 2015.
  * Halphen and Nguyen [1975] B. Halphen and Q.S. Nguyen.  Generalized standard materials.  _Journal de Mécanique_ , 14(1):39–63, 1975\.
  * Biot [1965] M.A. Biot.  _Mechanics of incremental deformations_.  Wiley, 1965.
  * Mori and Tanaka [1973] T. Mori and K. Tanaka.  Average stress in matrix and average elastic energy of materials with misfitting inclusions.  _Acta Metallurgica_ , 21(5):571–574, 1973.
  * Alessi et al. [2018b] R. Alessi, J.-J. Marigo, C. Maurini, and S. Vidoli.  Coupling damage and plasticity for a phase-field regularisation of brittle, cohesive and ductile fracture: one-dimensional examples.  _International Journal of Mechanical Sciences_ , 149:559–576, 2018b.
  * Samaniego et al. [2021] C. Samaniego, J. Ulloa, P. Rodríguez, G. Houzeaux, M. Vázquez, and E. Samaniego.  A phase-field model for ductile fracture with shear bands: A parallel implementation.  _International Journal of Mechanical Sciences_ , 200:106424, 2021.
  * Kuhn et al. [2015] C. Kuhn, A. Schlüter, and R. Müller.  On degradation functions in phase field fracture models.  _Computational Materials Science_ , 108:374–384, 2015.
  * Wu et al. [2020] J.-Y. Wu, V.P. Nguyen, C.T. Nguyen, D. Sutula, S. Bordas, and S. Sinaie.  Phase field modeling of fracture.  _Advances in Applied Mechanics_ , 53:1–183, 2020.
  * Hill [1948] R. Hill.  A variational principle of maximum plastic work in classical plasticity.  _The Quarterly Journal of Mechanics and Applied Mathematics_ , 1(1):18–28, 1948.
  * Moreau [1974] J.-J. Moreau.  On unilateral constraints, friction and plasticity.  In _New variational techniques in mathematical physics_ , pages 171–322. CIME, 1974.
  * Rockafellar [1970] R.T. Rockafellar.  _Convex analysis_ , volume 28.  Princeton University Press, 1970.
  * Han and Reddy [1999] W. Han and B.D. Reddy.  _Plasticity: mathematical theory and numerical analysis_ , volume 9.  Springer Science & Business Media, 1999.
  * Francfort [2018] G.A. Francfort.  Recovering convexity in non-associated plasticity.  _Comptes Rendus Mécanique_ , 346(3):198–205, 2018.
  * Laborde [1987] P. Laborde.  Analysis of the strain-stress relation in plasticity with non-associated laws.  _International Journal of Engineering Science_ , 25(6):655–666, 1987.
  * Bleyer and Alessi [2018] J. Bleyer and R. Alessi.  Phase-field modeling of anisotropic brittle fracture including several damage mechanisms.  _Computer Methods in Applied Mechanics and Engineering_ , 336:213–236, 2018.
  * Alessi [2016] R. Alessi.  Energetic formulation for rate-independent processes: remarks on discontinuous evolutions with a simple example.  _Acta Mechanica_ , 227(10):2805–2829, 2016.
  * Alessi and Bernardini [2015] R. Alessi and D. Bernardini.  Analysis of localization phenomena in shape memory alloys bars by a variational approach.  _International Journal of Solids and Structures_ , 73:113–133, 2015.
  * Rokoš et al. [2016] O. Rokoš, J. Zeman, and M. Jirásek.  Localization analysis of an energy-based fourth-order gradient plasticity model.  _European Journal of Mechanics-A/Solids_ , 55:256–277, 2016\.
  * Alessi et al. [2018c] R. Alessi, S. Vidoli, and L. De Lorenzis.  A phenomenological approach to fatigue with a variational phase-field model: The one-dimensional case.  _Engineering Fracture Mechanics_ , 190:53–73, 2018c.
  * Luege et al. [2018] M. Luege, A. Orlando, M.E. Almenar, and E.A. Pilotta.  An energetic formulation of a gradient damage model for concrete and its numerical implementation.  _International Journal of Solids and Structures_ , 155:160–184, 2018.
  * Lancioni and Alessi [2020] G. Lancioni and R. Alessi.  Modeling micro-cracking and failure in short fiber-reinforced composites.  _Journal of the Mechanics and Physics of Solids_ , 137:103854, 2020.
  * Ulloa et al. [2021b] J. Ulloa, J. Wambacq, R. Alessi, G. Degrande, and S. François.  Phase-field modeling of fatigue coupled to cyclic plasticity in an energetic formulation.  _Computer Methods in Applied Mechanics and Engineering_ , 373:113473, 2021b.
  * Simo and Hughes [1998] J.C. Simo and T.J.R. Hughes.  _Computational inelasticity_.  Springer, 1998.
  * de Souza Neto et al. [2011] E.A. de Souza Neto, D. Peric, and D.R.J. Owen.  _Computational methods for plasticity: theory and applications_.  John Wiley & Sons, 2011.
  * Borja [2013] R.I. Borja.  _Plasticity: modeling & computation_.  Springer Science & Business Media, 2013.
  * Sysala et al. [2016] S. Sysala, M. Cermak, T. Koudelka, J. Kruis, J. Zeman, and R. Blaheta.  Subdifferential-based implicit return-mapping operators in computational plasticity.  _ZAMM-Journal of Applied Mathematics and Mechanics/Zeitschrift für Angewandte Mathematik und Mechanik_ , 96(11):1318–1338, 2016.
  * Miehe [2011] C. Miehe.  A multi-field incremental variational framework for gradient-extended standard dissipative solids.  _Journal of the Mechanics and Physics of Solids_ , 59(4):898–923, 2011.
  * Wambacq et al. [2021] J. Wambacq, J. Ulloa, G. Lombaert, and S. François.  A dissipation-based path-following technique for the phase-field approach to brittle and ductile fracture.  _International Journal for Numerical Methods in Engineering_ , 2021\.
  * Ord et al. [1991] A. Ord, I. Vardoulakis, and R. Kajewski.  Shear band formation in gosford sandstone.  In _International Journal of Rock Mechanics and Mining Sciences & Geomechanics Abstracts_, volume 28, pages 397–409. Elsevier, 1991.
  * Labuz et al. [1996] J.F. Labuz, S.-T. Dai, and E. Papamichos.  Plane-strain compression of rock-like materials.  In _International Journal of Rock Mechanics and Mining Sciences & Geomechanics Abstracts_, volume 33, pages 573–584. Elsevier, 1996.
  * Fakhimi et al. [2002] A. Fakhimi, F. Carvalho, T. Ishida, and J.F. Labuz.  Simulation of failure around a circular opening in rock.  _International Journal of Rock Mechanics and Mining Sciences_ , 39(4):507–515, 2002.
  * Labuz et al. [2006] J.F. Labuz, J.J. Riedel, and S.-T. Dai.  Shear fracture in sandstone under plane-strain compression.  _Engineering Fracture Mechanics_ , 73(6):820–828, 2006.
  * Chang et al. [2002] S.-H. Chang, C.-I. Lee, and S. Jeon.  Measurement of rock fracture toughness under modes I and II and mixed-mode conditions by using disc-type specimens.  _Engineering Geology_ , 66(1-2):79–97, 2002.
  * Ayatollahi and Aliha [2007] M.R. Ayatollahi and M.R.M Aliha.  Wide range data for crack tip parameters in two disc-type specimens under mixed mode loading.  _Computational Materials Science_ , 38(4):660–670, 2007.
  * Ayatollahi and Aliha [2008] M.R. Ayatollahi and M.R.M. Aliha.  On the use of brazilian disc specimen for calculating mixed mode I–II fracture toughness of rock materials.  _Engineering Fracture Mechanics_ , 75(16):4631–4641, 2008.
  * Aliha et al. [2010] M.R.M. Aliha, M.R. Ayatollahi, D.J. Smith, and M.J. Pavier.  Geometry and size effects on fracture trajectory in a limestone rock under mixed mode loading.  _Engineering Fracture Mechanics_ , 77(11):2200–2212, 2010.
  * Haeri et al. [2014] H. Haeri, K. Shahriar, M.F. Marji, and P. Moarefvand.  Experimental and numerical study of crack propagation and coalescence in pre-cracked rock-like rock-like disks.  _International Journal of Rock Mechanics and Mining Sciences_ , 67:20–28, 2014.
  * Xiankai et al. [2018] B. Xiankai, T. Meng, and Z. Jinchang.  Study of mixed mode fracture toughness and fracture trajectories in gypsum interlayers in corrosive environment.  _Royal Society Open Science_ , 5(1):171374, 2018\.
  * Wang and Zhou [2020] L. Wang and X. Zhou.  Phase field model for simulating the fracture behaviors of some disc-type specimens.  _Engineering Fracture Mechanics_ , 226:106870, 2020.
  * Erdogan and Sih [1963] F. Erdogan and G.C. Sih.  On the crack extension in plates under plane loading and transverse shear.  _Journal of Basic Engineering_ , 85:519–525, 1963.
  * Smith et al. [2001] D.J. Smith, M.R. Ayatollahi, and M.J. Pavier.  The role of t-stress in brittle fracture for linear elastic materials under mixed-mode loading.  _Fatigue & Fracture of Engineering Materials & Structures_, 24(2):137–150, 2001.
  * Navidtehrani et al. [2021] Y. Navidtehrani, C. Betegón, and E. Martínez-Pañeda.  A simple and robust abaqus implementation of the phase field fracture method.  _Applications in Engineering Science_ , 6:100050, 2021.
  * Ramsey and Chester [2004] J.M. Ramsey and F.M. Chester.  Hybrid fracture and the transition from extension fracture to shear fracture.  _Nature_ , 428(6978):63–66, 2004.
  * Zhu and Shao [2015] Q.-Z. Zhu and J.-F. Shao.  A refined micromechanical damage–friction model with strength prediction for rock-like materials under compression.  _International Journal of Solids and Structures_ , 60:75–83, 2015.
  * Schofield and Wroth [1968] A. Schofield and P. Wroth.  _Critical state soil mechanics_.  McGraw-hill, 1968.
  * Salvati [2021] E. Salvati.  Residual stress as a fracture toughening mechanism: A phase-field study on a brittle material.  _Theoretical and Applied Fracture Mechanics_ , 114:103021, 2021.