Mesoscale Modelling of Concretes Subjected to Triaxial Loadings: Mechanical Properties and Fracture Behaviour

The mechanical properties and fracture behaviour of concretes under different triaxial stress states were investigated based on a 3D mesoscale model. The quasistatic triaxial loadings, namely, compression–compression–compression (C–C–C), compression–tension–tension (C–T–T) and compression–compression–tension (C–C–T), were simulated using an implicit solver. The mesoscopic modelling with good robustness gave reliable and detailed damage evolution processes under different triaxial stress states. The lateral tensile stress significantly influenced the multiaxial mechanical behaviour of the concretes, accelerating the concrete failure. With low lateral pressures or tensile stress, axial cleavage was the main failure mode of the specimens. Furthermore, the concretes presented shear failures under medium lateral pressures. The concretes experienced a transition from brittle fracture to plastic failure under high lateral pressures. The Ottosen parameters were modified by the gradient descent method and then the failure criterion of the concretes in the principal stress space was given. The failure criterion could describe the strength characteristics of concrete materials well by being fitted with experimental data under different triaxial stress states.


Introduction
In practical engineering structures, concretes are inevitably subject to biaxial and triaxial stress states [1]. In particular, for some massive concrete structures, such as water gates, dams, concrete offshore platforms, foundations of heavy machinery and nuclear reactors, concretes are in a complex stress state. The three principal stresses of concretes are unequal, and most of them are even different combinations of compressive stresses and tensile stresses. In order to accurately analyse and design these structures, it is necessary to study the mechanical responses and damage evolution characteristics of concretes under multiaxial compressive and tensile stress conditions [2]. Many experimental studies of the mechanical behaviour of concretes under multiaxial triaxial loadings have been carried out by researchers [3][4][5][6]. Su et al. studied the multiaxial mechanical behaviour of foamed concrete under triaxial confining pressures [7]. Javanmardi [8] calculated the strength relationship between the lateral confining stress and uniaxial stress by summarising the experimental data of extensive and systematic multiaxial tests and compared these with the proposed relationships by other researchers [9][10][11][12][13][14][15][16]. Glatz et al. developed a novel, high-temperature and high-pressure triaxial system that is suitable for X-ray computed tomography (CT) scanning to capture the role of effective stress on the fluid distribution of rocks under high triaxial pressures [17]. Watanabe et al. designed a new true triaxial cell to study the initiation and propagation of microfractures in granite under confining pressures [18]. Shang analysed the experimental results of normal air-entrained concretes under a tension-compression state and compared the mechanical behaviour of different types of concretes under the same loading conditions [19]. Wang and Song investigated the mechanical characteristics of mass concrete under a triaxial compression-tension-tension loading condition [1]. The specimens were tested in a specially designed, truly triaxial apparatus (designed by the State Key Laboratory of Coastal and Offshore Engineering), which was capable of developing three independent compressive or tensile forces. It can be noted that most of the experiments were done under triaxial confining pressures, but a few triaxial loading experiments were done under complex tensile and compressive stress states. Because the development of multiaxial tensile test system technology is difficult and has a high cost, and several key technical problems exist in the experiments, to date, only a few works have been dedicated to investigating the mechanical properties of concretes under multiaxial tensile stress states. According to the test results obtained by different researchers, it was shown that the mechanical behaviour of concretes under triaxial loads was affected by many factors, such as the loading rate, stiffness of the testing apparatus and the measuring and testing techniques, such as the adopted loading mode and the friction-reducing means along the loading direction. However, despite the difficult test system technology and the lack of uniform specifications, the mesoscopic modelling that demonstrated good robustness gives reliable and detailed damage evolution processes under complex multiaxial loadings. This is helpful towards further understanding the triaxial mechanisms of damage evolution and the mechanical behaviour of concretes.
As an acceptable, cost-effective and time-efficient alternative for experiments, the mesoscale concrete model with accurate constitutive behaviour provides an effective way to predict mechanical behaviour and understand the damage mechanisms of concretes under multiaxial stress conditions. In particular, with the developments in computer capabilities, mesoscale modelling should be a powerful tool for investigating the mechanical behaviour of heterogeneous materials, such as concretes, at the mesoscale level. Comby-Peyrot et al. [20] also suggested that modelling the heterogeneities of concretes was important for studying the damage initiation and progression and investigating the local deformation mechanisms in concrete. Karavelić et al. [21] used the Gaussian distribution to model the heterogeneity of material properties for the aggregate and cement phases and showed that the numerical model could simulate the salient characteristics of concrete-like responses well. With the application of mesoscale modelling in the study of the concrete structure [15,22], the behaviour of concrete material under multiaxial loadings can be well understood. Many researchers have studied the mechanical properties of concretes under uniaxial and biaxial loadings [21,23,24]. Caballero et al. [25] described a method for the 3D mesostructural analysis of concrete specimens under biaxial loading. The results indicated that both the cracking and failure modes and the envelope were in good agreement with the experimental results. However, due to the imperfections for some mesoscopic models and the lack of uniform specifications for the sample size and experimental apparatus, relatively few numerical simulations of triaxial loadings for concretes have been investigated. Currently, the triaxial compressive conditions are mainly studied, especially for high confining pressures. Tran et al. [26] studied the behaviour of concretes under high-confining pressures through a three-dimensional numerical model. The internal physical mechanism of the random cracks and mechanical behaviour of concretes could be well reflected by the mesoscale modelling. Dupray et al. [27] performed numerical analyses of the experiments of concretes under quasistatic triaxial compression at the mesoscopic scale. The numerical model provided a reproduction of the main characteristics of mechanical behaviour for concretes under high confining pressures, both qualitatively and quantitatively: mortar compaction, damage localisation and the influence of aggregates on the limit state. The results suggested that the numerical modelling was helpful for understanding the mesoscale phenomena and the influence of damage zones on macroscale behaviour.
As mentioned above, the study of the mechanical responses and damage evolution characteristics of concretes under multiaxial compressive and tensile stresses is of great significance for concrete structures in terms of rational design, safe operation and ma-terial savings. However, due to the heterogeneity of meso-concrete and the complexity of its internal structure, there are rare systematic and comprehensive numerical simulation studies of concretes under triaxial loadings, especially under the combination of unequal stresses of tension and compression. This study investigated the triaxial mechanical behaviour of 3D mesoscale concretes under various loading conditions, namely, compression-compression-compression (C-C-C), compression-tension-tension (C-T-T) and compression-compression-tension (C-C-T). Meanwhile, the damage evolution within the entire concrete specimens was further discussed. Moreover, rich simulation data were systematically obtained to identify the mechanical responses under different triaxial stress states. In this study, the gradient descent method was adopted to iteratively modify the Ottosen parameters to further modify the failure criteria in the stress space of the concretes. The good fitting between the experimental data and the failure surface verified the validity of triaxial numerical simulations in this study. On one hand, it was shown that the mechanical performance of concrete structures under various complex loadings could be predicted and reasonably checked using the modified Ottosen failure criteria. On the other hand, these simulations provided a foundation for establishing the quantitative relations between the macro-mechanical properties of concretes and the mesoscopic parameters of concrete components.

