CN112037297A - 一种基于网格重组的混合光传输建模方法 - Google Patents

Abstract

本发明公开了一种基于网格重组的混合光传输建模方法，根据生物组织的光学特性，将各组织划分为高散射组织区域和低散射组织区域；根据各组织所属区域将动物模型离散得到的由空间位置排序的节点和四面体重新排序，得到两个独立网格；基于光传输理论和界面耦合，把两个网格合成一个重组网格，建立混合光传输模型。光传输模型用于前向仿真，能够由光源信息和生物组织结构信息获得生物表面光分布，用于逆向仿真，能够由生物表面光分布和组织结构信息重建光源。本发明能够在精度和效率之间取得更好的平衡，适用于常见光学分子成像模态，包括生物发光断层成像、激发荧光断层成像，契伦科夫光学成像等，在医学分子影像、重建方法等领域有重要的应用价值。

Description

一种基于网格重组的混合光传输建模方法

技术领域

本发明属于分子影像领域，涉及一种基于网格重组的混合光传输建模方法，适用于常见的光学分子成像模态中的光传输建模，包括生物发光断层成像、激发荧光断层成像，契伦科夫光学成像、X射线发光断层成像等。

背景技术

光学分子影像是一个快速发展的生物医学成像领域，将光学过程与一定的分子性质相结合，利用特异性探针的光学性质，反映活体状态下分子水平变化，用于组织病理变化的早期检测。该成像方式具有无辐射、反馈快、灵敏度高、时空分辨率高、特异性高等优点，常被用于小动物研究中的分子、细胞和基因表达成像，促进药物开发、疾病研究和治疗干预。

常见的光学分子成像模态，根据成像过程是否需要外部能量激发，分为直接光学分子成像和间接光学分子成像。直接光学分子成像是无需外部能量激发的成像技术，主要通过探测特异性分子探针或放射性核素释放的光子实现成像，包括生物发光断层成像和契伦科夫光学成像。间接光学分子成像是需要外部能量激发的成像技术，如用激光作为外部能量源的荧光分子断层成像，用X射线作为激发源的X射线发光断层成像。

对于各模态光学分子断层成像中，光传输模型用于描述光在生物组织中的传输过程，其准确性和运算效率对光源的快速精确重建起着决定性作用。辐射传输方程(RTE)虽已经被较为广泛地用作描述可见光和近红外光子在生物组织中传播的标准模型。然而在实际应用中，对于复杂的生物组织，RTE的实现是极其复杂的，且需要大量的计算时间。一些RTE近似模型先后被提出，用于模拟光在非均匀介质中的传输，如扩散方程(DE)、简化球谐近似方程(SP_N)、离散序数方程(S_N)、球谐方程(P_N)和相位近似(PA)。然而，DE计算量小、求解速度快，但只适用于高散射特性的生物组织。高阶近似S_N、P_N和PA，可以精确地模拟光在生物组织内的传输，但较低的计算效率限制了它们在复杂非均匀性介质中的应用。SP_N在所有RTE高阶近似模型中，计算量相对较小，更适用于描述光在低散射特性介质中的传输过程。已有的基于扩散方程-简化球谐近似构建的混合光传输模型，求解过程较为复杂，计算量大。光传输模型计算精度和效率的缺陷，制约了光学分子断层成像中光源的快速精确重建。因此，目前需要一种新的光传输模型来实现计算精度和效率之间的平衡。

发明内容

针对现有技术中的缺陷和不足，本发明提供了一种基于网格重组的混合光传输建模方法，解决现有光传输模型模拟光在非均匀介质中的传输时，计算精度和效率之间难以达到有效平衡的问题。

为达到上述目的，本发明采取如下的技术方案：

一种基于网格重组的混合光传输建模方法，该方法根据动物组织的光学特性，将动物各组织划分为高散射组织区域和低散射组织区域；根据各组织所属区域将动物模型离散得到由空间位置排序的节点和四面体重新排序，得到两个独立网格；在这两个网格上，基于光传输理论，高散射组织区域采用扩散方程DE近似建模，低散射组织区域采用三阶简化球谐SP₃近似建模，两个区域有相应的相关系统矩阵；通过耦合这两个系统矩阵，把两个网格合成一个重组网格，建立混合光传输模型。

本发明还包括如下技术特征：

具体的，包括以下步骤：

步骤一，获取动物模型并进行预处理，得到动物模型各组织和初始离散网格；

