Hydraulic fracturing (HF) is a breakthrough technology to increase the permeability and productivity of the reservoir in oil/gas exploitation or enhanced geothermal system (EGS). The stimulation mechanism is not fully understood as it occurs in the deep subsurface (depth of 2-6 km). This research developed a novel experimental method to investigate the properties of consolidated porous materials and integrated it with fully coupled hydro-mechanical numerical models to study the elastic-plastic fracturing mechanisms during the hydraulic stimulation process. In this research, first, an alternative Gas Expansion Induced Water Intrusion Porosimetry (GEIWIP) method was experimented on different concrete samples with different porosities ranging from 10 to 20%. The GEIWIP method was identified as a useful technique to investigate the porosity and pore size distribution evolution due to the growth of the hydro-fracture/s through the samples. Next, a finite volume numerical model of a real rough fracture was simulated to study the heterogeneities effect on the fluid flow within the fracture. The results were validated against the experimental investigations, showing that the proposed constitutive model is reliable to simulate the fluid flow within the fracture. Another scenario was a finite element simulation of a fluid flow within the rock matrix in a fractured reservoir. The results were validated against the existing benchmark problems such as geothermal doublet and thermohaline aquifer-aquitard-aquifer and the corresponding numerical models were shown to be reliable in porous media modeling. Next, a fully coupled finite element analysis with the phase-field approach was carried out to capture fractures initiation and elastic-plastic propagations in a punch through shear (PTS) test on a mechanically loaded granite sample. The phase-field method combined with the finite element analysis showed the elastic-plastic fracture propagation within the rock matrix and the results were in agreement with the PTS tests. The sample failure was investigated after the sample underwent plastic deformation. Finally, the hydraulic fracturing phenomenon in a fractured deformable rock formation was modeled to investigate different poroelastic responses in different wellbore orientations. This modeling was done on one of the active hydraulic fracturing sites in Texas, USA.