新华夏构造体系结构面“米字型”分布与演化的数学模拟研究

吕古贤，王红才,2，韩 璐，张宝林，胡宝群，吕承训，马立成，焦建刚，毕珉峰

(1.中国地质科学院地质力学研究所，北京 100081；2.新构造运动与地质灾害重点实验室，北京 100081；3.中国科学院矿产资源研究重点实验室，中国科学院地质与地球物理研究所，北京 100029；4.东华理工大学 地球科学学院，江西 南昌 330013；5.中国地质调查局发展研究中心，北京 100037；6.长安大学 地球科学与资源学院，陕西 西安 710054；7.中国地质科学院矿产资源研究所，北京 100037)

0 引 言

新华夏构造体系是东亚地区特有的地质构造体系[1-6]。依据地质力学的观点，用新华夏构造体系分析方法，研究长江中下游铁铜矿带、胶东玲珑金矿田和焦家金矿田、豫西栾川南泥湖钼钨矿田、江西相山铀矿田、赣南银坑银多金属矿田和盘古山钨矿山等区域的构造控岩控矿和构造成岩成矿规律[7-18]，经过详细的构造变形岩相地质填图或编图、构造测量和应力场反演等研究工作[11,19-23]，发现上述地区构造在中—新生代受新华夏构造体系的应力场控制，构造体系的结构面呈现“米字型”分布样式。

新华夏构造体系的结构面在多个层次和不同规模上，较普遍呈现“米字型”样式。构造体系的结构面发育了较为规则的几个分布方向，具有相对稳定的力学性质。“米字型”样式系由NNE 25°方向挤压断裂和褶皱带、NNW 345°方向(大义山式)张扭断裂、NEE 75°方向(泰山式)压扭构造和NWW 300°方向(长江式)的横张构造组成[24-25]。

本文利用有限元法对新华夏系“米字型”构造样式进行应力场模拟，对比分析不同阶段结构面(带)的位移、主干拉应力、主干压应力、剪切应力、应变等参量的分布特征。实验揭示了“米字型”构造的形成过程，即是剪切变形、挤压变形和引张变形三个构造阶段的产物。对于构造体系方向和演化等理论难题的解决，这是一项重要的进展。

特别需要说明的是，在分析讨论中，为兼顾地质专家用语，采用了主干拉应力(对应于本文软件计算结果的主应力σ1，依照弹性力学规定，拉应力为正数，压应力为负数)、最小主压应力(对应于本文软件计算结果的主应力σ2)、最大主压应力(对应于主应力σ3)、剪切应力(剪应力τ)等名词，为避免混淆，文中的讨论分析均采用文字叙述。

1 新华夏构造体系“米字型”构造及分阶段形成的地质模型

1.1 新华夏构造体系“米字型”构造

野外地质调查和资料综合分析表明，新华夏构造体系的“米字型”结构(图1)的地质特征表现为以下4种构造。

图1 新华夏构造体系结构面组成的“米字型”示意图 Fig.1 Schematic diagram of the double-cross surface in the Neo-Cathaysian structural system1.新华夏挤压构造；2.“泰山式”压剪构造；3.“大义山式”张扭构造；4.“长江式”张性构造；5.等效作用力

(1)新华夏主干构造为北东25°～30°方向大型褶皱、隆起带和沉降带及伴生的压扭构造岩相带；以整体挤压变形为特征，生成时期跨度在印支—燕山运动及其以来，其中强烈活动时代集中于120～140 Ma。

(2)泰山式构造早期表现为北东东75°左右的扭性褶皱和断裂构造岩相带，压剪变形，活动时代在140～203 Ma；晚期表现为北东东75°左右的断裂构造岩相带，以花岗岩带、断陷盆地的线性展布为体现，压剪变形，活动时代在100～120 Ma。

(3)大义山式构造为北北西345°左右的扭性断裂构造岩相带，张剪变形，活动时代与泰山式和长江式构造活动相近；较晚期构造以岩浆活动为主，控制雁列构造和岩脉的分布。

(4)长江式构造早期显示为北西西向290°～300°的张性断裂构造-岩相带，张裂变形，具有规模性的岩浆侵入，活动时代在130～144 Ma；晚期长江式构造，北西西向295°～300°向张性断裂带，活动时期为100～120 Ma。