步骤二，根据动物模型各组织的光学参数，并将各组织分为高散射组织区域和低散射组织区域，其中高散射组织区域为Region1，低散射组织区域为Region2；

步骤三，将步骤一得到的初始离散网格中的节点和四面体重新排序，根据节点和四面体所属区域，属于Region1的节点和四面体组成网格1，其中包含的节点数表示为N₁；属于Region2的节点和四面体组成网格2，其中包含的节点数表示为N₂；

步骤四，在步骤三得到的网格1上，其所属区域作为高散射组织区域采用DE近似建模得到对应的系统矩阵M₁；在网格2上，其所属区域作为低散射组织区域采用SP₃近似建模得到对应的系统矩阵M₂；

步骤五，步骤四得到的系统矩阵M₁，维数为N₁×N₁，N₁为网格1包含的节点数；系统矩阵M₂，维数为2N₂×2N₂，N₂为网格2包含的节点数；

步骤六，将系统矩阵M₁和系统矩阵M₂联结组成混合系统矩阵M’，并对重复项耦合后得到耦合后的混合系统矩阵M；

步骤七，把网格1和网格2合成一个整体的重组网格，建立重组网格表面光子通量密度和光源功率密度的关系：

其中，M为步骤六获得的耦合后的混合系统矩阵，Φ为表面光子通量密度，S为光源功率密度，F是源权重矩阵；Φ₁、F₁和S₁对应于网格1；Φ₂、F₂和S₂对应于网格2；合成一个整体的重组网格后，上式表示为：

MΦ＝FS

经变换为

Φ＝M^-1FS

得到混合光传输模型；经过运算可由已知光源信息，获得动物模型表面光子通量密度分布；或已知动物模型表面光子通量密度分布信息，逆向重建出光源信息。

具体的，所述步骤一包括：

步骤1.1，由试验动物的CT成像切片数据中提取组织信息，获得动物模型，该动物模型包括肌肉、肺、心脏、肝脏、胃和肾脏各组织和光源；

步骤1.2，采用Amira软件对动物模型进行剖分，将动物模型离散获得初始离散网格，其中包括N个网格节点、X个内部四面体，该初始离散网格上，节点和四面体是按照其空间位置进行排列。

具体的，所述步骤二包括：

步骤2.1，根据光源发射波长获得各组织对应的光学特性参数，所述光学特性参数包括吸收系数μ_a，散射系数μ_s，各向异性系数g和相对折射率n，通过μ_s′＝μ_s(1-g)计算得到约化散射系数μ_s′；

步骤2.2，根据光学特性参数将各组织分组为两个区域：高散射组织区域Region1和低散射组织区域Region2，其中，当约化散射系数与吸收系数的比值大于给定阈值时，将该组织作为高散射组织区域；反之，该组织将作为低散射组织区域。

具体的，所述步骤六中，将系统矩阵M₁和系统矩阵M₂组成一个混合系统矩阵M’：

其中，由于网格1和网格2在边界上相连，边界上的光子通量密度是连续的，即

Φ₁′＝Φ₂′

Φ₁’为Region1边界上的光子通量密度，Φ₂’为Region2边界上的光子通量密度；

边界上的节点同时属于Region1和Region2，对应于空间同一点，在构成混合系统矩阵时产生了重复项，为避免边界节点的重复贡献，需对其进行耦合，包括以下步骤：

步骤6.1，首先找出网格1和网格2中，空间位置相同的节点，即边界上的重复节点；

步骤6.2，由上述重复节点，获得混合系统矩阵M’中来自于M₁和M₂部分对应的重复项，即M₁部分的一行对应M₂部分的两行；

步骤6.3，对重复项进行操作，将M₂中重复内容移动到M₁对应的行，然后将M₂对应项设为0；将整行的项设置为0时，这些共有节点与其自身的关系也被删除，为了恢复这种关系，将上述被设置为0的项中属于系统矩阵M’中对角线上的项设为1；在混合系统矩阵M’中，每个共有节点与其他共有节点的关系表示为：

其中，a对应Region1部分重复节点的项，b和c对应Region2部分重复节点的项；

通过上式对混合系统矩阵M’中对应项进行操作，完成对两个区域连接边界上重复节点的耦合，此过程中，M’中对应项发生变化，可由分块矩阵M₁₁、M₁₂、M₂₁、M₂₂表示，重复项耦合后得到耦合后的混合系统矩阵M，整个过程表示为：

