The analysis of fracture propagation in geomaterials (such as rock, concrete) still remains a challenge from the perspective of numerical modeling, which needs to capture micro crack (damage) inception and growth up to the point of macro fracture initiation and propagation, and the formation of fracture paths after the global softening in later simulation stages. State-of-the-Art finite element tools are either based on Continuum Damage Mechanics (CDM) or on Fracture Mechanics (FM). In CDM, fractures are considered as the ultimate stage of damage accumulation, at which the energy dissipated in the damage zone equals the energy that needs to be released to create new surfaces (
). In FM, discontinuous displacement fields propagate through a continuum according to physics-based criteria (
). FM theories are usually implemented in Cohesive Zone Models (CZM) or Extended Finite Element Methods (XFEM). Because CDM cannot be used to simulate fracture surface debonding at late crack propagation stages and since State-of-the-Art CZM and XFEM do not account for micro-crack propagation prior to fracture surface debonding, several researchers coupled CDM and XFEM (
) or CDM and CZM (
). However, even in these coupled models, there is no rigorous calibration of the critical damage threshold that marks the transition between micro-crack propagation and macro-fracture initiation. Moreover, the fracture path depends on the location of the CZ elements or on the interpolation order chosen in the XFEs. Thus in the present study, we propose a numerical method that couples a CDM model (for the bulk) to a Cohesive Zone Model (for the fracture). We present the theoretical framework and calibration method for the CDM model and for the CZM in Sections 2 and 3, respectively. In Section 4, we present the results obtained during simulations of biaxial compression tests in which CZs are placed at the boundaries of all the Finite Elements of the mesh.