[发明专利]一种求解薄基片应力的方法

权利要求书

1.一种求解薄基片应力的方法，其特征在于，所述方法包括如下步骤：

S1：获取待分析应力的薄基片的基础参数；所述基础参数包括直径、厚度、弹性模量、剪切模块、泊松比和密度；

S2：建立与所述待分析应力的薄基片具有相同基础参数的有限元模型；其中，所述有限元模型的构建方法包括如下步骤：

S21：定义材料属性和单元类型；详细过程如下：

首先定义单元编号选择壳单元的类型，然后依次设置第一材料层和第二材料层的材料的密度、弹性模型、泊松比和剪切模量；最后依次定义第一材料层和第二材料层的厚度；

S22：建立薄基片的有限元模型；详细过程如下：

首先将模型中的坐标系设置为柱面坐标系；在模型中长度等于1/3圆周的弧长的扇形区域内生成关键点，将生成的内部的关键点用直线连接，外部的关键点用圆弧线连接，使得扇形区域被分为多个四边形；然后选择圆周内的所有线，在每条线上生成等分关键点，将四边形进一步划分成更小的四边形；接着按照同样的方法将模型中其余2/3部分也生成关键点并进行连线，合并重复关键点；最后在模型中选择边缘的圆弧线生成圆面；

S23：对建立的有限元模型划分网格；详细过程如下：

首先将模型的网格类型选择为四边形网格划分，指定网格划分方法为映射网格划分；然后选择模型中除圆周以外的所有线，将选择的线组成一个线组，用线组分割上步骤生成的圆面，并将分割的面粘接在一起；接着删除之前的线组，选择模型中的所有线，重新将所选择的线组成一个新的线组；最后通过新的线组来控制网格密度，进行网格划分、生成单元和节点；

S24：为有限元模型中的节点添加约束；详细过程如下：

选择三个节点施加约束，对于第一个节点限制X、Y、Z方向的移动；对于第二个节点限制X，Z方向的移动；对第三个节点限制Z方向的移动；

S3：依次向有限元模型中的每个面施加单位应力载荷σ_u，然后对有限元模型进行求解；定义输出路径，并将模型中施加单位应力载荷的各个面的求解结果映射到路径上；输出模型中所有面的求解结果；具体过程如下：

S31：获取模型的总面数；

S32：在第一材料层上施加单位应力载荷，选择编号最小的面相关的壳单元，在壳单元上定义单位应力载荷的数值；

S33：对模型进行求解；

S34：读取求解的结果，并创建一个柱面坐标系，在柱面坐标系下创建一个圆形路径，设置映射项目，定相邻两点间分数段；

S35：定义一个数组，将分析计算结果中Z方向的位移映射到创建的路径上，获取路径中的数据并存放在数组中；

S36：在模型中从内到外依次建立多条路径，使得采样点分布均匀；同时将每条路径对应的数组赋值给一个新的数组，使得新的数组包含之前所有采样点的位置和位移信息；

S37：定义一个新的数组，将步骤S36中数组的采样点的X坐标、Y坐标和Z方向的位移，分别赋值到新定义的数组中，并对该数组中的数据进行输出；

S38：重复步骤S31-S37依次将单位应力载荷加载到模型中的所有面上，其中每次加载载荷时仅在模型其中一个面上加载，加载前删除前一个面上加载的单位应力载荷，依次获得每个加载的面中的数组的输出结果；

S4：根据有限元模型的求解结果，获取模型中施加单位应力载荷的面与其相邻面，进而构造包含各面之间中心点距离关系的矩阵L；具体过程如下：

S41：定义一个新的空数组，选择所有模型中的所有实体，获取模型的总面数；

S42：创建一个新的二维数组，将二维数组的行数设为步骤S41中的获取的模型的总面数，列数为4，二维数组的第一列为1至总面数；

S43：选择模型中最小编号的面，计算出最小编号面的几何数据，获取最小编号面中心点的X、Y、Z坐标，将坐标值赋值各二维数组的第2、3、4列；依次将所有施加载荷的面的中心点的X、Y、Z的坐标赋值给二维数组的各行；最后输出二维数组中的数据；

S44：对上步骤的输出结果进行数值分析，基于分析结果形成一个包含面的编号和面的中心点坐标的矩阵，获取矩阵行数；

S45：创建一个全0矩阵，所述矩阵的行数为模型的总面数，列数为5；

S46：把步骤S44中的矩阵的第一列元素赋值给步骤S45中创建的矩阵的第一列；并依次通过步骤S44中计算矩阵对应的模型的第一个面的中心点到所有面中心点的距离，将计算结果赋值给步骤S44中矩阵的第四列；

S47：对步骤S46中经过赋值的矩阵按第四列升序排列，然后把排列后的矩阵赋值各临时变量，形成一个新的矩阵；

S48：把步骤S47中创建的新矩阵的第2、3、4、5行的第一列元素分别赋值给步骤S45中矩阵的第一行的2、3、4、5列；并依次按照S46中的方法，计算剩余的面的中心点到所有面中心点的距离，完成步骤S45中矩阵所有行的赋值操作；

S49：首先，创建一个矩阵L，矩阵L为单位阵，其维数等于施加单位应力载荷的面数，主对角元素为1，其余元素为0；

其次，将创建的矩阵L中主对角线上元素全修改为4；

然后，对上步骤的矩阵L的第1行进行赋值：以步骤S48中处理后的矩阵第1行的2、3、4、5列元素值为索引，对矩阵L第1行中的相应列的元素均赋值为-1；

最后，按照相同的方式，依次完成矩阵L中所有行的赋值操作；

S5：对薄基片进行测量获取真实薄基片的面形，由多个真实薄基片测量点的位移组成的列向量构建矩阵b，由多个测量点分别在各个施加单位应力载荷的面作用下的位移组成矩阵记作W；具体过程如下：

S51：创建一个新矩阵，所述新矩阵的行数为薄基片面形获取过程中测量点的个数；列数为有限元模型中施加载荷的面的总数，依次将步骤S38中输出的结果中第3列数据按列导入到创建的新矩阵中，完成W矩阵的构造；

S52：使用传感器获得真实薄基片的面形，将测量坐标系与有限元坐标系相重合，使得测量坐标系中的点与有限元对应约束点的Z向位移分量为0，通过面形插补获得与有限元模型位移提取点相匹配的真实薄基片的测量点，并获取测量点的X，Y，Z坐标，从而得到b阵；

S6：采用正则化方法求解线性方程组，使其满足最小；绘制不同μ值时以正则化方法惩罚项为横坐标，以为纵坐标的取舍曲线；μ表示正则化方法惩罚项的系数，x为待求向量，是由α_j组成一个n维列向量；具体过程如下：

S61：定义变量μ为正则化方法惩罚项的系数，变量变化范围大于0；当给定一个μ，x满足下式：

x＝(W^TW+μL)^-1W^Tb，

上式中，^T表示矩阵的转置；W表示步骤S51中构造的矩阵；b表示步骤S52中得到的矩阵；

接着以和μ相对应的为横坐标，为纵坐标；绘制μ在大于0内散点拟合的取舍曲线；

S62：选择取舍曲线转折处的点所对应的μ，计算x，根据x计算的面形与测量面形的差值以及应力分布选取最后的μ值；

最后，根据如下公式计算薄基片的实际应力分布：

σ_j＝α_jσ_u,

上式中，σ_j表示薄基片第j个面的实际应力；σ_u表示有限元模型中各个面的单位应力载荷值。