图7.2　含节理方板

如图7.2 所示，宽为1 m、高为3 m 的平板，中间含有一条贯通节理，底部固定，上部受均布压力σ 和切向力τ 的作用。平板材料参数为弹性模量E=2000 MPa，泊松比为v=0.3，节理面摩擦系数为μ=0.3，凝聚力为c=0。将整个计算区域划分为21×31个背景积分网格，网格角点作为节点，采用4×4 高斯积分。将节理面划分为20 个区段，每个区段采用4 个高斯积分点。首先计算σ=50 MPa，τ=0 Pa 的情况。当惩罚因子En=Eτ的取值依次选取为E～1012E，计算的法向接触应力值与精确解（σ=50 MPa）均吻合得很好，最大相对误差仅为0.77%（当En=Eτ=E 时）；图7.3给出了沿节理面上的法向接触和切向接触应力分布（En=Eτ=109E 时）。可以看出，本文计算的切向接触应力远小于法向接触应力，其分布反对称于节理面。图7.4 和图7.5 分别给出了平板的变形图（放大30 倍）和位移云图，可以看出，计算结果与没有裂纹时的情况基本一致，这表明附加变量值基本都为零。然后计算σ=50 MPa，τ=5 MPa 的情况。图7.6 为沿节理面上的法向接触和切向接触应力分布。可以看出，切向接触应力为抛物线型分布，法向接触应力为线性分布。以上两种情况的计算结果与文献［133］十分接近。

图7.3　σ=50MPa，τ=0Pa节理面接触应力

图7.4　平板变形图

图7.5　平板位移云图

图7.6　σ=50MPa，τ=5MPa时节理面接触应力