← Back to branch summary

~ubuntu-branches/ubuntu/hardy/suitesparse/hardy

« back to all changes in this revision

Viewing changes to MESHND/meshnd.m

  • Committer: Bazaar Package Importer
  • Author(s): Rafael Laboissiere, Rafael Laboissiere, Ondrej Certik
  • Date: 2008-02-21 14:46:50 UTC
  • mfrom: (1.1.2 upstream) (5.1.1 hardy)
  • Revision ID: james.westby@ubuntu.com-20080221144650-tgeppgj0t7s759i8
Tags: 3.1.0-3
[ Rafael Laboissiere ]
* Upload to unstable

[ Ondrej Certik ]
* XS-DM-Upload-Allowed: yes field added
* Ondrej Certik added to uploaders

Show diffs side-by-side

added added

removed removed

Lines of Context:
MESHND/meshnd.m
 
1
function [G, p, pinv, Gnew] = meshnd (arg1,n,k)
 
2
%MESHND creation and nested dissection of a regular 2D or 3D mesh.
 
3
% [p G pinv Gnew] = meshnd (m,n) constructs a m-by-n 2D mesh G, and then finds
 
4
% a permuted mesh Gnew where Gnew = pinv(G) and G = p(Gnew).  meshnd(m,n,k)
 
5
% creates an m-by-n-by-k 3D mesh.
 
6
%
 
7
% [p G pinv Gnew] = meshnd (G) does not construct G, but uses the mesh G as
 
8
% given on input instead.
 
9
%
 
10
% Example:
 
11
% [G p pinv Gnew] = meshnd (4,5) ;
 
12
%
 
13
% returns
 
14
%    Gnew =
 
15
%        1     2    17     9    10
 
16
%        7     8    18    15    16
 
17
%        3     5    19    11    13
 
18
%        4     6    20    12    14
 
19
%    G =
 
20
%        1     2     3     4     5
 
21
%        6     7     8     9    10
 
22
%       11    12    13    14    15
 
23
%       16    17    18    19    20
 
24
%
 
25
% With no inputs, a few example meshes are generated and plotted.
 
26
%
 
27
% See also nested, numgrid.
 
28
 
 
29
% Copyright 2007, Timothy A. Davis, Univ. of Florida
 
30
 
 
31
% get the inputs and create the mesh if not provided on input
 
32
if (nargin == 0)
 
33
 
 
34
    % run a simple example
 
35
    meshnd_example ;
 
36
 
 
37
elseif (nargin == 1)
 
38
 
 
39
    % the mesh is provided on input
 
40
    G = arg1 ;
 
41
    [m n k] = size (G) ;
 
42
 
 
43
elseif (nargin == 2)
 
44
 
 
45
    % create the m-by-n-by-k mesh in "natural" (row-major) order.  This is how
 
46
    % a typical 2D mesh is ordered.  A column-major order would be better, since
 
47
    % in that case G(:) would equal 1:(m*n) ... but let's stick with tradition.
 
48
    m = arg1 ;
 
49
    k = 1 ;
 
50
    G = reshape (1:(m*n*k), n, m, k)' ;
 
51
 
 
52
elseif (nargin == 3)
 
53
 
 
54
    % create the m-by-n-by-k mesh in column-major order.  The first m-by-n-by-1
 
55
    % slice is in column-major order, followed by all the other slices 2 to k.
 
56
    m = arg1 ;
 
57
    G = reshape (1:(m*n*k), m, n, k) ;
 
58
 
 
59
else
 
60
 
 
61
    error ('Usage: [G p pinv Gnew] = meshnd(G), meshnd(m,n) or meshnd(m,n,k)') ;
 
62
 
 
63
end
 
64
 
 
65
if (nargout > 1)
 
66
    p = nd2 (G)' ;          % order the mesh
 
67
end
 
68
 
 
69
if (nargout > 2)
 
70
    pinv (p) = 1:(m*n*k) ;  % find the inverse permutation
 
71
end
 
72
 
 
73
if (nargout > 3)
 
74
    Gnew = pinv (G) ;       % find the permuted mesh
 
75
end
 
76
 
 
77
%-------------------------------------------------------------------------------
 
78
 
 
79
function p = nd2 (G)
 
80
%ND2 p = nd2 (G) permutes a 2D or 3D mesh G.
 
81
% Compare with nestdiss which uses p as a scalar offset and returns a modified
 
82
% mesh G that corresponds to Gnew in meshnd.  Here, the scalar offset p in
 
83
% nestdiss is not needed.  Instead, p is a permutation, and the modified mesh
 
84
% Gnew is not returned.
 
85
 
 
86
[m n k] = size (G) ;
 
87
 
 
88
if (max ([m n k]) <= 2)
 
89
 
 
90
    % G is small; do not cut it
 
91
    p = G (:) ;
 
92
 
 
93
elseif k >= max (m,n)
 
94
 
 
95
    % cut G along the middle slice, cutting k in half
 
96
    s = ceil (k/2) ;
 
97
    middle = G (:,:,s) ;
 
98
    p = [(nd2 (G (:,:,1:s-1))) ; (nd2 (G (:,:,s+1:k))) ; middle(:)] ;
 
99
 
 
100
elseif n >= max (m,k)
 
101
 
 
102
    % cut G along the middle column, cutting n in half
 
103
    s = ceil (n/2) ;
 
104
    middle = G (:,s,:) ;
 
105
    p = [(nd2 (G (:,1:s-1,:))) ; (nd2 (G (:,s+1:n,:))) ; middle(:)] ;
 
106
 
 
107
else   
 
108
 
 
109
    % cut G along the middle row, cutting m in half
 
110
    s = ceil (m/2) ;
 
111
    middle = G (s,:,:) ;
 
112
    p = [(nd2 (G (1:s-1,:,:))) ; (nd2 (G (s+1:m,:,:))) ; middle(:)] ;
 
113
 
 
114
end