1.2 新华夏系“米字型”构造的形成阶段

大量野外地质调查和观察结果表明[3,5,24]，新华夏系“米字型”构造的形成包括了如下几个阶段(图 2)：第一期，共轭构造发育阶段，包括NNW方向(大义山式)张扭断裂和NEE方向(泰山式)压扭构造。第二期，主压扭构造发育阶段，由NNE方向断裂带构成新华夏构造体系的主压扭带，形成 “3隆3拗格局”。第三期，横张构造发育阶段，发育NWW方向(长江式)的张裂构造。

图2 新华夏构造体系结构面 “米字型”的形成阶段图Fig.2 Schematic diagram showing the evolution stages of the double-cross surface of the Neo-Cathaysian structural system1.新华夏挤压构造;2.“泰山式”压剪构造；3.“大义山式”张扭构造；4.“长江式”张性构造；5.等效作用力；(a)NNW和NEE向共轭构造带；(b)NNE向主干构造带；(c)NWW向构造带

2 构造体系“米字型”构造的有限元模型

2.1 有限元方法(FEM)

结构体的应力应变分析法，包括有限元法、离散元、边界元、解析法等多种方法，其中有限元法有其不可比拟的优势。由于野外地质问题、岩石的力学性质非常复杂，在利用有限元法进行模拟分析时，需要做适当简化，有效解决问题，又避免了模型计算过于复杂而增加巨大困难[26-29]。

国内外许多文献中，多数应力场的数学模拟仅仅限于模拟岩石的弹性变形。尽管也可以按粘-弹性材料进行计算，但粘-弹性问题是采用流变方程进行计算的[30]，计算结果与时间紧密相关，存在巨大的伸缩性(或说具有不确定性)。因而，国内外绝大多数所谓的粘-弹性构造应力场有限元数学模拟，不做粘性问题计算，而只是按弹性问题处理[19,31-33]。显然，应力场数值模拟比物理实验模拟具有明显的灵活性优势[34]。总体上，地质问题模拟过程，可大致分为四个重要阶段：(1)客观建立地质模型和获取参数；(2)将地质模型简化为合理的数学力学模型；(3)划分单元网格进行数值模拟；(4)结果后处理分析[19-20,31,35-41]。

采用广为流行应用的ansys等应力应变模拟软件。

2.2 有限元模型

构造体系属于新华夏构造应力场的产物，本文选择代表性金属热液矿田(例如银坑矿田和玲珑矿田等)[7-8,11,13]，进行有限元法对新华夏系构造不同构造带(面)进行应力和应变场及其分布的模拟分析。

需要指出的是，尽管新华夏系构造不同构造带(面)涉及到不同时期形成演化、遭受了不同构造应力环境作用过程、且各分支具有不同的力学性质差异性，不具备严格的数学力学可比性，但一个非常重要的事实是：在野外地质工作中，我们会同时面对许多不同地质时期形成的构造带(面)，它们共同存在、相互叠加、许多时候难以区分。更为重要的是，它们要被放在一起加以分析、对比。因此，将它们置于相同或相近的背景下进行对比研究，寻找其中的规律性，以便为地质找矿服务。

为了便于对比分析，本文将新华夏系构造不同构造带(面)置于相同级别的应力环境，也即具有相同的莫尔圆空间和相近的边界约束下，模拟反演计算，并综合分析讨论。

实验观察到，自然界中岩石材料受力后的变形性质，不仅仅是弹性的，还表现出粘性、塑性，以及强度峰值之后的应变软化性质等[41]。为了方便起见，本文把岩石材料视为连续弹性体进行处理。有关地质和岩石力学参数，参考了中国大陆主要岩石类型物性参数库[42]。围岩和构造带赋予不同的岩石力学参数，以区分其力学性质。由于实际地质问题非常复杂，本文进行简化，以理想化的模型开展模拟分析。具体参数确定为：围岩E=5×103MPa，ν=0.25，ρ=2 060 kg/m3；构造带E=1×103MPa，ν=0.30，ρ=1 700 kg/m3；本模拟反演过程中，ANSYS软件采用弹性力学假定，规定拉应力为正、压应力为负。

