Drexel dragonThe Math ForumDonate to the Math Forum



Search All of the Math Forum:

Views expressed in these public forums are not endorsed by Drexel University or The Math Forum.


Math Forum » Discussions » Software » comp.soft-sys.matlab

Topic: refining a Quadrilateral mesh.
Replies: 3   Last Post: Jan 31, 2011 11:00 AM

Advanced Search

Back to Topic List Back to Topic List Jump to Tree View Jump to Tree View   Messages: [ Previous | Next ]
Kwesi Moore

Posts: 19
Registered: 11/27/09
refining a Quadrilateral mesh.
Posted: Jan 7, 2011 1:16 PM
  Click to see the message monospaced in plain text Plain Text   Click to reply to this topic Reply

I have a quadrilateral mesh.i want to create refinements within the mesh.
This is what i have done so far....

element=[1 2 4 5;2 3 5 6;4 5 7 8;5 6 8 9];
coordinates=[-1 0 1 -1 0 1 -1 0 1;-1 -1 -1 0 0 0 1 1 1]';
dirichlet=[]1 2 1 4;2 3 4 7]';
neumann=[3 6 7 8;6 9 8 9]';
nodes2edge=[1 2 1 4 3 2 5 4 5 7 6 8;4 5 2 5 6 3 6 7 8 8 9 9]';

function [coordinates,NewElement,dirichlet,neumann,ParentElement]=...
refinement(coordinates, elements, dirichlet,neumann,nodes2edge)

% Compute number of edges, nodes and elements
no_of_edges = full(max(max(nodes2edge)));
no_of_nodes = size(coordinates,1);
no_of_elements = size(elements,1);

% Create New Arrays
Node2Edge = zeros(no_of_edges,2);
% Get the indices of the nodes
[a,b,c] = find(nodes2edge);
Node2Edge(c,:) = [a,b];

% Coordinates
coordinates(no_of_nodes+1:no_of_nodes+no_of_edges,:)=...
(coordinates(Node2Edge(:,1),:)+coordinates(Node2Edge(:,2),:))/2;

% Elements
NewElement = zeros(4*size(elements,1),3);
ParentElement = zeros(size(NewElement,1),1);

% Element Refinements
ell = 0;

for j = 1:no_of_elements
% describes the order of the faces >>left (face 1) right south north.
I = [nodes2edge(elements(j,1),elements(j,3));
nodes2edge(elements(j,2),elements(j,4));
nodes2edge(elements(j,1),elements(j,2));
nodes2edge(elements(j,3),elements(j,4))];
NewNode=(I+no_of_nodes)';
NewElement(ell+1:ell+5,:)=[NewNode([1 2 3 4]);
elements(j,1), NewNode(1), NewNode(2), NewNode(3);
NewNode(2), elements(j,2), NewNode(3), NewNode(4);
NewNode(3), NewNode(4), elements(j,3),NewNode(1);
NewNode(4), NewNode(1), NewNode(2), elements(j,4)];

ParentElement(ell+1:ell+5)=j;
ell=ell+5;
end

% Dirichlet and Neumann
tmp=0;
dirichletnr=diag(nodes2edge(dirichlet(:,1),dirichlet(:,2)));
for k=1:size(dirichlet,1)
dirichlet=[dirichlet(1:k-1+tmp,:);
dirichlet(k+tmp,1),dirichletnr(k)+no_of_nodes;
dirichletnr(k)+no_of_nodes,dirichlet(k+tmp,2);
dirichlet(k+1+tmp:size(dirichlet,1),:)];
tmp=tmp+1;
end

tmp=0;
if ~isempty(neumann)
neumannnr=diag(nodes2edge(neumann(:,1),neumann(:,2)));
end
for k=1:size(neumann,1)
neumann = [neumann(1:k-1+tmp,:);
neumann(k+tmp,1),neumannnr(k)+no_of_nodes;
neumannnr(k)+no_of_nodes,neumann(k+tmp,2);
neumann(k+1+tmp:size(neumann,1),:)];
tmp=tmp+1;
end



Point your RSS reader here for a feed of the latest messages in this topic.

[Privacy Policy] [Terms of Use]

© Drexel University 1994-2014. All Rights Reserved.
The Math Forum is a research and educational enterprise of the Drexel University School of Education.