其中，M₁为Region1对应的系统矩阵，M₂为Region2对应的系统矩阵，M’为混合系统矩阵，M为耦合后的混合系统矩阵；M₁₁、M₁₂、M₂₁、M₂₂为分块矩阵。

本发明与现有技术相比，有益的技术效果是：

本发明通过网格重组和系统矩阵的联结，使混合光传输模型能够一步求解，联结后混合系统矩阵M的维数为N₁+2N₂；与传统混合模型中，不同的近似方程需要分步求解，混合系统矩阵的维数为有限元体网格总节点数的两倍相比，本发明求解步骤简单，得到的混合系统矩阵维数较低，运算量相对较小，计算效率得到提高。

根据高散射组织区域和低散射组织区域分别采用DE近似模型和SP₃近似模型；为每个组织选择合适的光传输模型，提高了计算精度。

附图说明

图1为本发明总体流程图；

图2为网格重组过程的示意图；

图3为一个高度简化的包含2对共有节点的混合网格例子；

图4各模型仿真的表面光子通量密度通过Tecplot软件的展示图；其中图(a)为分子光学模拟环境(MOSE)平台仿真的结果；图(b)为DE近似模型的运行结果；图(c)为SP₃近似模型的运行结果；图(d)为本发明提出的混合模型的运行结果。

图5数字鼠实验中各个光传输模型获得的表面光子通量密度比较结果；其中图(a)为各采样点对应光子通量密度的比较；图(b)为各模型对应的平均光子通量密度误差的比较。

具体实施方式

以下对本发明的具体实施方式进行详细说明。应当理解的是，此处所描述的具体实施方式仅用于说明和解释本发明，并不用于限制本发明。

本发明提出了一种基于网格重组的混合光传输建模方法，该方法根据动物组织的光学特性，将动物各组织划分为高散射组织区域和低散射组织区域；根据各组织所属区域将动物模型离散得到由空间位置排序的节点和四面体重新排序，得到两个独立网格；在这两个网格上，基于光传输理论，高散射组织区域采用扩散方程DE近似建模，低散射组织区域采用三阶简化球谐SP₃近似建模，两个区域有相应的相关系统矩阵；通过耦合这两个系统矩阵，把两个网格合成一个重组网格，建立混合光传输模型。

包括以下步骤：

步骤一，获取动物模型并进行预处理，得到动物模型各组织和初始离散网格；

步骤1.1，由试验动物的CT成像切片数据中提取组织信息，获得动物模型，该动物模型包括肌肉、肺、心脏、肝脏、胃和肾脏各组织和光源；

步骤1.2，采用Amira软件对动物模型进行剖分，将动物模型离散获得初始离散网格，其中包括N个网格节点、X个内部四面体，该初始离散网格上，节点和四面体是按照其空间位置进行排列。

步骤二，根据动物模型各组织的光学参数，并将各组织分为高散射组织区域和低散射组织区域，其中高散射组织区域为Region1，低散射组织区域为Region2；

步骤2.1，根据光源发射波长获得各组织对应的光学特性参数，所述光学特性参数包括吸收系数μ_a，散射系数μ_s，各向异性系数g和相对折射率n，通过μ_s′＝μ_s(1-g)计算得到约化散射系数μ_s′；

步骤2.2，根据光学特性参数将各组织分组为两个区域：高散射组织区域Region1和低散射组织区域Region2，其中，当约化散射系数与吸收系数的比值大于给定阈值时，将该组织作为高散射组织区域；反之，该组织将作为低散射组织区域。

步骤三，将步骤一得到的初始离散网格中的节点和四面体重新排序，根据节点和四面体所属区域，属于Region1的节点和四面体组成网格1，其中包含的节点数表示为N₁；属于Region2的节点和四面体组成网格2，其中包含的节点数表示为N₂；

步骤四，在步骤三得到的网格1上，其所属区域作为高散射组织区域采用DE近似建模；在网格2上，其所属区域作为低散射组织区域采用SP₃近似建模；

步骤四中，DE近似建模为：

其中，

为梯度，Φ(r)是光密度，S(r)是光源的功率密度，μ_a为吸收系数，μ_s为散射系数，g为各向异性系数，An为边界

处的折射失配系数；ν(r)是边界上的外向单位法向量；

所述SP₃近似建模为：

其中

