佛山网站建设3lue,常州做网站哪家好,嘉兴市城市建设门户网站,鹤壁公司做网站我的主页#xff1a;
技术邻#xff1a;小铭的ABAQUS学习的技术邻主页博客园 : HF_SO4的主页哔哩哔哩#xff1a;小铭的ABAQUS学习的个人空间csdn#xff1a;qgm1702
博客园文章链接#xff1a;
https://www.cnblogs.com/aksoam/p/17590039.html 分别是向量判别法…我的主页
技术邻小铭的ABAQUS学习的技术邻主页博客园 : HF_SO4的主页哔哩哔哩小铭的ABAQUS学习的个人空间csdnqgm1702
博客园文章链接
https://www.cnblogs.com/aksoam/p/17590039.html 分别是向量判别法算法来自他人论文、体积判别法code 是我从网上找的。 方法一 向量判别法
方法来自一会议论文《判断点与多面体空间位置关系的一个新算法_石露》2008年,知网、万方、百度学术有收录。 优点
适合任意多面体计算简单速度快
算法原理 Matlab 实现
主函数为InPolyhedronByVJM(Nodes,P),当前仅支持ABAQUS 四面体和五面体单元其余有需要的可自行在switch语句添加函数。 输入
Nodes包含节点顺序和坐标的N x 3的矩阵P1x3的行向量
这是我根据ABAQUS单元规则编写的函数因此其他软件的网格信息文件还需要重写。
function flag InPolyhedronByVJM(Nodes,P)
% 根据向量判别法判断点和多面体的位置关系.
% 算法依据石露,白冰,李小春. 判断点与多面体空问位置关系的一个新算法[C].
% 输入
% Nodes: n x 3 matrix
% P : 1 x 3 row vector
% 单元节点顺序需要满足ABAQUS的约定。
% Nodes是一个n x 3的数值矩阵每一行表示多面体的一个节点的空间坐标。
% 点在多面体内部的充分必要条件每个face上任一点到点P的向量和该face法向量的数量积都
% 小于等于0否则在体外。
NumOfNode size(Nodes, 1);
switch NumOfNodecase 4flag InTetrahedron(Nodes, P);disp(tetra)case 5flag InPyramid(Nodes, P);disp(parymid)
end
end针对C3D4类单元的判别函数
function flag InTetrahedron(Nodes, P)
% tetra elem have 4 face.
% ABAQUS rule about node ordering and face numbering on element
FaceIDX [1 2 3;1 4 2;2 4 3;3 4 1];
for faceId 1 : 1 : size(FaceIDX, 1)% judge face i: node 1-node 2-node 3 face% get face normal vector(outside)n -1 .* GetNormVector(Nodes(FaceIDX(faceId, 1), :), Nodes(FaceIDX(faceId, 2), :), Nodes(FaceIDX(faceId, 3), :));% calculate dot product of P and normal vector nN1P P - Nodes(FaceIDX(faceId, 1), :);f dot(N1P, n);if f 0flag 0;returnend
end
flag1;
end针对C3D5类单元的判别函数
function flag InPyramid(Nodes, P)
% Pyramid elem have 5 face.
% ABAQUS rule about node ordering and face numbering on element
FaceIDX [1 2 3 4;1 5 2 0;2 5 3 0;3 5 4 0;4 5 1 0];
for faceId 1 : 1 : size(FaceIDX, 1)% judge face i: node 1-node 2-node 3 face% get face normal vector(outside)n -1 .* GetNormVector(Nodes(FaceIDX(faceId, 1), :), Nodes(FaceIDX(faceId, 2), :), Nodes(FaceIDX(faceId, 3), :));% calculate dot product of P and normal vector nN1P P - Nodes(FaceIDX(faceId, 1), :);f dot(N1P, n);if f 0flag 0;returnend
end
flag 1;
end求面法向向量的函数
function NormalVector GetNormVector(p1, p2, p3)
% function return a Normal Vector,base on RightHand Rule, according to three point (row vector)
% NormalVectorcross product of (p2-p1) and (p3-p1)
% check
if (~isrow(p1)) || (~isrow(p2)) || (~isrow(p3))return
end
p1p2 p2 - p1;
p1p3 p3 - p1;
NormalVector cross(p1p2, p1p3);
end附abqus四面体和五面体单元的节点约定 方法二体积判别法 这个代码是我在网上找的出处已经忘了。 Matlab实现(仅限四面体)
function inflag inpolyhedronByVolCal(points_mat, p_detected)
% input:
% points_set : 4 points x y z coordinate,integrated in a 4X3 matrix
% p_detected : a point needed to detect ,1X3 matrix contain x y z
% coordinate
% return : inflag: 0 or 1
D0 det([points_mat(1, :) 1;points_mat(2, :) 1;points_mat(3, :) 1;points_mat(4, :) 1]);
D1 det([p_detected 1;points_mat(2, :) 1;points_mat(3, :) 1;points_mat(4, :) 1]);
D2 det([points_mat(1, :) 1;p_detected 1;points_mat(3, :) 1;points_mat(4, :) 1]);
D3 det([points_mat(1, :) 1;points_mat(2, :) 1;p_detected 1;points_mat(4, :) 1]);
D4 det([points_mat(1, :) 1;points_mat(2, :) 1;points_mat(3, :) 1;p_detected 1]);if (D0 0 D1 0 D2 0 D3 0 D4 0) || ((D0 0 D1 0 D2 0 D3 0 D4 0))inflag 1;
elseinflag 0;
end
end