Plot lines from FEM (PDE) stress vectors
5 views (last 30 days)
I would like to calculate the vectors of the principal stresses in each node of the mesh after an FEM analysis with Matlab (with the Partial Differential Equation Toolbox) of an imported 3D component (.stl) and plot them in my component. This also works very well, just takes a relatively long time because I do this via a loop (surely it can be done better). I calculate the vectors from the double lists of the different stress components generated during the FEM analysis.
Now, at each node (xyz coordinates) in the component, the direction of the stress is available in the form of a vector (u, v, w) - plotted with quiver3.
%Here would be the code for the FEM analysis
%As relevant results, I get several list (doubles) with 1x x-thousand entries about the stress components at each node (in this example, the first 500 in the list)
for i = 1 : 500 % for i continuous from 1 to X (here just 500 because 5000 takes to long)
S = [qx(i) txy(i) txz(i); txy(i) qy(i) tyz(i); txz(i) tyz(i) qz(i)]; % matrix of all stresses
E = eye(3); % matrix of units
V, D, W] = eig(S); % D = matrix with principal stresses on the diagonals
HS = [D(3, 3) ; D(2, 2); D(1, 1)]; % vector with principal stresses
V1 = [W(1, 3); W(2, 3); W(3, 3)]; % vector of eigendirection 1
q = quiver3(msh.Nodes(1,i), msh.Nodes(2,i), msh.Nodes(3,i), W(1,1), W(2,1), W(3,1)); % msh.Nodes -> coordinates of the nodes
set(q,'AutoScale','on', 'AutoScaleFactor', 5,);
Now to my question: I would like to display these vectors as a kind of streamlines. Lines that follow the approximate directions of the vectors. Preferably in variable quantity.
First I tried this with the streamline command, but it doesn't work because I don't have a real vector field, but only plot single vectors at certain xyz coordinates.
I am very grateful for any help.
darova on 28 Jul 2021
% generata data
n = 20; % number of points
x = rand(20,1);
y = rand(20,1);
u = rand(20,1);
v = rand(20,1);
% calculate vector field (interpolate the data)
[x1,y1] = meshgrid(0:0.1:1);
u1 = griddata(x,y,u,x1,y1);
v1 = griddata(x,y,v,x1,y1);
x00 = 0:.1:1;