[发明专利]一种三维弹性模量成像方法

权利要求书

1.一种三维弹性模量成像方法，其特征在于，包括以下步骤：

(A)建立对象初始状态的三维有限元模型；

(B)向对象施加系统载荷；

(C)测量对象的有限元节点在施加系统载荷后的三维运动信息，所 述的三维运动信息包括位移、速度；

(D)使用基于状态空间的H∞滤波估计方法，得到如下连续状态估 计方程(20)～(23)，利用方程(20)～(23)估计出对象弹性模量的三 维分布，实现对象的三维弹性模量成像；

x·^θ·^=FAC(x(t))00x^(t)θ^+BCw(t)0+Σ(t)-1DT(y-Dx^(t))---(20)]]>

Σ·=-Σ(t)FAC(x(t))00-FT0ACT(x(t))0Σ(t)+I00-γ-2Q-Σ(t)I000Σ(t)---(21)]]>

x^(t+Δt)θ^=x^(t)θ^+x·^θ·^Δt---(22)]]>

Σ(t+Δt)=Σ(t)+Σ·Δt---(23)]]>

其中：

F=0I0-αI;]]>

AC=00-M-1G1-βM-1G2;]]>

BC=000M-1;]]>

x^(t)=U^(t)U·^T;]]>

U为对象的有限元节点的位移向量，其一阶微分和二阶微分分别 表示速度向量和加速度向量；

I为与U维数相应的单位矩阵；

α和β为瑞利阻尼常数，可取0.01～0.05；

M为对象的有限元质量矩阵，是对象材料密度的函数；

G₁和G₂是计算过程中的过渡矩阵；

t为时间变量；

Δt为估计时间间隔；

为对象的弹性模量向量估计值，θ初值根据对象的经验知识预先设 定，要求与对象的真实弹性模量接近；

w(t)为控制向量，w(t)＝[0 R]^T，R为系统载荷向量；

y(t)为步骤(C)中测量对象得到的三维运动信息；

D为测量矩阵，D的行数为三维运动信息可以测量得到的有限元节点 数，D的列数为x(t)的维数，D的每行中，索引位置对应的有限元节点的 三维运动信息如果可以测量得到，该索引位置的元素设为1，否则元素设 为0；

为x(t)的导数，表示在具体的t时刻的x(t)的变化量；

为θ的修正值；

γ为估计的质量的约束值，影响估计值精确度，取1；

Σ=Σ1Σ2Σ2TΣ3,]]>其中∑₁、∑₂、∑₃的初值∑₁(0)、∑₂(0)、∑₃(0)分别为

∑₁(0)＝Q₁、∑₂(0)＝0、∑₃(0)＝Q₀；

为∑的导数，表示在具体的t时刻的∑的变化量；

Q、Q₀、Q₁为权重矩阵；

所述的过渡矩阵G₁和G₂的构造方法如下：

(a)初始化两个N×N_e的空矩阵G₁和G₂，N为系统节点变量数，在 三维空间中，每个节点有三个方向的自由度，即系统总节点数×3，N_e为系 统总单元数；初始化一个N×N的空矩阵K_g′，作为不含弹性模量E的系统 临时刚度矩阵；

(b)对于一个编号为j的单元，构造不包含该单元弹性模量E_j的刚 度矩阵K_j′；

(c)按照系统生成网格时单元和节点的编号规则，把K_j′中的元素插 入系统临时刚度矩阵K_g′中相应的位置中，并将其单元内下标改成相应的 系统下标；

(d)K_g′与U相乘得到的向量，即是矩阵G₁的第j列；K_g′与相乘 得到的向量，即是矩阵G₂的第j列；

(e)返回到第(b)步直到所有节点遍历，完成G1和G2的构造。