为光子通量合成矩，μ_ai＝μ_a+μ_s(1-gⁱ)(i＝1,2,3)为第i阶吸收系数，S(r)是光源的功率密度，μ_a为吸收系数，μ_s为散射系数，g为各向异性系数，ν(r)是边界上的外向单位法向量，A_i,B_i,C_i,D_i(i＝1,2),J_i(i＝0,1,2,3)都为系数；

步骤五，由DE方程，获得Region1对应的系统矩阵M₁，其维数为N₁×N₁；由SP₃方程，获得Region2对应的系统矩阵M₂，其维数为2N₂×2N₂；

步骤六，将系统矩阵M₁和系统矩阵M₂联结组成混合系统矩阵M’，并对重复项耦合后得到耦合后的混合系统矩阵M；

具体的，步骤六中，将系统矩阵M₁和系统矩阵M₂组成混合系统矩阵M’：

其中，由于网格1和网格2在边界上相连，边界上的光子通量密度是连续的，即

Φ₁′＝Φ₂′

Φ₁’为Region1边界上的光子通量密度，Φ₂’为Region2边界上的光子通量密度；

边界上的节点同时属于Region1和Region2，对应于空间同一点，在构成混合系统矩阵时产生了重复项，为避免边界节点的重复贡献，需对其进行耦合，包括以下步骤：

步骤6.1，首先找出网格1和网格2中，空间位置相同的节点，即边界上的重复节点；

步骤6.2，由上述重复节点，获得系统矩阵M’中来自于M₁和M₂部分对应的重复项，即M₁部分的一行对应M₂部分的两行；

步骤6.3，对重复项进行操作，将M₂中重复内容移动到M₁对应的行，然后将M₂对应项设为0；将整行的项设置为0时，这些共有节点与其自身的关系也被删除，为了恢复这种关系，将上述被设置为0的项中属于混合系统矩阵M’中对角线上的项设为1；在混合系统矩阵M’中，每个共有节点与其他共有节点的关系表示为：

其中，a对应Region1部分重复节点的项，b和c对应Region2部分重复节点的项；

通过上式对混合系统矩阵M’中对应项进行操作，完成对两个区域连接边界上重复节点的耦合，此过程中，M’中对应项发生变化，可由分块矩阵M₁₁、M₁₂、M₂₁、M₂₂表示，重复项耦合后得到耦合后的混合系统矩阵M，整个过程表示为：

其中，M₁为Region1对应的系统矩阵，M₂为Region2对应的系统矩阵，M’为混合系统矩阵，M为耦合后的混合系统矩阵；M₁₁、M₁₂、M₂₁、M₂₂为分块矩阵。

步骤七，把网格1和网格2合成一个整体的重组网格，建立重组网格表面光子通量密度和光源功率密度的关系：

其中，M为步骤六获得的混合系统矩阵，Φ为表面光子通量密度，S为光源功率密度，F是源权重矩阵；Φ₁、F₁和S₁对应于网格1；Φ₂、F₂和S₂对应于网格2；合成一个整体的重组网格后，上式表示为：

MΦ＝FS

经变换为

Φ＝M^-1FS

得到混合光传输模型；经过运算可由已知光源信息，获得动物模型表面光子通量密度分布；或已知动物模型表面光子通量密度分布信息，逆向重建出光源信息。

实施例1：

本实施例采用常用的数字小鼠模型，仅选取小鼠躯干部分，高度35mm，作为研究区域。选用波长为610nm作为本发明中生物发光断层成像的发射波长。光源设置为一个半径为0.8mm、高为1.6mm的圆柱形，放置在肝脏中，中心点坐标(19mm；7mm；17.5mm)，光源强度设为1000nw。

将小鼠躯干模型肌肉、心脏、肝脏、肺、胃、肾脏、光源等7个组织，离散为12967个网格节点、67102个内部四面体。并根据光学特性参数分组为高散射特性区域和低散射特性区域，其中肝脏和光源作为低散射特性区域(Region2)，其他组织作为高散射特性区域(Region1)，如附图1所示。

Region1对应的网格1包含11665个网格节点、55275个内部四面体，作为高散射特性区域采用扩散方程近似(DE)建模，系统矩阵M₁的维数为11665×11665。Region2对应的网格2包含2776个网格节点、11827个内部四面体。作为低散射特性区域采用三阶简化球谐近似(SP₃)建模，系统矩阵M₂的维数为5552×5552。通过对节点空间位置的比较得到两个网格中包含1474个重复节点，对重复节点耦合，获得混合系统矩阵M，维数为17217×17217。把两个网格合成一个重组网格，获得混合光传输模型。