Model Generation
The convex polyhedral aggregates inside concrete specimens were generated through Voronoi diagrams. Each Voronoi cell V i consists of each point P in the region. Based on the 3D Voronoi model, a random scaling factor is introduced to generate the random convex aggregates via a shrinking process, and the ITZ (interfacial transition zone) is generated through an extension process. The mesoscale concrete models are efficiently generated.
To avoid malformed cells, the minimum distance between any two nuclei should be controlled based on Equation (1): where δ 0 is the distance between two nuclei and the value of K = 0.2 [28] was adopted to achieve complete randomness regarding the shape of the aggregates. MATLAB (Version R2014a, MathWorks, Natick, MA, USA) [29] was used to establish the initial model of aggregates according to the 3D Voronoi diagram, which is plotted in Figure 1. A random shrinking process was performed to separate adjacent cells. For each cell, the vector → v = (x, y, z) from a vertex p i to its corresponding seed S i can be given  and shrunk by the factor q (0 < q < 1). The updated vertex p i is calculated using Equation (2): The shrinking scheme of a single polyhedral cell is shown in Figure 2. In order to achieve the grading requirements, various shrinking factors are employed according to the Fuller curve [30], which is expressed as Equation (3): where B(d) is the total percentage of aggregates that pass the sieve with an aperture diameter d, d max is the maximum diameter of aggregates and n represents an exponent where d max and d min represent the maximum and minimum diameters of the aggregates, V agg [d i , d i+1 ] is the aggregates' volume within the aggregates' size range [d i , d i+1 ] and V represents the specimen volume. For the grading range [d i , d i+1 ], the random shrinking factor q can be given as: where ω represents a random variable ranging from 0 to 1. The shrinking scheme of a single polyhedral cell is shown in Figure 2. In order to achieve the grading requirements, various shrinking factors are employed according to the Fuller curve [30], which is expressed as Equation (3): where B(d) is the total percentage of aggregates that pass the sieve with an aperture diameter d, max d is the maximum diameter of aggregates and n represents an exponen ranging from 0.45 to 0.7, but is generally considered to be 0.5. The aggregate volume within each grading range [d i , d i 1 ] can be described as:    Figure 3 shows the modelling results of the random aggregates with the grading sizes. Figure 4 gives the diagram of the ITZ-generating process. Compared with coarse aggregates, the characteristic length of ITZ was very small. Taking into account both the accuracy and computational efficiency, an appropriate thickness of 0.1 mm was employed through an extending process [31,32] and the "free-falling" procedure was adopted to improve the aggregate content [33].

Damaged Plastic Constitutive Model
The heterogeneous constitution, random distribution of phases and a microcrack patterns all contribute to the complex mechanical behaviour of concre The damage plasticity model was used to characterise the irreversible deformations and the stiffness degradation behaviour of the concretes. The mo been well validated for the damage evolution analysis of mesoscopic concretes Under multiaxial loadings, the effective plastic strain rates of tension and compres calculated as follows:

Damaged Plastic Constitutive Model
The heterogeneous constitution, random distribution of phases microcrack patterns all contribute to the complex mechanical behaviour of The damage plasticity model was used to characterise the irreve deformations and the stiffness degradation behaviour of the concretes. T been well validated for the damage evolution analysis of mesoscopic con Under multiaxial loadings, the effective plastic strain rates of tension and co calculated as follows:

Damaged Plastic Constitutive Model
The heterogeneous constitution, random distribution of phases and arbitrary microcrack patterns all contribute to the complex mechanical behaviour of concretes [34]. The damage plasticity model was used to characterise the irreversible plastic deformations and the stiffness degradation behaviour of the concretes. The model has been well validated for the damage evolution analysis of mesoscopic concretes [35][36][37][38]. Under multiaxial loadings, the effective plastic strain rates of tension and compression are calculated as follows: . ε pl t = r σ ˆ. ε pl max (6) .
whereˆ. ε pl max andˆ. ε pl min are the maximum and minimum eigenvalues of the plastic strain rate tensor . ε pl ; r σ represents a stress weight factor.
The yield condition proposed by Lee and Fenves [39,40] was used in the plastic damage concrete model to account for the different strength evolutions under different loadings. For effective stresses, the form of the yield function is: where α and γ are dimensionless material constants, p is the effective hydrostatic pressure, q is the Mises equivalent effective stress, andσ max is the algebraical maximum eigenvalue of the effective stress tensor. Furthermore: where σ t and σ c are the effective tensile and compressive cohesion stresses, respectively. The coefficient α can be calculated using the initial equibiaxial yield stress σ b0 and the uniaxial compressive yield stress σ c0 . Whenσ max < 0, the coefficient γ enters the yield function only for the triaxial compressive stress states. γ can be determined by comparing the yield conditions along the tensile meridians (TMs) and compressive meridians (CMs). The corresponding yield conditions are given as Equations (12) and (13). As a typical concrete value, a value of γ = 3 gives K c = 2 3 .
Typical yield surfaces are shown in Figure 5 in the deviatoric plane and in Figure 6 for plane stress conditions. For the damaged plastic constitutive model, the flow potential G was chosen to be the Drucker-Prager hyperbolic function: where ψ is the dilation angle that is measured on the failure curve in the p-q plane at high confining pressures, ω represents the eccentricity and σ t0 is the uniaxial tensile stress at failure.