边界载荷确定：基于大量野外地质观察，认为新华夏构造体系不同阶段的主压应力方向为北西—南东到北西西—南东东方向的渐变[1]，因此，考虑3～5 km深度范围的平面应变地质模型，第一期模型在NW 45°边界AB上施加北西—南东向挤压应力10 MPa载荷，第二期模型在NW 60°边界AB上施加北西西—南东东向主压应力10 MPa挤压应力载荷，第三期受到北西西—南东东向挤压应力，在NW 60°边界AB上施加10 MPa挤压应力载荷。

在模型上施加边界约束条件的设定：第一期，在NW 45°边界AB的对面边界CD边上施加滑轮约束并原点固定，其它边界为自由边界。第二期，在NW 60°边界AB的对面边界CD边上施加滑轮约束并角点固定，其它边界为自由边界。第三期与第二期相同[35]。

同时，在建立模型时，为了避免边界效应影响导致的计算结果偏差，在划分网格时，将模型的外围进行适当扩边。在此基础上建立有限元模型，划分二维网格节点1 458个，采用四边形及三角形混合剖分，共划分989个等参元，如图3所示。

图3 应力场有限元模型Fig.3 Stress field finite element model (a)第一期，NNW向构造和NEE向构造带共轭发育；(b)第二期，NNE向主干构造带；(c)第三期，NWW向“长江式”构造带

3 有限元法模拟结果与讨论

3.1 分阶段模拟结果

各阶段的局部位移、应变和应力具有明显不同的变化(图4、图5、图6)。

3.1.1 第一期剪切共轭阶段

受到NW向挤压应力的作用，区内应变最大的地方主要在模型的西北侧，NW向断裂的北端附近。断裂带内位移比围岩内量值较大，因此，位移等值线在断裂带内有比较明显的转折。围岩内应变较均匀、且量值较小，而断裂带内发生的应变量值大，断裂带内主要呈现为NW向挤压和NE向拉伸作用并具有较强烈的应变(图4(a))。例如银坑矿田断裂带内发育的断层泥、断层角砾和透镜体等，都是强烈变形的结果。

最小主压应力，在区内呈现NW向和SE向量值较大，向中心逐渐减小的趋势，而在断裂的两端部位置呈现高应力区(应力集中)，应力最大值约为5.58 MPa，最小值为0.054 6 MPa。NNW向断裂带内应力量值较小，NEE向断裂带内应力量值较大，最大约为0.833 MPa。

最大主压应力分布(图4(c))，最大为19.822 MPa，最小为1.289 MPa，围岩总体较大，断裂带内总体较小，为2～8 MPa，NEE向断裂带内明显高于NNW向断裂。在断裂带两端的围岩内最大主应力均较高。说明在新华夏应力环境下，围岩受到的挤压应力最大，断裂带内受到的挤压应力比围岩小，而NNW向断裂带受到的挤压应力又比NEE向构造带小。而且高应力区主要分布在断裂带的两端，而并非断裂带的交汇部位。

区内剪切应力最大为8.121 MPa，最小为0.796 MPa(图4(d))，在围岩中表现为北东侧和南西侧较高，向中间降低的趋势，断裂带内则为北西侧和南东侧较高，中心较低的趋势，且NNW向断裂带受到的剪切应力总体比NEE向断裂带低。在断裂的两端，有较高的剪切应力分布。

图4 新华夏构造体系第一期构造带应变应力分布图Fig.4 Strain and stress distribution of stage-1 Neo-Cathaysian structural system

由于模型整体受到压应力作用为挤压应力环境，因此，拉张应力总体较小，在NEE向断裂带的两端，有较高的拉张应力，最大值为3.35 MPa，NNW向断裂带内较弱并由南到北呈现增强，为0.2～2.0 MPa。

野外地质工作表明，NNW方向(大义山式)张扭断裂和NEE方向(泰山式)压扭构造附近，发育着大量基性-中性岩石及相关矿床，表明其对该类矿床有明显的控制作用[43-44]。

3.1.2 第二期挤压阶段构造

从位移场看，模型总体受到NWW向挤压应力的作用，在围岩中呈现比较规律的NWW向变形，在断裂带内呈现明显的NWW向压缩变形和NNE向拉伸变形，相比之下，围岩中变形(应变)量值大，且以NWW向压缩为主(图5(a))。