图2为网格重组过程的示意图；其中图(a)为动物模型离散得到的初始网格，包含18个节点，按照其空间位置进行排列；图(b)和图(c)分别为重组后的网格1和网格2，节点根据所属区域重新排序，网格1包含14个节点，网格2包含10个节点。在网格1和网格2的边界上，同属于两个网格的为6个重复节点。

图3为一个高度简化的包含2对共有节点的混合网格例子；进一步解释混合系统矩阵中重复节点对应项的耦合过程：以一组基本网格为例，节点3、4与节点5、6为重复节点；所得混合系统矩阵，用对应的节点编号标记每一行，如附图1中的界面耦合过程；例如，节点3和节点5对应于空间中的同一点，应当被耦合；将混合系统矩阵中节点5的项移动到第3行，然后将节点5中的每一项设为零；为了恢复共有节点与其自身的关系，将节点5中被设为零的项属于混合系统矩阵对角线元素的项设置为1；为了恢复每个共有节点与其他共有节点的关系，节点5的项中连接b和a的元素设置为1，连接b和c的元素设置为2/3。

以分子光学模拟环境(MOSE)平台的仿真结果作为标准，采用单一的DE近似模型和SP₃近似模型与本发明提出的基于网格重组的混合光传输模型(HSDM)进行对比。分别采用不同近似建立的模型进行前向仿真，附图4为采用MOSE、DE近似模型、SP₃近似模型和HSDM计算的表面光子通量通过Tecplot软件的展示图，其中图(a)为分子光学模拟环境(MOSE)平台仿真的结果；图(b)为DE近似模型的运行结果；图(c)为SP₃近似模型的运行结果；图(d)为本发明提出的混合模型的运行结果。由于数值方法计算的结果和MOSE平台仿真的结果单位量纲不同，为了更明确地说明比较，所有结果都进行归一化处理。

为了进一步评估实验结果，附图5展示了实验中各个光传输模型获得的表面光子通量密度在x＝19mm处各采样点的结果比较，并计算DE近似模型、SP₃近似模型和HSDM的计算结果与MOSE仿真结果的平均相对误差ARE：

其中Simulation_i为各模型通过MATLAB仿真得到的表面光子通量密度，Mose_i为通过MOSE平台得到的表面光子通量密度，N为采样点总数。

图5(a)为各采样点对应光子通量密度的比较；图(b)为各模型对应的平均光子通量密度误差的比较，从附图5可以明显看出，与DE近似模型和SP₃近似模型相比，HSDM的结果更接近于MOSE的结果。DE、SP₃和HSDM模型的ARE分别为0.3011，0.2479和0.1606，HSDM的平均相对误差最小。

本实施例中，数字小鼠实验得到了表面光子通量密度，同时记录了各模型运行消耗的时间。DE近似模型平均用时103.57s，SP₃近似模型平均用时561.89s，HSDM平均用时247.81s，计算比较得HSDM比SP₃近似模型用时减少55.9％。

实验结果表明，与DE近似模型和SP₃近似模型相比，本发明提出的基于网格重组的混合光传输模型在精度和效率之间取得了更好的平衡，相对更适合作为光学分子断层成像的光传输模型。

Claims (5)

1.一种基于网格重组的混合光传输建模方法，其特征在于，该方法根据动物组织的光学特性，将动物各组织划分为高散射组织区域和低散射组织区域；根据各组织所属区域将动物模型离散得到由空间位置排序的节点和四面体重新排序，得到两个独立网格；在这两个网格上，基于光传输理论，高散射组织区域采用扩散方程DE近似建模，低散射组织区域采用三阶简化球谐SP₃近似建模，两个区域有相应的相关系统矩阵；通过耦合这两个系统矩阵，把两个网格合成一个重组网格，建立混合光传输模型。

2.如权利要求1所述的基于网格重组的混合光传输建模方法，其特征在于，包括以下步骤：

步骤一，获取动物模型并进行预处理，得到动物模型各组织和初始离散网格；

步骤二，根据动物模型各组织的光学参数，并将各组织分为高散射组织区域和低散射组织区域，其中高散射组织区域为Region1，低散射组织区域为Region2；

步骤三，将步骤一得到的初始离散网格中的节点和四面体重新排序，根据节点和四面体所属区域，属于Region1的节点和四面体组成网格1，其中包含的节点数表示为N₁；属于Region2的节点和四面体组成网格2，其中包含的节点数表示为N₂；