Mechanical Behaviour and Failure Analysis of Concretes under Triaxial Loadings
In this section, the results of the triaxial mechanical behaviour investigation o mesoscopic concretes under various loadings, such as C-C-C, C-T-T and C-C-T, are presented. Then, the damage evolution within the entire concrete specimens is further discussed to clarify the damage mechanism of concretes under triaxial loadings. In this study, the specimen considered was a 25 mm side length cube. The aggregates and mortar phases adopted the mesh size set to 1.0 mm as an average size using tetrahedral elements For the wedge element used in the ITZ, the element size depended on the contiguous aggregate element, while the size in the thickness direction was an invariable 0.1 mm. An implicit solver was employed to investigate the triaxial loadings using the nonlinear finite element (FE) code ABAQUS/standard (Version 6.14, Dassault Systèmes, Paris, France

Mechanical Behaviour and Failure Analysis of Concretes under Triaxial Loa
In this section, the results of the triaxial mechanical behaviour investig mesoscopic concretes under various loadings, such as C-C-C, C-T-T and C-C presented. Then, the damage evolution within the entire concrete specimens is discussed to clarify the damage mechanism of concretes under triaxial loadings study, the specimen considered was a 25 mm side length cube. The aggregates and phases adopted the mesh size set to 1.0 mm as an average size using tetrahedral e For the wedge element used in the ITZ, the element size depended on the con aggregate element, while the size in the thickness direction was an invariable 0.1 Figure 6. Illustration of the yield surface due to the plane stress, as defined by the concrete damaged plasticity (CDP) model.

Mechanical Behaviour and Failure Analysis of Concretes under Triaxial Loadings
In this section, the results of the triaxial mechanical behaviour investigation of mesoscopic concretes under various loadings, such as C-C-C, C-T-T and C-C-T, are presented. Then, the damage evolution within the entire concrete specimens is further discussed to clarify the damage mechanism of concretes under triaxial loadings. In this study, the specimen considered was a 25 mm side length cube. The aggregates and mortar phases adopted the mesh size set to 1.0 mm as an average size using tetrahedral elements. For the wedge element used in the ITZ, the element size depended on the contiguous aggregate element, while the size in the thickness direction was an invariable 0.1 mm. An implicit solver was employed to investigate the triaxial loadings using the nonlinear finite element (FE) code ABAQUS/standard (Version 6.14, Dassault Systèmes, Paris, France) [41], where a small time increment had to be used to reduce the influence of the loading kinematics as far as possible. This model is able to precisely control the components, including the geometrical characteristics and the mechanical parameters, and has the advantage of stable convergence by using the implicit algorithm during the simulational processes. The volume fraction of aggregates in the cubic concrete specimens was 37.38%. The aggregate size segment is shown in Table 1. The main computational parameters of the concrete mesoscale contents adopted in the current simulations are given in Table 2. Considering the mortar matrix with a Young's modulus of 25 GPa, the standard strength of mortar was assumed to be 35 MPa [42]. Due to the limitations of existing measuring techniques, it is hard to clarify the physical and mechanical properties of the ITZ with significant inherent defects in its internal structure. In this study, the ITZ strength was calculated based on the reduction in mortar strength, which can be expressed as the ratio of the ITZ strength to the mortar strength. The ratios adopted by different researchers range from 0.5 to 0.9 [34,43]. One study found that it is reasonable for the equivalent ITZ strength to be 70% of the mortar strength [44]. For the current model of the concretes, the compressive strength of the equivalent ITZ was around 20 MPa, with the Young's modulus set to 18 GPa. The aggregates' properties could depend significantly on the natural types. For crushed stones, the Young's modulus of the aggregates is between 40 and 60 GPa [45]. In the concrete damaged plasticity (CDP) model, the eccentricity was determined using the default value for concretes provided by ABAQUS [41], namely, 0.1. The viscosity coefficient was 1.0 × 10 −5 [46].

Model Validation under Triaxial Confining Pressures
Under triaxial stress states, the overall strength of concretes of the same specimen size is improved with an increased uniaxial compressive strength f c . Therefore, in order to eliminate the influence of different values of f c on the macrostrength, the stress ratios of the axial stress of σ divided by the uniaxial compressive strength of f c can be used to investigate the triaxial mechanical behaviour of concretes.
For numerical simulations of concretes under triaxial loadings, the common loading method is shown in Figure 7, which represents the skeleton diagram of concrete specimens under triaxial loadings. The compressive stresses or tensile stresses are applied in the horizontal direction according to the stress ratio σ 1 /σ 3 (along the directions of the X-axis (σ 1 ) and Y-axis (σ 3 )), and the displacement is applied along the vertical direction (the direction of the Z-axis (σ 2 )). The mesh convergence study was conducted before the loading analysis [32]. An element length of 1 mm was employed in the following analysis. With equal confining stresses (σ1/σ3 = 1) along the horizontal direction, Figure 8 illustrates the stress-strain curves of experimental tests and numerical simulations under different compressive stresses. Experimental data were obtained from [15], where the confining pressure ranged from 0 to 30 MPa. The confining pressure was applied to the concrete specimens by changing the value of the hydrostatic pressure. The stress-strain curves changed with the increasing confining pressure. With low confining pressures, the stress-strain curves had clear peak stress points and then experienced the falling stage. Beyond the peak stresses, the softening stages increased with the increasing confining pressures and the stress during the softening stage tended to be stable. Under the confining pressure of 9 MPa, there were no significant descending stages in the curves. When the confining pressure continued to increase to 30 MPa, which was close to the uniaxial compressive strength, the maximum stress of the softening stage could exceed the peak stress, as shown in Figure 8. At that time, the curve went up continuously as the strain increased, but the growth rate of the stress, which is represented by the slope of the curve, decreased gradually and approached 0. Under low confining pressures, the concretes exhibited a brittle fracture feature. With sufficiently large confining pressures, the brittle behaviour of the concretes would disappear, preventing With equal confining stresses (σ 1 /σ 3 = 1) along the horizontal direction, Figure 8 illustrates the stress-strain curves of experimental tests and numerical simulations under different compressive stresses. Experimental data were obtained from [15], where the confining pressure ranged from 0 to 30 MPa. The confining pressure was applied to the concrete specimens by changing the value of the hydrostatic pressure. With equal confining stresses (σ1/σ3 = 1) along the horizontal direction, F illustrates the stress-strain curves of experimental tests and numerical simulations different compressive stresses. Experimental data were obtained from [15], wh confining pressure ranged from 0 to 30 MPa. The confining pressure was applied concrete specimens by changing the value of the hydrostatic pressure. The stress-strain curves changed with the increasing confining pressure. W confining pressures, the stress-strain curves had clear peak stress points an experienced the falling stage. Beyond the peak stresses, the softening stages inc with the increasing confining pressures and the stress during the softening stage to be stable. Under the confining pressure of 9 MPa, there were no significant desc stages in the curves. When the confining pressure continued to increase to 30 MPa was close to the uniaxial compressive strength, the maximum stress of the softenin could exceed the peak stress, as shown in Figure 8. At that time, the curve w continuously as the strain increased, but the growth rate of the stress, which is repre by the slope of the curve, decreased gradually and approached 0. Under low co pressures, the concretes exhibited a brittle fracture feature. With sufficiently confining pressures, the brittle behaviour of the concretes would disappear, prev The stress-strain curves changed with the increasing confining pressure. With low confining pressures, the stress-strain curves had clear peak stress points and then experienced the falling stage. Beyond the peak stresses, the softening stages increased with the increasing confining pressures and the stress during the softening stage tended to be stable. Under the confining pressure of 9 MPa, there were no significant descending stages in the curves. When the confining pressure continued to increase to 30 MPa, which was close to the uniaxial compressive strength, the maximum stress of the softening stage could exceed the peak stress, as shown in Figure 8. At that time, the curve went up continuously as the strain increased, but the growth rate of the stress, which is represented by the slope of the curve, decreased gradually and approached 0. Under low confining pressures, the concretes exhibited a brittle fracture feature. With sufficiently large confining pressures, the brittle behaviour of the concretes would disappear, preventing cracks from propagating. In these conditions, the failure of the concretes was driven by the consolidation of the microcracks and the collapse of the micropores, leading to a mechanical response of the brittle materials at the macrolevel. Correspondingly, the numerical simulation results show that the trends of the curves agreed well with the experiments below the confining pressure of 9 MPa. The deviation between the simulation and experimental results fluctuated within a reasonable range at the points of peak stress with increasing confining stress. However, the rising stage of the numerical simulation was slightly different from the experiment under the confining pressure of 30 MPa, which can be explained as follows. Concrete materials are extremely sensitive to confining pressures, which have a great influence on the yield strength of concretes [15]. The increasing confining pressure leads to a change in the principal failure mechanism and substantial flaws or microcracks are inevitably generated inside the initial concretes [47]. The collapse of pores and closure of microcracks during the compacting process result in significantly large strains, up to 3% [5]. In triaxial compressive tests, the mechanical response of concretes presented a clear peak, followed by a relatively smooth decline curve under low confining pressures. Then, with the monotonic increase in confinement, the response was gradually transformed to be ductile, resulting in a plateau. In this study, as a highly nonhomogeneous artificial composite material, the mesoscopic model of the concretes comprised aggregate, mortar and ITZs. The CDP model based on an isotropic damage hypothesis was employed as the constitutive laws of mortar and the ITZ to investigate the inelastic behaviour of concretes. Up to the confining pressure of 30 MPa, the axial strain corresponding to the maximum stress was less than 1%, which was slightly less than that in experiments. The tendencies of the curves and the maximum stresses were close to the experimental results. Consequently, the numerical simulations agreed well with the experiments, which made clear that the mesoscopic model in this study can be used for the mechanical analysis of concretes under triaxial complex stress conditions.

Compression-Compression-Compression Loading
Based on the uniaxial numerical results, the uniaxial compressive strength of the concretes' f c was 32.3 MPa and the uniaxial tensile strength f t was 2.8 MPa. For the sake of description, the stress ratio of each axial stress to the uniaxial strength, namely, σ 1 /f c , σ 2 /f c and σ 3 /f c , was adopted to represent the triaxial numerical simulation results. In order to ensure that the specimens were not damaged prematurely, the imposed stresses should not be higher than the uniaxial compressive and tensile strengths. The compressive stress ratios of σ 1 /σ 3 under various C-C-C loadings are shown in Table 3. There were 25 kinds of compressive loading combinations, where the compressive stresses along the X-axis and Y-axis varied from 0.083f c to f c , while the stress ratios along the X and Y directions varied from 0.083 to 12. Table 3. Compressive stress ratios of the X-axis to the Y-axis directions. The stress-strain curves under different compressive stresses of σ 1 when the compressive stress of σ 3 (Y) was 0.25f c are shown in Figure 9. Compared to the uniaxial compression, the whole strength of the concretes under triaxial confining pressures increased to 2f c -3f c and then experienced a small drop after reaching the peak. Keeping the compressive stress along the Y-axis constant, the peak stress increased with the increasing compressive stress along the X-axis. Furthermore, with the increasing stress of σ 1 (X), the softening descending stage was increasingly shorter and the stress platform stage became higher. Because of the strong constraint, the stress was even strengthened with the increasing strain. The tendency of the curves is shown in Figure 9, in which the peak stress was not necessarily the maximum of the curves. After the softening stage, the stress increased to a certain extent and even exceeded the peak stress, as shown at point "Q" in Figure 9. However, the internal structure of the concretes may have begun to be damaged during the process of the stresses increasing after the peak stress for real concretes. For safety reasons, the peak stress was still considered as its compressive strength in the subsequent analysis. After the compressive stress along the X-direction reached 0.5f c , the compressive strength was no longer significantly increased.    Figure 10 shows the processes of the damage evolution and failure modes of the concrete specimens under different compressive σ1/σ3 stress ratios while keeping the stress along the Y-axis constant. It can be seen that no obvious damage zones were formed at the peak stress points. After the peak stress, the damage evolution processes along the Y-axis were basically identical. The damage initiated at the edges of specimen surfaces perpendicular to the Y-axis, then developed obliquely in the Y plane and finally formed intersecting cracks. The damage pattern of the X-Y plane was closely related to the change in the compressive stress ratios. When σ1/σ3 was less than 1, the generated single cracks were perpendicular to the Z-direction. Single cracks formed into grid cracks as σ1/σ3 became equal to one and then grid cracks turned into single cracks after σ1/σ3 was greater than 1. It can be seen that with the increasing stress ratio, the failure mode changed from single cracks to crossing cracks and then to single cracks again.  Figure 10 shows the processes of the damage evolution and failure modes of the concrete specimens under different compressive σ 1 /σ 3 stress ratios while keeping the stress along the Y-axis constant. It can be seen that no obvious damage zones were formed at the peak stress points. After the peak stress, the damage evolution processes along the Y-axis were basically identical. The damage initiated at the edges of specimen surfaces perpendicular to the Y-axis, then developed obliquely in the Y plane and finally formed intersecting cracks. The damage pattern of the X-Y plane was closely related to the change in the compressive stress ratios. When σ 1 /σ 3 was less than 1, the generated single cracks were perpendicular to the Z-direction. Single cracks formed into grid cracks as σ 1 /σ 3 became equal to one and then grid cracks turned into single cracks after σ 1 /σ 3 was greater than 1. It can be seen that with the increasing stress ratio, the failure mode changed from single cracks to crossing cracks and then to single cracks again.
The damage modes of the specimens under different compressive stresses with a strain of 0.01 when the stress ratio σ 1 /σ 3 was 1 are shown in Figure 11. As can be seen, the damage patterns were basically identical and oblique latticed damage zones were formed on the three loading surfaces along the X-, Y-and Z-directions. The main difference between the concretes under different compressive stresses was that the weak lateral restraint was conducive to the damage development with low confining pressure. Basically, the cracks propagated through the whole concrete specimens. With the increasing confining pressure, the lateral restraint was strengthened and the process of damage evolution slowed down. The damage modes of the specimens under different compressive stresses with a strain of 0.01 when the stress ratio σ1/σ3 was 1 are shown in Figure 11. As can be seen, the damage patterns were basically identical and oblique latticed damage zones were formed on the three loading surfaces along the X-, Y-and Z-directions. The main difference between the concretes under different compressive stresses was that the weak lateral restraint was conducive to the damage development with low confining pressure. Basically, the cracks propagated through the whole concrete specimens. With the increasing confining pressure, the lateral restraint was strengthened and the process of damage evolution slowed down.  The damage modes of the specimens under different compressive stresses with a strain of 0.01 when the stress ratio σ1/σ3 was 1 are shown in Figure 11. As can be seen, the damage patterns were basically identical and oblique latticed damage zones were formed on the three loading surfaces along the X-, Y-and Z-directions. The main difference between the concretes under different compressive stresses was that the weak lateral restraint was conducive to the damage development with low confining pressure. Basically, the cracks propagated through the whole concrete specimens. With the increasing confining pressure, the lateral restraint was strengthened and the process of damage evolution slowed down.  The results of the 25 groups of numerical simulations are shown in Figure 12, representing the strength of the concretes under triaxial compression. The horizontal axis in Figure 12a is the absolute value of the stress ratio σ 1 /f c and the horizontal axis in Figure 12b is the absolute value of the stress ratio σ 1 /σ 3 . It can be noted that the compressive stress of σ 1 along the X-axis had an insignificant influence on the growth of concrete compressive strength with low compressive σ 3 (Y) stresses. However, with the increasing compressive σ 3 (Y) stress, the overall strength of the concrete specimens was improved significantly.
From the curves, it can be seen that the effects of the lateral compressive stresses on the overall strength could be categorised into one of three situations. The strength of the concretes tended to be relatively constant as the lateral compressive stress was less than 10% of the uniaxial strength. On the other hand, when the lateral compressive stress increased to 50% of the uniaxial strength, the strength of the concretes reached the stress maximum. The strength of the concretes tended to increase linearly when the lateral compressive stress increased to the uniaxial strength, as shown in Figure 12b. in Figure 12a is the absolute value of the stress ratio σ1/fc and the horizontal axis in Figu 12b is the absolute value of the stress ratio σ1/σ3. It can be noted that the compressive str of σ1 along the X-axis had an insignificant influence on the growth of concrete compress strength with low compressive σ3(Y) stresses. However, with the increasing compress σ3(Y) stress, the overall strength of the concrete specimens was improved significant From the curves, it can be seen that the effects of the lateral compressive stresses on overall strength could be categorised into one of three situations. The strength of concretes tended to be relatively constant as the lateral compressive stress was less th 10% of the uniaxial strength. On the other hand, when the lateral compressive str increased to 50% of the uniaxial strength, the strength of the concretes reached the str maximum. The strength of the concretes tended to increase linearly when the late compressive stress increased to the uniaxial strength, as shown in Figure 12b.

Compression-Tension-Tension Loading
Compared to the triaxial compressive loading, the loading modes of C-T-T involv changing the stresses along the X-and Y-axes, with constant stress along the Z-axis. Ta 4 shows various stress of σ1(X)/σ3(Y) ratios under C-T-T loadings. The tensile str ranged from 0.02fc to 0.1fc and the stress ratio σ1/σ3 varied from 0.2 to 5. Keeping the tensile stress along the X-direction constant (0.02fc), the dama evolution process of the concrete specimens under different σ1/σ3 stress ratios is shown Figure 13. It can be observed that there was no significant difference between the dama situations of concrete specimens at the peak stresses, as shown in Figure 13a. During following process of damage evolution, there were obvious damage bands along

Compression-Tension-Tension Loading
Compared to the triaxial compressive loading, the loading modes of C-T-T involved changing the stresses along the X-and Y-axes, with constant stress along the Z-axis. Table 4 shows various stress of σ 1 (X)/σ 3 (Y) ratios under C-T-T loadings. The tensile stress ranged from 0.02f c to 0.1f c and the stress ratio σ 1 /σ 3 varied from 0.2 to 5. Table 4. Tensile stress ratios in the X-and Y-axis directions.  Keeping the tensile stress along the X-direction constant (0.02f c ), the damage evolution process of the concrete specimens under different σ 1 /σ 3 stress ratios is shown in Figure 13. It can be observed that there was no significant difference between the damage situations of concrete specimens at the peak stresses, as shown in Figure 13a. During the following process of damage evolution, there were obvious damage bands along the direction of compressive stresses in the Y-Z plane. Furthermore, the damage in the X-Z plane was relatively delayed when the tensile stress along the X-direction was less than that along the Y-direction, as shown in Figure 13b. Then, the damage patterns began to change with the increase of stress ratio σ 3 /σ 1. When the stress ratio σ 1 (X)/σ 3 (Y) was larger than 1 (σ 1 /σ 3 > 1), the damage distribution was uniform and the intersecting cracks developed along the X-axis. The final failure modes under different σ 1 /σ 3 stress ratios can be seen in Figure 13c. With a relatively small stress ratio of σ 3 /σ 1 , the interaction of tensile stresses resulted in damage belts in the X-Z plane developing diagonally. Furthermore, in this circumstance, the X-Z plane was the main failure surface, where the specimen split along the direction of the compressive stress. With the increase of the stress ratio σ 3 /σ 1 , the tensile stresses along the X-and Y-directions interacted with each other and the specimens showed tensile failure with intersecting cracks in each plane. direction of compressive stresses in the Y-Z plane. Furthermore, the damage in the X-Z plane was relatively delayed when the tensile stress along the X-direction was less than that along the Y-direction, as shown in Figure 13b. Then, the damage patterns began to change with the increase of stress ratio σ3/σ1. When the stress ratio σ1(X)/σ3(Y) was larger than 1 (σ1/σ3 > 1), the damage distribution was uniform and the intersecting cracks developed along the X-axis. The final failure modes under different σ1/σ3 stress ratios can be seen in Figure 13c. With a relatively small stress ratio of σ3/σ1, the interaction of tensile stresses resulted in damage belts in the X-Z plane developing diagonally. Furthermore, in this circumstance, the X-Z plane was the main failure surface, where the specimen split along the direction of the compressive stress. With the increase of the stress ratio σ3/σ1, the tensile stresses along the X-and Y-directions interacted with each other and the specimens showed tensile failure with intersecting cracks in each plane. Under an equal stress condition (σ1(X)/σ3(Y) = 1), the internal damage pattern of the concrete specimens showed cracks that developed along the boundary of the aggregates in the direction of the compressive stress (the Z-axis), as shown in Figure 14. Figure 15 shows the damage model inside the concretes along the X-and Y-directions when the stress ratio σ3(Y)/σ1(X) was equal to 5. Obviously, when the applied stresses were not equal (σ1(X) ≠ σ3(Y)), the specimens showed spallation failure and the cracks were perpendicular to the Y-axis (Figure 15a). The damage belts in the X-Z plane turned into damage zones propagating throughout the whole plane, but the cracks in the Y-Z plane were perpendicular to the Y-direction and propagated along the aggregate boundaries, as shown in Figure 15b,c. Under an equal stress condition (σ 1 (X)/σ 3 (Y) = 1), the internal damage pattern of the concrete specimens showed cracks that developed along the boundary of the aggregates in the direction of the compressive stress (the Z-axis), as shown in Figure 14. Figure 15 shows the damage model inside the concretes along the X-and Y-directions when the stress ratio σ 3 (Y)/σ 1 (X) was equal to 5. Obviously, when the applied stresses were not equal (σ 1 (X) = σ 3 (Y)), the specimens showed spallation failure and the cracks were perpendicular to the Yaxis (Figure 15a). The damage belts in the X-Z plane turned into damage zones propagating throughout the whole plane, but the cracks in the Y-Z plane were perpendicular to the Y-direction and propagated along the aggregate boundaries, as shown in Figure 15b,c. tensile stress σ3(Y) was lower than 0.06fc and the tensile stress ratio was less than 1, overall strength increased linearly with the stress ratio. After the stress ratio of σ exceeded 2, the growth of the overall strength slowed down with the increasing st ratio σ3/fc. When the tensile stress was relatively high (>0.08fc), the overall strength longer increased with the continuously increasing stress ratio σ3(Y)/σ1(X) as it approac its extreme strength.  increased, the overall strength of the specimens also improved. In particular, when the tensile stress σ3(Y) was lower than 0.06fc and the tensile stress ratio was less than 1, the overall strength increased linearly with the stress ratio. After the stress ratio of σ3/σ1 exceeded 2, the growth of the overall strength slowed down with the increasing stress ratio σ3/fc. When the tensile stress was relatively high (>0.08fc), the overall strength no longer increased with the continuously increasing stress ratio σ3(Y)/σ1(X) as it approached its extreme strength.  According to the 25 groups of triaxial loading conditions shown in Table 4, Figure 16 shows the corresponding numerical simulation results. The curves in Figure 16 illustrate the compressive strengths of the concretes under C-T-T loadings, taking the σ 3 /f c and σ 3 /σ 1 stress ratios as the respective abscissas. With low tensile stresses in the X-direction, the increasing tensile stress in the Y-direction had a significant influence on the overall strength of the concretes, as represented in Figure 16a. It can be seen that the relationship between the overall strength and the compressive stresses along the Y-direction could be approximated as linear. Furthermore, the slopes of the curves under different tensile stresses of σ 1 (X) varied. The whole strength decreased with the increase of σ 3 (Y). In particular, when the concretes were under lower tensile stress for σ 1 (Y), the overall strength witnessed a more significant decline. The trend of the curves in Figure 16b was completely opposite to that in Figure 16a. As the lateral tensile stress ratio σ 3 (Y)/σ 1 (X) increased, the overall strength of the specimens also improved. In particular, when the tensile stress σ 3 (Y) was lower than 0.06f c and the tensile stress ratio was less than 1, the overall strength increased linearly with the stress ratio. After the stress ratio of σ 3 /σ 1 exceeded 2, the growth of the overall strength slowed down with the increasing stress ratio σ 3 /f c . When the tensile stress was relatively high (>0.08f c ), the overall strength no longer increased with the continuously increasing stress ratio σ 3 (Y)/σ 1 (X) as it approached its extreme strength.

Compression-Compression-Tension Loading
The stress ratios applied to the concretes under triaxial C-C-T loadings are illustrated in Table 5. The X-and Z-axes were the directions of compressive stresses ranging from 0.08fc to fc and the Y-axis was the direction of the tensile stress ranging from 0.021fc to 0.1fc. The absolute values of the stress ratios of the compressive stress to the tensile stress were calculated to be 0.8-48. Table 5. Compressive/tensile stress ratio of the X-and Y-axis directions.  Figure 17a,b shows the failure patterns of the concretes when σ3 = 0.0210fc and σ2 was −0.08fc and −0.25fc, respectively. With low compressive stresses, the overall trend of the cracks was perpendicular to the Z-axis. Furthermore, due to the influence of the compressive stresses acting on the specimen, these cracks intersected with each other. In addition, there were also damage zones in the X-Y planes. When the compressive stress continuously increased to 0.1fc, the failure modes of the concretes made a difference, as represented in Figure 17b. The failure mode changed from a longitudinal splitting failure to a layered splitting failure. Furthermore, the direction of the main cracks perpendicular to the Y-axis turned towards being perpendicular to the Z-axis.

Compression-Compression-Tension Loading
The stress ratios applied to the concretes under triaxial C-C-T loadings are illustrated in Table 5. The X-and Z-axes were the directions of compressive stresses ranging from 0.08f c to f c and the Y-axis was the direction of the tensile stress ranging from 0.021f c to 0.1f c .
The absolute values of the stress ratios of the compressive stress to the tensile stress were calculated to be 0.8-48. Table 5. Compressive/tensile stress ratio of the X-and Y-axis directions.  Figure 17a,b shows the failure patterns of the concretes when σ 3 = 0.0210f c and σ 2 was −0.08f c and −0.25f c , respectively. With low compressive stresses, the overall trend of the cracks was perpendicular to the Z-axis. Furthermore, due to the influence of the compressive stresses acting on the specimen, these cracks intersected with each other. In addition, there were also damage zones in the X-Y planes. When the compressive stress continuously increased to 0.1f c , the failure modes of the concretes made a difference, as represented in Figure 17b. The failure mode changed from a longitudinal splitting failure to a layered splitting failure. Furthermore, the direction of the main cracks perpendicular to the Y-axis turned towards being perpendicular to the Z-axis.
The damage patterns of concrete specimens under different tensile stresses along the Y-direction are shown in Figure 18, where the compressive stress was kept constant along the X-direction at −0.5f c. The compressive and tensile damage zones and the failure modes of the concrete specimens were compared. With low tensile stresses, the failure modes were mainly manifested as compression failures, which indicates that the compressive stress was the main cause of the failure. With the increasing tensile stresses, the tensile damage and compressive damage were evenly distributed. This situation meant that at that time, the compressive and tensile stresses worked together to determine the damage development of the specimens. When the tensile stresses increased to the uniaxial tensile strength, there was almost no compressive damage on the surfaces of the specimens and the failure mode was completely caused by the tensile stresses. It can be observed that the stress ratios of the compressive stress to the tensile stress influenced not only the distribution of the damage but also the failure types. The damage patterns of concrete specimens under different tensile stresses along the Y-direction are shown in Figure 18, where the compressive stress was kept constant along the X-direction at −0.5fc. The compressive and tensile damage zones and the failure modes of the concrete specimens were compared. With low tensile stresses, the failure modes were mainly manifested as compression failures, which indicates that the compressive stress was the main cause of the failure. With the increasing tensile stresses, the tensile damage and compressive damage were evenly distributed. This situation meant that at that time, the compressive and tensile stresses worked together to determine the damage development of the specimens. When the tensile stresses increased to the uniaxial tensile strength, there was almost no compressive damage on the surfaces of the specimens and the failure mode was completely caused by the tensile stresses. It can be observed that the stress ratios of the compressive stress to the tensile stress influenced not only the distribution of the damage but also the failure types. Figure 19 shows the compressive strengths of the concretes under the C-C-T loadings, taking the stress ratios of σ1(X)/fc and σ1(X)/σ3(Y) as the respective abscissas. When the tensile stress was less than 0.0208fc, the overall strength under different compressive stresses was above fc. Furthermore, the overall strength was below fc as the tensile stress was greater than 0.1fc, as shown in Figure 19b. When the tensile stress was close to the median value of the given parameters, which was 0.0625fc, the strength varied between 0.8fc and 1.2fc with the increase of the compressive stresses. Obviously, with low compressive stresses less than 0.5fc, the strength of the concrete specimens increased with the compressive stresses, as shown in Figure 19a. However, the overall strength tended to be constant after the compressive stress exceeded 0.5fc. Then, the whole strength was no longer continuously improved with the increasing lateral compressive stresses.   Figure 19 shows the compressive strengths of the concretes under the C-C-T loadings, taking the stress ratios of σ 1 (X)/f c and σ 1 (X)/σ 3 (Y) as the respective abscissas. When the tensile stress was less than 0.0208f c , the overall strength under different compressive stresses was above f c . Furthermore, the overall strength was below f c as the tensile stress was greater than 0.1f c , as shown in Figure 19b. When the tensile stress was close to the median value of the given parameters, which was 0.0625f c , the strength varied between 0.8f c and 1.2f c with the increase of the compressive stresses. Obviously, with low compressive stresses less than 0.5f c , the strength of the concrete specimens increased with the compressive stresses, as shown in Figure 19a. However, the overall strength tended to be constant after the compressive stress exceeded 0.5f c . Then, the whole strength was no longer continuously improved with the increasing lateral compressive stresses.

Modification of Ottosen Parameters
According to the numerical results discussed above, the strengths of the concretes under multiaxial loadings were a function of the triaxial stress states and could not be predicted only by the independent limits of pure tension, compression and shear stresses Therefore, in order to correctly estimate the concretes' strengths, the failure criterion of concrete materials under the combined complex stress states was expressed using a formula. By combining the failure criterion with the constitutive relationship of concrete materials, the mechanical behaviour of concretes can be predicted and then it can be used in engineering analysis and optimisation design under complex multiaxial stress states.
The gradient descent method is widely used in machine learning; it is a common firstorder optimisation method which is one of the most simple and classical methods that are used to solve unconstrained optimisation problems. The main purpose of this method is to find the minimum value of the objective function through iteration. The first order

Modification of Ottosen Parameters
According to the numerical results discussed above, the strengths of the concretes under multiaxial loadings were a function of the triaxial stress states and could not be predicted only by the independent limits of pure tension, compression and shear stresses. Therefore, in order to correctly estimate the concretes' strengths, the failure criterion of concrete materials under the combined complex stress states was expressed using a formula. By combining the failure criterion with the constitutive relationship of concrete materials, the mechanical behaviour of concretes can be predicted and then it can be used in engineering analysis and optimisation design under complex multiaxial stress states.
The gradient descent method is widely used in machine learning; it is a common first-order optimisation method which is one of the most simple and classical methods that are used to solve unconstrained optimisation problems. The main purpose of this method is to find the minimum value of the objective function through iteration. The first order means using only the first derivative of the target function and not its higher-order derivatives. Based on the numerical simulation data obtained in Section 3, this study adopted the gradient descent method and took the Ottosen parameter obtained from the classical experimental results as the initial value [48].
The machine learning method was used to update the Ottosen parameters to obtain the yield surface with higher precision. The new partial plane yield curve further validated the numerical simulations of the concretes under triaxial loading that were undertaken in this study. The specific methods and steps are described below.

Calculating the Invariants
Based on the 65 sets of triaxial numerical simulation results obtained in Section 3, the invariants of the stress and deviatoric stress tensors were calculated according to the principal stresses (σ 1 , σ 2 , σ 3 ) of the above simulations: where I 1 , I 2 and I 3 are the first, second and third invariants of the stress tensor, respectively. Furthermore, where J 2 and J 3 represent the second and third invariants of the deviatoric stress tensor. In order to meet the geometric requirements of the yield surface of concrete materials, Ottosen [49] proposed a yield function that is denoted by three stress invariants I 1 , J 2 and θ: where a, b, k 1 and k 2 are constants. The yield function of the concrete materials was obtained in the Haigh-Westergaard coordinate system, which is defined by the cylindrical coordinates ξ, ρ and θ representing the hydrostatic component, the deviatoric component and the Lode angle component, respectively. These coordinates are functions of the invariants I 1 , J 2 and J 3 expressed in terms of the principal stress tensors σ 1 , σ 2 and σ 3 as follows:

Fitting Parameters
The gradient of each parameter (a, b, k 1 and k 2 ) is calculated using the chain rule [50]: In order to make the fitted parameters more accurate and narrow the range of the iteration values, the gradient descent method needs to give initial values to the unknown parameters. The initial parameter values used for the calculations in this study were derived from [49]. The four initial parameters (a, b, k 1 and k 2 ) were determined based on the experimental results of four kinds of concrete tests, including the uniaxial compressive strength of f c , the uniaxial tensile strength of f t , the biaxial compressive strength of f bc and the triaxial stress state at the compression meridian (θ = 60 • ). Then, the four initial parameters (a, b, k 1 and k 2 ) and the three stress invariants (I 1 , J 2 and cos3θ) obtained from the results of the numerical simulations in Section 4.1 were substituted into Equations (22) and (23). In machine learning, the optimal weight parameters can be found based on a certain index, which is used to indicate the current state. The index is defined as the loss function, which represents the unfitting level of the data. In this study, the loss function was defined as follows [51]: The loss values of all groups were calculated, and then the absolute values of all losses were averaged, which can be expressed as |loss|.
The learning rate of the machine learning lr was set to 0.001 [52]. The updating formula [50] for the parameters is given in Equation (33): The updated parameters were plugged into Equations (22) and (23) again and the yield function f(I 1 , J 2 , θ) was iterated 100,000 times to obtain the final values of the updated parameters. The initial and updated parameters, as well as the corresponding loss values, are listed in Table 6. Based on using the original data as the numerical results of the triaxial stress states, it was found that the accuracy of the new parameters obtained using the machine learning method were 4.25 times higher than those of the original parameters (see Table 6).

Validation of the Modified Ottosen Parameters Based on the Triaxial Experimental Tests
The modified Ottosen parameters can provide valuable information, not only for modification purposes of the failure criteria but also for predicting the macroscopic strength of concretes under multiaxial stresses. In order to verify the validity and accuracy of the modified failure criterion, 97 sets of data obtained from the triaxial experimental tests [53] were substituted into the failure criterion. The first step in such a validation was the calculation of invariants based on the experimental data and then the coordinate points (ξ, ρ, θ) were plotted on the yield surface curve on the deviatoric plane and the failure surface in the principal stress space.
The new parameters obtained in Section 4.2 were put into the Haigh-Westergaard stress space coordinate system: When the yield curve in the partial plane was plotted using the coordinate points (ξ, ρ, θ), the angle θ of the partial plane was divided equally by the angle interval of π 120 . The value of ξ ranged from 0.13 to −10.5 at intervals of 0.3. Then, the calculated coordinate points of (ξ, ρ, θ) were substituted into Equation (34) and the values of ρ in the implicit function were solved.
The coordinate points of (ξ, ρ, θ), which were obtained using calculations based on the triaxial stresses (σ 1 , σ 2 , σ 3 ) from the experiments, were plotted on the yield surface curve on the deviatoric plane, which was modified by the Ottosen parameters, as shown in Figure 20. Furthermore, the transformation relationship of the coordinates between the Haigh-Westergaard stress space and the principal stress space is as follows:

Conclusions
This study focused on the mechanical behaviour of concretes under triaxial loading conditions and was devoted to characterising the damage evolution of the triaxial behaviour of concrete specimens under various triaxial stress states. Triaxial stresses (including tension and compression) were proportionally applied to the concrete specimens to perform numerical experiments. The following conclusions can be drawn: The coordinates of (ξ, ρ, θ) above were transformed into the coordinates of (σ 1 , σ 2 , σ 3 ) and these points were plotted on the failure surface in principal stress space, as shown in Figure 21. It can be observed that the experimental data points were basically located on the partial plane curve, which verified the validity of the simulation results of the concrete specimens under triaxial loading. At the same time, the modified parameters were also shown to be reasonable and reliable.

Conclusions
This study focused on the mechanical behaviour of concretes under triaxial loading conditions and was devoted to characterising the damage evolution of the triaxial

Conclusions
This study focused on the mechanical behaviour of concretes under triaxial loading conditions and was devoted to characterising the damage evolution of the triaxial behaviour of concrete specimens under various triaxial stress states. Triaxial stresses (including tension and compression) were proportionally applied to the concrete specimens to perform numerical experiments. The following conclusions can be drawn: (1) The mesoscopic modelling with good robustness gave reliable and detailed damage evolution processes under different triaxial stress states. Rich simulation data were systematically obtained to identify the mechanical responses and damage mechanisms of the concretes. The results show that the multiaxial mechanical behaviour of the concretes was significantly influenced by the level of confining pressures and the lateral tensile stress. Under the combination of tensile and compressive stress, the lateral tensile stresses rapidly aggravated the decline of concrete strength. (2) The triaxial stress ratios applied to the concrete specimens significantly influenced the damage evolution and failure patterns. Under C-C-C loadings, the direction of the crack propagation mainly depended on the larger lateral compressive stress. Under C-T-T loadings, the initiation and progression of damage belts directly depended on the stress ratios and the intersecting cracks were distributed evenly inside the concretes where the ratio was closer to one. Under C-C-T loadings, the compressive damage inside the concretes gradually transformed into tensile damage with the increasing lateral tensile stress. (3) The Ottosen parameters were modified using the gradient descent method. The fitting results between the experimental data and the failure surface verified the validity of the triaxial numerical simulations in this study. The validation results indicated that the gradient descent method with good data processing ability was suitable for analysing parameters quickly. Furthermore, the modified Ottosen failure criteria can be used to analyse the mechanical performance of concrete structures under further complex multiaxial stress states.
Finally, further work is needed, where it is likely the influence of mesoscopic components on the macroscopic strength and damage evolution law can be investigated when concretes are under uniaxial and multiaxial loadings, respectively. In the mesoscopic structures of concretes, the identification of the mesoscopic parameters of concrete components is still an open topic. Further investigations on how to establish the quantitative relations between the macromechanical properties of concretes and the mesoscopic parameters should be carried out.

Conflicts of Interest:
The authors declare no conflict of interest.