最小压应力在围岩内较为规律，呈现明显的从北西侧和南东侧向断裂带降低的趋势，最大值约为3.241 MPa，最小值为0.045 MPa。而在断裂内的分布规律与围岩明显不同，在断裂带内表现为中心最高并向两端递减的趋势(图5(b))。

图5 新华夏构造体系第二期构造带应变分布图Fig.5 Strain and stress distribution of stage-2 Neo-Cathaysian structural system

最大压应力分布(图5(c))，模型内最大压应力总体呈现从断裂的北东边部和南西部两端向中心位置递减趋势，最大值为20.097 MPa，最小值为5.776 MPa，断裂带两端的围岩位置也出现最大压应力值，断裂带内，8.4～11.2 MPa，相比最大压应力在断裂带内稍低于围岩。

区内剪切应力(图5(d)、(e))最大为7.895 MPa，最小为2.719 MPa，总体呈现出从区域北西侧和南东侧向断裂带逐渐增加趋势，围岩中断裂带两端附近位置出现剪切应力高值区。断裂带内，剪切应力总体较小，为2.40～3.57 MPa，但是，断裂带内的剪切应力为4.80～6.12 MPa,相比之下断裂带内的剪切应力大于围岩中，以剪切破坏为主导。

模型内拉张应力量值总体较小，介于0～4.23 MPa之间。拉张应力主要发育在围岩中，拉张应力最大为4.23 MPa，拉张应力等值线近似沿断裂带平行方向条带状展分布，离开断裂带的方向上迅速减弱并在2倍断裂带范围内接近消失。而在断裂带内，拉张应力极微弱，可以忽略，这归因于本阶段为挤压环境阶段主要产生挤压应力。

f野外地质工作发现，主压扭构造发育阶段，由NNE方向断裂带构成构造体系的主压扭带，形态上有“3隆3拗格局”，它对于沉积型的油气、煤炭和矿床具有重要的控制作用。

3.1.3 第三期“长江式”构造带引张阶段

总体模型受到NWW向挤压，发生NWW向压缩变形，同时，在与压缩变形垂直的NNE向出现不同规模的拉张变形，拉张变形以断裂带内较为明显且分布较为均匀，围岩内拉张应变量值较小且分布不均匀为特征(图6(a))。

应力场总体呈现NWW—SEE的挤压应力主导，有不规则变化，断裂带内应力较均匀，且发育NNW向微弱张应力。模型内最小压应力 (图6(c))，量值上较小，为0.036～3.179 MPa，呈现从模型的北西侧中部和南东侧中部断裂带端部向围岩内逐渐增大的趋势，围岩中心部位趋于最大且分布较均匀。在断裂内的最小压应力分布规律与围岩基本一致。

图6 新华夏构造体系第三期构造带应变分布图Fig.6 Strain and stress distribution of stage-3 Neo-Cathaysian structural system

最大压应力分布(图6(c))，区域围岩内最大压应力比断裂带内大，最大压应力的峰值最高为2.4～16.5 MPa，且在断裂带端部两侧的围岩内出现高峰值并向周围内迅速降低。同时，最大压应力σ3的最小值也分布在断裂带延伸方向围岩内、向周围呈近似环带状分布并减弱。在断层内,最大压应力明显小于围岩，为2.3～4.0 MPa，断层中心部位为局部低值区，向两端逐渐增大，具有近似对称分布。

区内剪切应力(图6(d))介于1.5～6.8 MPa之间，总体上剪切应力在围岩内量值较大而在断裂带内较小。围岩内，剪切应力高值区出现在断裂带的端部附近，高达6.8 MPa。最小值出现在断裂带内。在断裂带内中心地带为低值区。

模型内拉张应力量值总体较小，介于0.036～2.280 MPa之间。围岩内量值普遍较小，在断裂带端部延伸处存在，<0.027 MPa；在断裂带内靠近端部出现拉张应力，量值较小，<0.24 MPa。

结合野外地质工作，发现NWW方向(长江式)的张裂构造对于局域型地层、岩浆岩和相关矿床有控制作用[43-44]。

3.2 新华夏系应力应变场的分析探讨

新华夏系“米字型”构造的上述数学模拟，也是构造体系的分阶段演化和形成过程的反演，这展示了各阶段构造带及其地质力学特性的内在协调性。基于上述新华夏系应力应变场模拟分析，可得到以下几点规律性认识：