步骤四，在步骤三得到的网格1上，其所属区域作为高散射组织区域采用DE近似建模得到对应的系统矩阵M₁；在网格2上，其所属区域作为低散射组织区域采用SP₃近似建模得到对应的系统矩阵M₂；

步骤五，步骤四得到的系统矩阵M₁，维数为N₁×N₁，N₁为网格1包含的节点数；系统矩阵M₂，维数为2N₂×2N₂，N₂为网格2包含的节点数；

步骤六，将系统矩阵M₁和系统矩阵M₂联结组成混合系统矩阵M’，并对重复项耦合后得到耦合后的混合系统矩阵M；

步骤七，把网格1和网格2合成一个整体的重组网格，建立重组网格表面光子通量密度和光源功率密度的关系：

其中，M为步骤六获得的耦合后的混合系统矩阵，Φ为表面光子通量密度，S为光源功率密度，F是源权重矩阵；Φ₁、F₁和S₁对应于网格1；Φ₂、F₂和S₂对应于网格2；合成一个整体的重组网格后，上式表示为：

MΦ＝FS

经变换为

Φ＝M^-1FS

得到混合光传输模型；经过运算可由已知光源信息，获得动物模型表面光子通量密度分布；或已知动物模型表面光子通量密度分布信息，逆向重建出光源信息。

3.如权利要求2所述的基于网格重组的混合光传输建模方法，其特征在于，所述步骤一包括：

步骤1.1，由试验动物的CT成像切片数据中提取组织信息，获得动物模型，该动物模型包括肌肉、肺、心脏、肝脏、胃和肾脏各组织和光源；

步骤1.2，采用Amira软件对动物模型进行剖分，将动物模型离散获得初始离散网格，其中包括N个网格节点、X个内部四面体，该初始离散网格上，节点和四面体是按照其空间位置进行排列。

4.如权利要求2所述的基于网格重组的混合光传输建模方法，其特征在于，所述步骤二包括：

步骤2.1，根据光源发射波长获得各组织对应的光学特性参数，所述光学特性参数包括吸收系数μ_a，散射系数μ_s，各向异性系数g和相对折射率n，通过μ_s′＝μ_s(1-g)计算得到约化散射系数μ_s′；

步骤2.2，根据光学特性参数将各组织分组为两个区域：高散射组织区域Region1和低散射组织区域Region2，其中，当约化散射系数与吸收系数的比值大于给定阈值时，将该组织作为高散射组织区域；反之，该组织将作为低散射组织区域。

5.如权利要求2所述的基于网格重组的混合光传输建模方法，其特征在于，所述步骤六中，将系统矩阵M₁和系统矩阵M₂组成一个混合系统矩阵M’：

其中，由于网格1和网格2在边界上相连，边界上的光子通量密度是连续的，即

Φ₁’＝Φ₂’

Φ₁’为Region1边界上的光子通量密度，Φ₂’为Region2边界上的光子通量密度；

边界上的节点同时属于Region1和Region2，对应于空间同一点，在构成混合系统矩阵时产生了重复项，为避免边界节点的重复贡献，需对其进行耦合，包括以下步骤：

步骤6.1，首先找出网格1和网格2中，空间位置相同的节点，即边界上的重复节点；

步骤6.2，由上述重复节点，获得混合系统矩阵M’中来自于M₁和M₂部分对应的重复项，即M₁部分的一行对应M₂部分的两行；

步骤6.3，对重复项进行操作，将M₂中重复内容移动到M₁对应的行，然后将M₂对应项设为0；将整行的项设置为0时，这些共有节点与其自身的关系也被删除，为了恢复这种关系，将上述被设置为0的项中属于系统矩阵M’中对角线上的项设为1；在混合系统矩阵M’中，每个共有节点与其他共有节点的关系表示为：

其中，a对应Region1部分重复节点的项，b和c对应Region2部分重复节点的项；

通过上式对混合系统矩阵M’中对应项进行操作，完成对两个区域连接边界上重复节点的耦合，此过程中，M’中对应项发生变化，可由分块矩阵M₁₁、M₁₂、M₂₁、M₂₂表示，重复项耦合后得到耦合后的混合系统矩阵M，整个过程表示为：

其中，M₁为Region1对应的系统矩阵，M₂为Region2对应的系统矩阵，M’为混合系统矩阵，M为耦合后的混合系统矩阵；M₁₁、M₁₂、M₂₁、M₂₂为分块矩阵。

Images