1 function gridp = gridpart(grid,eind)
2 %
function gridp = gridpart(grid,eind)
3 %
function extracting a part of a .triagrid or .rectgrid defined by the given
4 % element indices in the vector
'eind'.
6 % @note The neighbour information of the
new resulting boundaries is set to
10 % distance-information in the
new boundary elements are simply copied. I.e.
11 % these fields
do not completely meet the definition in the constructor. They
12 % might be chosen slightly different, such that the
'gridp' would be really
13 % identical to the result generated by the constructor on the subset of points.
16 % eind: vector of cell indices which shall be extracted from the grid.
19 % gridp: the partial grid of type .gridbase with extracted cells
'eind'.
21 % Bernard Haasdonk 15.5.2007
24 gridp.nelements = length(eind);
26 % generate elementid translation map: T \mapsto T_{local}
27 new_el_id = zeros(1,grid.nelements);
28 new_el_id(eind) = 1:length(eind);
30 % generate vertex translation map: V \mapsto V_{local}
31 new_vertex_id = zeros(1,grid.nvertices);
32 mask = zeros(1,grid.nvertices);
33 mask(grid.VI(eind,:))= 1;
35 new_vertex_id(vind) = 1:length(vind);
37 % set number of vertices, area vector, and inverse area vector
38 gridp.nvertices = max(new_vertex_id);
39 gridp.A = gridp.A(eind);
40 gridp.Ainv = gridp.Ainv(eind);
42 % set vertex coordinates
43 gridp.X = grid.X(vind);
44 gridp.Y = grid.Y(vind);
47 gridp.VI = grid.VI(eind,:); %: VI(i,j) is the global vertex index of j-th
48 gridp.VI = new_vertex_id(gridp.VI);
50 % set coordinates of elements' barycenters
51 gridp.CX = gridp.CX(eind);
52 gridp.CY = gridp.CY(eind);
54 % ATTENTION! Set neighbour indices
55 gridp.NBI = grid.NBI(eind,:);
56 i = find(gridp.NBI>0);
57 gridp.NBI(i) = new_el_id(gridp.NBI(i));
58 i = find(gridp.NBI == 0);
63 gridp.INB = grid.INB(eind,:);
65 gridp.EL = grid.EL(eind,:);
66 gridp.DC = grid.DC(eind,:);
67 gridp.NX = grid.NX(eind,:);
68 gridp.NY = grid.NY(eind,:);
69 gridp.ECX = grid.ECX(eind,:);
70 gridp.ECY = grid.ECY(eind,:);
71 gridp.SX = grid.SX(eind,:);
72 gridp.SY = grid.SY(eind,:);
73 gridp.ESX = grid.ESX(eind,:);
74 gridp.ESY = grid.ESY(eind,:);
75 gridp.DS = grid.DS(eind,:);
alpha
geometry bound (simultaneously satisfying and
Base class for all grid classes.
hmin
minimal element-diameter