(1)三阶段相比，构造带内主干拉应力以第一期NNW向构造带内为最大，第三期NWW向“长江式”构造带次之，第二期NNE向和第一期NEE向构造带内均最小，非常微弱，接近于零。

(2)三阶段相比，第一期NEE向构造带中等，最大压应力以第二期NNE向构造带为量值最大，第一期NNW向构造带和第三期NWW向“长江式”构造带为最小。在银坑矿田内即是如此[11,45]：NNE向主干推覆构造表现为强烈的挤压性质，不仅发生强烈的挤压推覆作用，而且在断裂带附近有明显的断层角砾和断层泥发育，局部地区可以见到透镜体发育。NWW向“长江式”构造带为区内主要控矿构造，矿脉和岩脉多显示出张性特征，说明该组构造受到的主干应力为张性，挤压应力较小。

(3)新华夏构造体系三期构造带，最小主压应力相比，第二期NNE>第三期NWW>第一期。第一期构造带最小主压应力，量值最低，为0.054 6～0.833 0 MPa，且在北西向断裂带量值较小、在NEE向带较大。第二期构造带内的最小主压应力，量值最高，分布为0.045～3.241 MPa。第三期NWW构造带内最小主应力，量值居中间，0.036～3.179 MPa。

(4)三期构造带均受到均匀分布的近似NW—SE向挤压应力作用，而第三期NWW“长江式”构造带则在挤压应力作用的同时还有NE—SW向拉张应力作用，产生较大的NE向拉张应变。表明构造体系的基本构造应力场的性质是复杂的，在外力、约束、岩石材料多种因素作用下宏观挤压环境内也会出现引张变形导致拉张构造，有利于成矿热液运移聚集。

(5)以上三期阶段构造的应变应力大小、方向及分布在围岩和断裂带内对比分析表明，新华夏系的构造应力应变场有其特有规律性，各演化阶段构造带的位移、应力和应变的变化具有显著的标志性。

4 结 论

(1)本文解析了新华夏构造体系的“米字型”构造样式：新华夏构造体系，普遍发育着“米字型”构造；“米字型”构造显示出分阶段的发育演化特征；依据野外区域构造应力特征、应力边界变化及其空间的关联性，“米字型”构造各组成部分应该归为同一构造体系分阶段或分期次的产物。

(2)通过计算机有限元法反演，获得了“米字型”构造应力和应变的变化特征：新华夏构造体系三期应变特征相比表现为，第一阶段发育NNW向和NEE向共轭应变带、第二阶段发育NNE向挤压应变构造带、第三阶段发育NWW向张剪应变带。新华夏构造体系三期构造带，最小压应力量值大小相比，第二期>第三期>第一期。第一期最小主应力在NNW向断裂带，量值较小，在NEE向带量值较大，为0.054 6～0.833 0 MPa，为三期中最小。第二期构造带最小压应力值分布为0.045～3.241 MPa，为三期中最大。第三期构造带内最小主应力，0.036～3.179 MPa，居三期中间。三阶段相比，构造带内主干拉应力以第一期NNW向构造带内为最大，第三期NWW向“长江式”构造带次之，第二期NNE向和第一期NEE向构造带内均微弱、几乎接近为零，表现出，由老到新，最高下降至微弱之后再回升的趋势。最大压应力第一期NEE向构造带中等，第二期NNE向构造带为最大，第三期NWW向“长江式”构造带和NNW向构造带为最小。表现出时间由老到新，遵循先增强至最高值再下降至最小的趋势。

(3)本研究模拟了新华夏构造体系结构面或构造带的分布特征，反演了三期演化特点，揭示了构造带应变应力的大小和方向，演示了它们在围岩和断裂带的变化规律。研究证明，构造体系的应力和应变场具备特定的分布演化规律；构造体系的构造带成生演化的位移、应力和应变具有标志性意义。构造体系结构面(带)的分布变化不仅能指导应力-应变特征研究，而且有助于准确预测成矿方向和成矿阶段与发现成矿部位。

致谢：感谢王连捷研究员和刘瑞珣教授对于“构造附加静水压力”研究的长期支持和具体指教，由于王连捷研究员的具体指导才完成本模拟的设计、计算和分析。