FEM grid generation with inclined sides

Question

Left Terry on 6 Feb 2022

0
Link

Direct link to this question

https://in.mathworks.com/matlabcentral/answers/1643955-fem-grid-generation-with-inclined-sides

Edited: Left Terry on 6 Feb 2022

Hello. The following code works properly only for x_i = 0, y_i = 0, theta_x = 0 and theta_y = 0. For any other values of these four variables i do not get the desired grid. Any suggestion or any help would be really useful. Also any suggestion in order to warn the user in the case that enters a non numerical input would be also helpful. Thank you in advance. Please keep in mind that i am a beginner.

function [offsets,lengths,elements,angles] = grid_t(x_i,y_i,Lx,Ly,nex,ney,theta_x,theta_y)
% grid
% This programm creates a grid according to the users inputs which are the x- and y- axes offsets (x_i,y_i) (given in meters) and lengths (Lx,Ly) (also given in meters),
% the number of elements on each direction (nex,ney) and the angles between the grid and the axes (theta_x,theta_y). The angles must be given in degrees.
% Theta_x is measured from the horizontal to the nearest side of the grid. Theta_y is measured from the vertical axis to the nearest side of the grid.
% E.g. typing grid_t(1,3,5,2,6,4,20,10) in the command window the resulting grid will start from (x_i,y_i)=(1,3) with x-axis length of 5 units, y-axis length of 2 units,
% number of elements in the x-direction will be equal to 6 and in the y-direction will be equal to 4 and the angles form x-axis and y-axis will be 20 degrees and 10 degrees respectively.
if Lx > 0 && Ly > 0 && nex > 0 && ney > 0 && mod(nex,1) == 0 && mod(ney,1) == 0 % Lengths, Lx and Ly, must be positive numbers and the number of elements, nex and ney, must be positive and integers.
theta_x = pi*theta_x/180; % Converting degrees to radians
theta_y = pi*theta_y/180; % Converting degrees to radians
x0 = Ly*sin(theta_y); % Maximum displecement in the x-axis
y0 = Lx*sin(theta_x); % Maximum displecement in the y-axis
x_f = x_i+Lx*cos(theta_x)+x0; % Last point in x-direction
y_f = y_i+Ly*cos(theta_y)+y0; % Last point in y-direction
elements      = zeros(nex*ney,4);         % Creates a nex*ney by 4 matrix of all zeros. This will be the elements matrix. The total number of elements is nex*ney.
nodes_coords  = zeros((nex+1)*(ney+1),2); % Creates a (nex+1)*(ney+1) by 2 matrix. This will be the x- and y- coordinates matrix for the nodes. The total number of  nodes is (nex+1)*(ney+1).
x_coords      = zeros(nex+1,1);           % Creates a nex+1 by 1 matrix (column vector) of all zeros. This will be the x-coordinates matrix for the nodes. The total number of  nodes in the x-direction is nex+1
y_coords      = zeros(ney+1,1);           % Creates a ney+1 by 1 matrix (column vector) of all zeros. This will be the y-coordinates matrix for the nodes. The total number of  nodes in the y-direction is ney+1
 
x_coords(:,1) = linspace(x_i,x_f,nex+1) % Generates nex+1 points between x_i (initial point in x-direction) and x_f (last point in x-direction). The points are displayd as a column vector.
y_coords(:,1) = linspace(y_i,y_f,ney+1) % Generates ney+1 points between y_i (initial point in x-direction) and y_f (last point in x-direction). The points are displayd as a column vector.
axis([x_i-1 x_f+1 y_i-1 y_f+1])
% This section creates the nodes' x- and y- coordinates matrix. It is a (nex+1)*(ney+1) by 2 matrix.
for j = 1:ney+1
for i = 1:nex+1
    
    nodes_coords((nex+1)*(j-1)+i,1) = x_coords(i); % This line creates the y-coordinates of the nodes and places them in the 1st column of the above-mentioned matrix.
    nodes_coords((nex+1)*(j-1)+i,2) = y_coords(j); % This line creates the x-coordinates of the nodes and places them in the 2nd column of the above-mentioned matrix.
       
end
end
for j = 1:ney
for i = 1:nex
    
    elements(i+nex*(j-1),:) = [i+(j-1)*(nex+1) i+1+(j-1)*(nex+1) i+2+nex+(j-1)*(nex+1) i+1+nex+(j-1)*(nex+1)]; % Creates the elements matrix
    
    % Graphing the horizontal lines
    
    line([nodes_coords(elements(i,1),1) nodes_coords(elements(i,2),1)], [nodes_coords(elements(i,1),2)+(j-1)*Ly/ney nodes_coords(elements(i,2),2)+(j-1)*Ly/ney],'Marker','.');
    
    if j == ney
        
       line([nodes_coords(elements(i,1),1) nodes_coords(elements(i,2),1)], [nodes_coords(elements(i,1),2)+j*Ly/ney nodes_coords(elements(i,2),2)+j*Ly/ney],'Marker','.');
       
    end
    
    % Graphing the vertical lines
    
    line([nodes_coords(elements(i,1),2)+(i-1)*Lx/nex nodes_coords(elements(i,2),2)+(i-1)*Lx/nex], [y_coords(j) y_coords(j+1)],'Marker','.');
    
    if i == nex
      
       line([nodes_coords(elements(i,1),2)+(i)*Lx/nex nodes_coords(elements(i,2),2)+(i)*Lx/nex], [y_coords(j) y_coords(j+1)],'Marker','.');
       
    end
    end
end
else
    disp("Give appropriate inputs")
    
end