# How to make every element in the form of isoparametric, Lagrangian quadratic element?

2 views (last 30 days)
Lokesh Mangal on 17 Jan 2020
Answered: KSSV on 17 Jan 2020
%% Code
clc ; clear;
% Dimensions of the plate
L = 1 ; % Length of the Plate along y-axes
B = 1 ; % Breadth of the Plate along x-axes
% Number of Elements required
Nx = 4 ; % Number of Elements along X-axes
Ny = 4 ; % Number of Elements along Y-axes
%----------------------------------------
nel = Nx*Ny ; % Total Number of Elements in the Mesh
nnel = 9 ; % Number of nodes per Element
% Number of points on the Length and Breadth
npx = 2*Nx+1 ;
npy = 2*Ny+1 ;
nnode = npx*npy ;
nx = linspace(0,L,npx) ;
ny = linspace(0,B,npy) ;
[xx yy] = meshgrid(nx,ny) ;
% To get the Nodal Connectivity Matrix
coordinates = [xx(:) yy(:)] ;
NodeNo = 1:nnode ;
nodes = zeros(nel,nnel) ;
NodeNo = reshape(NodeNo,npx,npy);
%%
nodes(:,1) = reshape(NodeNo(1:2:npx-2,1:2:npy-2),nel,1);
nodes(:,2) = reshape(NodeNo(3:2:npx,1:2:npy-2),nel,1);
nodes(:,3) = reshape(NodeNo(3:2:npx,3:2:npy),nel,1);
nodes(:,4) = reshape(NodeNo(1:2:npx-2,3:2:npy),nel,1);
nodes(:,5) = reshape(NodeNo(2:2:npx-1,1:2:npy-2),nel,1);
nodes(:,6) = reshape(NodeNo(3:2:npx,2:2:npy-1),nel,1);
nodes(:,7) = reshape(NodeNo(2:2:npx-1,3:2:npy),nel,1);
nodes(:,8) = reshape(NodeNo(1:2:npx-2,2:2:npy-1),nel,1);
nodes(:,9) = reshape(NodeNo(2:2:npx-1,2:2:npy-1),nel,1);
%%
X = zeros(nnel,nel) ;
Y = zeros(nnel,nel) ;
%%
for iel = 1:nel
Y(:,iel) = coordinates(nodes(iel,:),1) ;
X(:,iel) = coordinates(nodes(iel,:),2) ;
end
%%
x1 = X(1,:);
x2 = X(2,:);
x3 = X(3,:);
x4 = X(4,:);
x5 = X(5,:);
x6 = X(6,:);
x7 = X(7,:);
x8 = X(8,:);
x9 = X(9,:);
y1 = Y(1,:);
y2 = Y(2,:);
y3 = Y(3,:);
y4 = Y(4,:);
y5 = Y(5,:);
y6 = Y(6,:);
y7 = Y(7,:);
y8 = Y(8,:);
y9 = Y(9,:);
%%
figure
plot(x1,y1,'b.','MarkerSize',10)
xlim([-(B/10) B+(B/10)])
ylim([-(L/10) L+(L/10)])
hold on
plot(x2,y2,'b.','MarkerSize',10)
plot(x3,y3,'b.','MarkerSize',10)
plot(x4,y4,'b.','MarkerSize',10)
plot(x5,y5,'b.','MarkerSize',10)
plot(x6,y6,'b.','MarkerSize',10)
plot(x7,y7,'b.','MarkerSize',10)
plot(x8,y8,'b.','MarkerSize',10)
plot(x9,y9,'b.','MarkerSize',10)
hold off
%%
for i = 0 : (B/Nx) : B
line([i i],[0 L],'LineStyle','-','Color','black') % To draw vertical lines
for j = 0 : (L/Ny) : L
line([0,B],[j j],'LineStyle','-','Color','black') % To draw horizontal lines
end
end

KSSV on 17 Jan 2020
nodes = [nodes(:,1) nodes(:,5) nodes(:,2) nodes(:,6) nodes(:,3) nodes(:,7) nodes(:,4) nodes(:,8) nodes(:,9)] ;