/[escript]/branches/refine/buckley/src/dump.cc
ViewVC logotype

Contents of /branches/refine/buckley/src/dump.cc

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3668 - (show annotations)
Wed Nov 16 01:49:46 2011 UTC (7 years, 10 months ago) by jfenwick
File size: 6984 byte(s)
Stage 1 rename
1 #include <iostream>
2 #include "OctTree.h"
3
4 using namespace std;
5 using namespace buckley;
6
7 void f(const OctCell& c, void* v)
8 {
9
10
11 static int i=1;
12 cout << i++ << " " << (c.centre[0]-c.sides[0]/2) << ' ' << (c.centre[1]-c.sides[1]/2) << ' ' << (c.centre[2]-c.sides[2]/2) << '\n';
13 cout << i++ << " " << (c.centre[0]+c.sides[0]/2) << ' ' << (c.centre[1]-c.sides[1]/2) << ' ' << (c.centre[2]-c.sides[2]/2) << '\n';
14 cout << i++ << " " << (c.centre[0]+c.sides[0]/2) << ' ' << (c.centre[1]+c.sides[1]/2) << ' ' << (c.centre[2]-c.sides[2]/2) << '\n';
15 cout << i++ << " " << (c.centre[0]-c.sides[0]/2) << ' ' << (c.centre[1]+c.sides[1]/2) << ' ' << (c.centre[2]-c.sides[2]/2) << '\n';
16
17 cout << i++ << " " << (c.centre[0]-c.sides[0]/2) << ' ' << (c.centre[1]-c.sides[1]/2) << ' ' << (c.centre[2]+c.sides[2]/2) << '\n';
18 cout << i++ << " " << (c.centre[0]+c.sides[0]/2) << ' ' << (c.centre[1]-c.sides[1]/2) << ' ' << (c.centre[2]+c.sides[2]/2) << '\n';
19 cout << i++ << " " << (c.centre[0]+c.sides[0]/2) << ' ' << (c.centre[1]+c.sides[1]/2) << ' ' << (c.centre[2]+c.sides[2]/2) << '\n';
20 cout << i++ << " " << (c.centre[0]-c.sides[0]/2) << ' ' << (c.centre[1]+c.sides[1]/2) << ' ' << (c.centre[2]+c.sides[2]/2) << '\n';
21 }
22
23 void cen_pts(const OctCell& c, void* v)
24 {
25 int& i=*(reinterpret_cast<int*>(v));
26 cout << i++ << ' ' << c.centre[0] << ' ' << c.centre[1] << ' ' << c.centre[2] << endl;
27 }
28
29 void cen_elts(const OctCell& c, void* v)
30 {
31 int& i=*(reinterpret_cast<int*>(v));
32 cout << i << " 15 3 0 "<<i <<" 0 " << i << endl; // 15 is the code for a single point element
33 i++;
34 }
35
36 void g(const OctCell& c, void* v)
37 {
38 static int e=1;
39 static int i=1;
40
41 cout << e << " 5 0 "; // element number
42 cout << i << ' ' << (i+1) << ' ' << (i+2) << ' ' << (i+3) << ' ';
43 cout << (i+4) << ' ' << (i+5) << ' ' << (i+6) << ' ' << (i+7) << endl;
44 e++;
45 i+=8;
46 }
47
48 void countleaves(const OctCell& c, void* v)
49 {
50 (*reinterpret_cast<int*>(v))++;
51 }
52
53 void dumpGrid(OctTree& ot)
54 {
55 cout << "$MeshFormat\n2.2 0 8\n$EndMeshFormat\n$Nodes\n";
56 int c=0;
57 ot.walkLeaves(countleaves, &c);
58 cout << c*8 << endl;
59 int i=1;
60 ot.walkLeaves(f,0);
61 cout << "$EndNodes\n";
62 cout << "$Elements\n";
63 cout << c << endl;
64 ot.walkLeaves(g, 0);
65 i=1;
66 cout << "$EndElements\n";
67
68 }
69
70 void dumpCen(OctTree& ot)
71 {
72 int c=0;
73 ot.walkLeaves(countleaves, &c);
74 cout << c << endl;
75 int i=1;
76 ot.walkLeaves(cen_pts, &i);
77 cout << "$EndNodes\n";
78 cout << "$Elements\n";
79 cout << c << endl;
80 i=1;
81 ot.walkLeaves(cen_elts, &i);
82 cout << "$EndElements\n";
83
84 }
85
86 const int HANG_NODE=1;
87
88
89 void countHang(const OctCell& c, void* v)
90 {
91 for (int i=0;i<8;++i)
92 {
93 if (c.leafinfo->pmap[i]==HANG_NODE)
94 {
95 (*reinterpret_cast<int*>(v))++;
96 }
97 }
98 }
99
100 void hangPoints(const OctCell& c, void* v)
101 {
102
103
104
105 int& k=*(reinterpret_cast<int*>(v));
106 double corners[8][3]; // perhaps not efficient but this is for debug only
107 corners[0][0]=(c.centre[0]-c.sides[0]/2); corners[0][1]=(c.centre[1]-c.sides[1]/2); corners[0][2]=(c.centre[2]-c.sides[2]/2);
108 corners[1][0]=(c.centre[0]+c.sides[0]/2); corners[1][1]=(c.centre[1]-c.sides[1]/2); corners[1][2]=(c.centre[2]-c.sides[2]/2);
109 corners[2][0]=(c.centre[0]+c.sides[0]/2); corners[2][1]=(c.centre[1]+c.sides[1]/2); corners[2][2]=(c.centre[2]-c.sides[2]/2);
110 corners[3][0]=(c.centre[0]-c.sides[0]/2); corners[3][1]=(c.centre[1]+c.sides[1]/2); corners[3][2]=(c.centre[2]-c.sides[2]/2);
111
112 corners[4][0]=(c.centre[0]-c.sides[0]/2); corners[4][1]=(c.centre[1]-c.sides[1]/2); corners[4][2]=(c.centre[2]+c.sides[2]/2);
113 corners[5][0]=(c.centre[0]+c.sides[0]/2); corners[5][1]=(c.centre[1]-c.sides[1]/2); corners[5][2]=(c.centre[2]+c.sides[2]/2);
114 corners[6][0]=(c.centre[0]+c.sides[0]/2); corners[6][1]=(c.centre[1]+c.sides[1]/2); corners[6][2]=(c.centre[2]+c.sides[2]/2);
115 corners[7][0]=(c.centre[0]-c.sides[0]/2); corners[7][1]=(c.centre[1]+c.sides[1]/2); corners[7][2]=(c.centre[2]+c.sides[2]/2);
116
117 for (int i=0;i<8;++i)
118 {
119 if (c.leafinfo->pmap[i]==HANG_NODE)
120 {
121 cout << k++ << ' ' << corners[i][0] << ' ' << corners[i][1] << ' ' << corners[i][2] << endl;
122 }
123 }
124 }
125
126 void dumpHang(OctTree& ot)
127 {
128 cout << "$MeshFormat\n2.2 0 8\n$EndMeshFormat\n$Nodes\n";
129 int c=0;
130 ot.walkLeaves(countHang, &c);
131 cout << c << endl;
132 int i=1;
133 ot.walkLeaves(hangPoints, &i);
134 cout << "$EndNodes\n";
135 cout << "$Elements\n";
136 cout << c << endl;
137 for (int k=1;k<=c;++k)
138 {
139 cout << k << " 15 3 0 "<<k <<" 0 " << k << endl; // 15 is the code for a single point element
140 }
141 cout << "$EndElements\n";
142 }
143
144
145 void printLInfo(const OctCell& c, void* v)
146 {
147 if (c.leaf)
148 {
149 cout << c.leafinfo << " belongs to " << &c << " " << c.centre[0] << "," << c.centre[1] << "," << c.centre[2] << "," << endl;
150 }
151 }
152
153 void neigh(const OctCell& c, void* v)
154 {
155 if (c.leaf)
156 {
157 cout << &c << endl;
158 for (int i=0;i<6;++i)
159 {
160 cout << " " << i << ": ";
161 if (c.leafinfo->next[i])
162 {
163 cout << c.leafinfo->next[i];
164 if (c.leafinfo->next[i]->owner)
165 {
166 cout << " at depth " << c.leafinfo->next[i]->owner->depth <<endl;
167 }
168 else
169 {
170 cout << " NULL owner\n";
171 }
172 }
173 else
174 {
175 cout << "null\n";
176 }
177 }
178 }
179 }
180
181 void maintest()
182 {
183 int c=0;
184 // cout << "$MeshFormat\n2.2 0 8\n$EndMeshFormat\n$Nodes\n";
185 OctTree ot(1,1,1);
186 ot.allSplit(6);
187
188
189 // ot.allCollapse(7);
190 // ot.collapsePoint(0.120, 0.120, 0.120,1);
191
192 // ot.collapsePoint(0.5, 0.5, 0.5, 1);
193
194 // ot.collapsePoint(0.51, 0.51, 0.51, 1);
195 // ot.collapsePoint(0.49, 0.49, 0.49, 2);
196
197 //return 0;
198 // ot.debug();
199
200
201
202 ot.collapsePoint(0.49, 0.49, 0.49, 2);
203
204
205 // dumpGrid(ot);
206 //return 0;
207
208 // ot.walkLeaves(printLInfo, 0);
209 // cout << "----------\n";
210 // ot.walkLeaves(neigh, 0);
211 // cout << "\n\n";
212
213
214
215
216 ot.collapsePoint(0.1, 0.49, 0.49, 3);
217
218
219 // ot.walkLeaves(printLInfo, 0);
220 /* cout << "----------\n";
221 ot.walkLeaves(neigh, 0);
222 cout << "\n\n"; */
223
224
225 ot.collapsePoint(0.3, 0.49, 0.49, 3);
226 ot.collapsePoint(0.7, 0.49, 0.49, 3);
227 ot.collapsePoint(0.99, 0.49, 0.49,3);
228
229 ot.collapsePoint(0.1, 0.99, 0.99,3);
230 ot.collapsePoint(0.49, 0.99, 0.99,3);
231 ot.collapsePoint(0.99, 0.99, 0.99,3);
232
233
234 // ot.splitPoint(0.512, 0.126,0.99, 10);
235
236 // ot.splitPoint(0.2, 0.8,0.01, 30);
237 cerr << "Assigning IDs\n";
238 ot.assignIDs();
239 cerr << "Done assigning IDs\n";
240 dumpGrid(ot);
241 //dumpCen(ot);
242
243
244 }
245
246
247
248
249 int main()
250 {
251 // maintest();
252 OctTree ot(1,1,1);
253 ot.allSplit(3);
254
255 ot.collapsePoint(0,0,0,1);
256 ot.collapsePoint(1,1,1,1);
257 ot.assignIDs();
258 dumpHang(ot);
259 //dumpGrid(ot);
260 //dumpCen(ot);
261
262 }

  ViewVC Help
Powered by ViewVC 1.1.26