DGtal  0.9.2
testEmbedding.cpp
1 #include "DGtal/math/linalg/EigenSupport.h"
2 #include "DGtal/dec/DiscreteExteriorCalculus.h"
3 #include "DGtal/dec/DiscreteExteriorCalculusSolver.h"
4 #include "DGtal/dec/DiscreteExteriorCalculusFactory.h"
5 
6 #define NOVIEWER // comment to enable 3d viewers
7 
8 #if !defined(NOVIEWER)
9 #include "DGtal/io/viewers/Viewer3D.h"
10 #endif
11 #include "DGtal/io/boards/Board2D.h"
12 #include "DGtal/base/Common.h"
13 #include "DGtal/helpers/StdDefs.h"
14 
15 using namespace DGtal;
16 using namespace std;
17 using namespace Eigen;
18 
19 template <typename OperatorAA, typename OperatorBB>
20 bool
21 equal(const OperatorAA& aa, const OperatorBB& bb)
22 {
23  return MatrixXd(aa.myContainer) == MatrixXd(bb.myContainer);
24 }
25 
26 int main(int , char** )
27 {
28 #if !defined(NOVIEWER)
30 #endif
34 
35 #if !defined(NOVIEWER)
36  QApplication app(argc, argv);
37  Z3i::KSpace kspace_3d;
38 #endif
39  Z2i::KSpace kspace_2d;
40 
41 #if !defined(NOVIEWER)
42  Viewer viewer1(kspace_3d);
43  viewer1.show();
44  viewer1.setWindowTitle("embedding_1d_calculus_3d");
45  viewer1.camera()->setPosition( Vec(2,2,2) );
46  viewer1.camera()->setUpVector( Vec(0,0,1), false );
47  viewer1.camera()->lookAt( Vec(0,0,0) );
48 
49  Viewer viewer2(kspace_3d);
50  viewer2.show();
51  viewer2.setWindowTitle("embedding_2d_calculus_3d");
52  viewer2.camera()->setPosition( Vec(2,2,2) );
53  viewer2.camera()->setUpVector( Vec(0,0,2), false );
54  viewer2.camera()->lookAt( Vec(0,0,0) );
55 
56  Viewer viewer3(kspace_3d);
57  viewer3.show();
58  viewer3.setWindowTitle("embedding_2d_calculus_3d_no_border");
59  viewer3.camera()->setPosition( Vec(2,2,2) );
60  viewer3.camera()->setUpVector( Vec(0,0,2), false );
61  viewer3.camera()->lookAt( Vec(0,0,0) );
62 #endif
63 
64  {
65  trace.beginBlock("1d manifold embedding");
66 
72 
74  typedef std::set<Calculus1D::SCell> SCells1D;
75  SCells1D ncells_1d_factory;
77  SCells1D cells_1d_manual;
78  Calculus1D calculus_1d_manual;
79  for (int kk=0; kk<31; kk++)
80  {
81  Calculus1D::KSpace::Point point;
82  point[0] = kk;
83  const Calculus1D::SCell cell = calculus_1d_manual.myKSpace.sCell(point);
84  calculus_1d_manual.insertSCell(cell, 1, kk == 0 || kk == 30 ? 1/2. : 1);
85  cells_1d_manual.insert(cell);
86  if (kk%2 != 0) ncells_1d_factory.insert(cell);
87  }
88  calculus_1d_manual.updateIndexes();
89  trace.info() << "calculus_1d_manual=" << calculus_1d_manual << endl;
90 
91  {
92  Board2D board;
93  board << Z2i::Domain(Z2i::Point(-1,-1), Z2i::Point(15,0));
94  board.setFillColor( DGtal::Color(128,128,128) );
95  for (Calculus1D::ConstIterator ii=calculus_1d_manual.begin(), ie=calculus_1d_manual.end(); ii!=ie; ii++)
96  {
97  const Z2i::Point point(calculus_1d_manual.myKSpace.uKCoord(ii->first, 0), 0);
98  const Z2i::KSpace::SCell cell = kspace_2d.sCell(point);
100  }
101  board.saveSVG("embedding_1d_calculus_1d.svg");
102  }
103 
105  const Calculus1D calculus_1d_factory = CalculusFactory::createFromNSCells<1>(ncells_1d_factory.begin(), ncells_1d_factory.end(), true);
107  trace.info() << "calculus_1d_factory=" << calculus_1d_factory << endl;
108 
109  Calculus2D calculus_2d_manual;
110  {
111  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(6,0)), 1, 1/2. );
112  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(6,1), Calculus2D::KSpace::POS) );
113  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(6,2)) );
114  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(7,2), Calculus2D::KSpace::POS) );
115  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(8,2)) );
116  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(8,1), Calculus2D::KSpace::NEG) );
117  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(8,0)) );
118  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(8,-1), Calculus2D::KSpace::NEG) );
119  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(8,-2)) );
120  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(7,-2), Calculus2D::KSpace::NEG) );
121  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(6,-2)) );
122  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(5,-2), Calculus2D::KSpace::NEG) );
123  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(4,-2)) );
124  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(3,-2), Calculus2D::KSpace::NEG) );
125  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(2,-2)) );
126  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(1,-2), Calculus2D::KSpace::NEG) );
127  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(0,-2)) );
128  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(-1,-2), Calculus2D::KSpace::NEG) );
129  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(-2,-2)) );
130  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(-2,-1), Calculus2D::KSpace::POS) );
131  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(-2,0)) );
132  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(-2,1), Calculus2D::KSpace::POS) );
133  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(-2,2)) );
134  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(-1,2), Calculus2D::KSpace::POS) );
135  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(0,2)) );
136  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(1,2), Calculus2D::KSpace::POS) );
137  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(2,2)) );
138  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(2,1), Calculus2D::KSpace::NEG) );
139  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(2,0)) );
140  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(1,0), Calculus2D::KSpace::NEG) );
141  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(0,0)), 1, 1/2.);
142  calculus_2d_manual.updateIndexes();
143  }
144  trace.info() << "calculus_2d_manual=" << calculus_2d_manual << endl;
145 
146  {
147  Board2D board;
148  board << Z2i::Domain(Z2i::Point(-2,-2), Z2i::Point(4,1));
149  board << calculus_2d_manual;
150  board.saveSVG("embedding_1d_calculus_2d.svg");
151  }
152 
154  typedef std::list<Calculus2D::SCell> SCells2D;
155  SCells2D ncells_2d_factory;
157 
158  {
159  ncells_2d_factory.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(6,1), Calculus2D::KSpace::POS) );
160  ncells_2d_factory.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(7,2), Calculus2D::KSpace::POS) );
161  ncells_2d_factory.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(8,1), Calculus2D::KSpace::NEG) );
162  ncells_2d_factory.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(8,-1), Calculus2D::KSpace::NEG) );
163  ncells_2d_factory.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(7,-2), Calculus2D::KSpace::NEG) );
164  ncells_2d_factory.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(5,-2), Calculus2D::KSpace::NEG) );
165  ncells_2d_factory.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(3,-2), Calculus2D::KSpace::NEG) );
166  ncells_2d_factory.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(1,-2), Calculus2D::KSpace::NEG) );
167  ncells_2d_factory.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(-1,-2), Calculus2D::KSpace::NEG) );
168  ncells_2d_factory.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(-2,-1), Calculus2D::KSpace::POS) );
169  ncells_2d_factory.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(-2,1), Calculus2D::KSpace::POS) );
170  ncells_2d_factory.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(-1,2), Calculus2D::KSpace::POS) );
171  ncells_2d_factory.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(1,2), Calculus2D::KSpace::POS) );
172  ncells_2d_factory.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(2,1), Calculus2D::KSpace::NEG) );
173  ncells_2d_factory.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(1,0), Calculus2D::KSpace::NEG) );
174  }
176  const Calculus2D calculus_2d_factory = CalculusFactory::createFromNSCells<1>(ncells_2d_factory.begin(), ncells_2d_factory.end(), true);
178  trace.info() << "calculus_2d_factory=" << calculus_2d_factory << endl;
179 
180  SCells2D cells_2d_manual;
181  {
182  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(6,0)) );
183  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(6,1), Calculus2D::KSpace::POS) );
184  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(6,2)) );
185  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(7,2), Calculus2D::KSpace::POS) );
186  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(8,2)) );
187  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(8,1), Calculus2D::KSpace::NEG) );
188  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(8,0)) );
189  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(8,-1), Calculus2D::KSpace::NEG) );
190  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(8,-2)) );
191  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(7,-2), Calculus2D::KSpace::NEG) );
192  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(6,-2)) );
193  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(5,-2), Calculus2D::KSpace::NEG) );
194  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(4,-2)) );
195  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(3,-2), Calculus2D::KSpace::NEG) );
196  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(2,-2)) );
197  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(1,-2), Calculus2D::KSpace::NEG) );
198  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(0,-2)) );
199  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(-1,-2), Calculus2D::KSpace::NEG) );
200  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(-2,-2)) );
201  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(-2,-1), Calculus2D::KSpace::POS) );
202  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(-2,0)) );
203  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(-2,1), Calculus2D::KSpace::POS) );
204  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(-2,2)) );
205  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(-1,2), Calculus2D::KSpace::POS) );
206  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(0,2)) );
207  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(1,2), Calculus2D::KSpace::POS) );
208  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(2,2)) );
209  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(2,1), Calculus2D::KSpace::NEG) );
210  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(2,0)) );
211  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(1,0), Calculus2D::KSpace::NEG) );
212  cells_2d_manual.push_back( calculus_2d_manual.myKSpace.sCell(Z2i::Point(0,0)) );
213  }
214 
215  Calculus3D calculus_3d_manual;
216  {
217  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(0,0,0)), 1, 1/2. );
218  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(1,0,0), Calculus3D::KSpace::POS) );
219  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,0,0)) );
220  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(3,0,0), Calculus3D::KSpace::POS) );
221  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(4,0,0)) );
222  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(4,1,0), Calculus3D::KSpace::POS) );
223  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(4,2,0)) );
224  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(4,3,0), Calculus3D::KSpace::POS) );
225  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(4,4,0)) );
226  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(3,4,0), Calculus3D::KSpace::NEG) );
227  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,4,0)) );
228  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(1,4,0), Calculus3D::KSpace::NEG) );
229  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(0,4,0)) );
230  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(0,3,0), Calculus3D::KSpace::NEG) );
231  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(0,2,0)) );
232  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(1,2,0), Calculus3D::KSpace::POS) );
233  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,2,0)) );
234  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,2,1), Calculus3D::KSpace::POS) );
235  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,2,2)) );
236  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,3,2), Calculus3D::KSpace::POS) );
237  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,4,2)) );
238  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,5,2), Calculus3D::KSpace::POS) );
239  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,6,2)) );
240  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,6,1), Calculus3D::KSpace::NEG) );
241  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,6,0)) );
242  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,6,-1), Calculus3D::KSpace::NEG) );
243  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,6,-2)) );
244  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,5,-2), Calculus3D::KSpace::NEG) );
245  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,4,-2)) );
246  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,3,-2), Calculus3D::KSpace::NEG) );
247  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,2,-2)), 1, 1/2. );
248  calculus_3d_manual.updateIndexes();
249  }
250  trace.info() << "calculus_3d_manual=" << calculus_3d_manual << endl;
251 
252 #if !defined(NOVIEWER)
254  viewer1 << Viewer::updateDisplay;
255 #endif
256 
258  typedef std::vector<Calculus3D::SCell> SCells3D;
259  SCells3D ncells_3d_factory;
261  {
262  ncells_3d_factory.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(1,0,0), Calculus3D::KSpace::POS) );
263  ncells_3d_factory.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(3,0,0), Calculus3D::KSpace::POS) );
264  ncells_3d_factory.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(4,1,0), Calculus3D::KSpace::POS) );
265  ncells_3d_factory.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(4,3,0), Calculus3D::KSpace::POS) );
266  ncells_3d_factory.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(3,4,0), Calculus3D::KSpace::NEG) );
267  ncells_3d_factory.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(1,4,0), Calculus3D::KSpace::NEG) );
268  ncells_3d_factory.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(0,3,0), Calculus3D::KSpace::NEG) );
269  ncells_3d_factory.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(1,2,0), Calculus3D::KSpace::POS) );
270  ncells_3d_factory.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,2,1), Calculus3D::KSpace::POS) );
271  ncells_3d_factory.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,3,2), Calculus3D::KSpace::POS) );
272  ncells_3d_factory.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,5,2), Calculus3D::KSpace::POS) );
273  ncells_3d_factory.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,6,1), Calculus3D::KSpace::NEG) );
274  ncells_3d_factory.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,6,-1), Calculus3D::KSpace::NEG) );
275  ncells_3d_factory.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,5,-2), Calculus3D::KSpace::NEG) );
276  ncells_3d_factory.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,3,-2), Calculus3D::KSpace::NEG) );
277  }
279  const Calculus3D calculus_3d_factory = CalculusFactory::createFromNSCells<1>(ncells_3d_factory.begin(), ncells_3d_factory.end(), true);
281  trace.info() << "calculus_3d_factory=" << calculus_3d_factory << endl;
282 
283  SCells3D cells_3d_manual;
284  {
285  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(0,0,0)) );
286  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(1,0,0), Calculus3D::KSpace::POS) );
287  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,0,0)) );
288  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(3,0,0), Calculus3D::KSpace::POS) );
289  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(4,0,0)) );
290  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(4,1,0), Calculus3D::KSpace::POS) );
291  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(4,2,0)) );
292  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(4,3,0), Calculus3D::KSpace::POS) );
293  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(4,4,0)) );
294  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(3,4,0), Calculus3D::KSpace::NEG) );
295  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,4,0)) );
296  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(1,4,0), Calculus3D::KSpace::NEG) );
297  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(0,4,0)) );
298  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(0,3,0), Calculus3D::KSpace::NEG) );
299  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(0,2,0)) );
300  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(1,2,0), Calculus3D::KSpace::POS) );
301  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,2,0)) );
302  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,2,1), Calculus3D::KSpace::POS) );
303  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,2,2)) );
304  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,3,2), Calculus3D::KSpace::POS) );
305  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,4,2)) );
306  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,5,2), Calculus3D::KSpace::POS) );
307  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,6,2)) );
308  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,6,1), Calculus3D::KSpace::NEG) );
309  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,6,0)) );
310  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,6,-1), Calculus3D::KSpace::NEG) );
311  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,6,-2)) );
312  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,5,-2), Calculus3D::KSpace::NEG) );
313  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,4,-2)) );
314  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,3,-2), Calculus3D::KSpace::NEG) );
315  cells_3d_manual.push_back( calculus_3d_manual.myKSpace.sCell(Z3i::Point(2,2,-2)) );
316  }
317  trace.beginBlock("checking operators");
318 
319  const Calculus1D::PrimalIdentity0 reorder_0cell_1d = calculus_1d_manual.reorder<0, PRIMAL>(cells_1d_manual.begin(), cells_1d_manual.end());
320  const Calculus2D::PrimalIdentity0 reorder_0cell_2d = calculus_2d_manual.reorder<0, PRIMAL>(cells_2d_manual.begin(), cells_2d_manual.end());
321  const Calculus3D::PrimalIdentity0 reorder_0cell_3d = calculus_3d_manual.reorder<0, PRIMAL>(cells_3d_manual.begin(), cells_3d_manual.end());
322  const Calculus1D::PrimalIdentity1 reorder_1cell_1d = calculus_1d_manual.reorder<1, PRIMAL>(cells_1d_manual.begin(), cells_1d_manual.end());
323  const Calculus2D::PrimalIdentity1 reorder_1cell_2d = calculus_2d_manual.reorder<1, PRIMAL>(cells_2d_manual.begin(), cells_2d_manual.end());
324  const Calculus3D::PrimalIdentity1 reorder_1cell_3d = calculus_3d_manual.reorder<1, PRIMAL>(cells_3d_manual.begin(), cells_3d_manual.end());
325  const Calculus1D::DualIdentity0 reorder_0cellp_1d = calculus_1d_manual.reorder<0, DUAL>(cells_1d_manual.begin(), cells_1d_manual.end());
326  const Calculus2D::DualIdentity0 reorder_0cellp_2d = calculus_2d_manual.reorder<0, DUAL>(cells_2d_manual.begin(), cells_2d_manual.end());
327  const Calculus3D::DualIdentity0 reorder_0cellp_3d = calculus_3d_manual.reorder<0, DUAL>(cells_3d_manual.begin(), cells_3d_manual.end());
328  const Calculus1D::DualIdentity1 reorder_1cellp_1d = calculus_1d_manual.reorder<1, DUAL>(cells_1d_manual.begin(), cells_1d_manual.end());
329  const Calculus2D::DualIdentity1 reorder_1cellp_2d = calculus_2d_manual.reorder<1, DUAL>(cells_2d_manual.begin(), cells_2d_manual.end());
330  const Calculus3D::DualIdentity1 reorder_1cellp_3d = calculus_3d_manual.reorder<1, DUAL>(cells_3d_manual.begin(), cells_3d_manual.end());
331 
332  { // testing primal operators
333  const Calculus1D::PrimalIdentity0 primal_laplace_1d = reorder_0cell_1d * calculus_1d_manual.laplace<PRIMAL>() * reorder_0cell_1d.transpose();
334  const Calculus2D::PrimalIdentity0 primal_laplace_2d = reorder_0cell_2d * calculus_2d_manual.laplace<PRIMAL>() * reorder_0cell_2d.transpose();
335  const Calculus3D::PrimalIdentity0 primal_laplace_3d = reorder_0cell_3d * calculus_3d_manual.laplace<PRIMAL>() * reorder_0cell_3d.transpose();
336  trace.info() << "primal_laplace_1d=" << primal_laplace_1d << endl;
337  trace.info() << "primal_laplace_2d=" << primal_laplace_2d << endl;
338  trace.info() << "primal_laplace_3d=" << primal_laplace_3d << endl;
339  trace.info() << "primal_laplace_container=" << endl << MatrixXd(primal_laplace_1d.myContainer) << endl;
340  FATAL_ERROR( equal(
341  reorder_1cellp_1d * calculus_1d_manual.hodge<0, PRIMAL>() * reorder_0cell_1d.transpose(),
342  reorder_1cellp_2d * calculus_2d_manual.hodge<0, PRIMAL>() * reorder_0cell_2d.transpose()
343  ) );
344  FATAL_ERROR( equal(
345  reorder_1cellp_1d * calculus_1d_manual.hodge<0, PRIMAL>() * reorder_0cell_1d.transpose(),
346  reorder_1cellp_3d * calculus_3d_manual.hodge<0, PRIMAL>() * reorder_0cell_3d.transpose()
347  ) );
348  FATAL_ERROR( equal(
349  reorder_0cellp_1d * calculus_1d_manual.hodge<1, PRIMAL>() * reorder_1cell_1d.transpose(),
350  reorder_0cellp_2d * calculus_2d_manual.hodge<1, PRIMAL>() * reorder_1cell_2d.transpose()
351  ) );
352  FATAL_ERROR( equal(
353  reorder_0cellp_1d * calculus_1d_manual.hodge<1, PRIMAL>() * reorder_1cell_1d.transpose(),
354  reorder_0cellp_3d * calculus_3d_manual.hodge<1, PRIMAL>() * reorder_1cell_3d.transpose()
355  ) );
356  FATAL_ERROR( equal(
357  reorder_1cell_1d * calculus_1d_manual.derivative<0, PRIMAL>() * reorder_0cell_1d.transpose(),
358  reorder_1cell_2d * calculus_2d_manual.derivative<0, PRIMAL>() * reorder_0cell_2d.transpose()
359  ) );
360  FATAL_ERROR( equal(
361  reorder_1cell_1d * calculus_1d_manual.derivative<0, PRIMAL>() * reorder_0cell_1d.transpose(),
362  reorder_1cell_3d * calculus_3d_manual.derivative<0, PRIMAL>() * reorder_0cell_3d.transpose()
363  ) );
364  FATAL_ERROR( equal(primal_laplace_1d, primal_laplace_2d) );
365  FATAL_ERROR( equal(primal_laplace_1d, primal_laplace_3d) );
366  }
367 
368  { // testing dual operators
369  const Calculus1D::DualIdentity0 dual_laplace_1d = reorder_0cellp_1d * calculus_1d_manual.laplace<DUAL>() * reorder_0cellp_1d.transpose();
370  const Calculus2D::DualIdentity0 dual_laplace_2d = reorder_0cellp_2d * calculus_2d_manual.laplace<DUAL>() * reorder_0cellp_2d.transpose();
371  const Calculus3D::DualIdentity0 dual_laplace_3d = reorder_0cellp_3d * calculus_3d_manual.laplace<DUAL>() * reorder_0cellp_3d.transpose();
372  trace.info() << "dual_laplace_1d=" << dual_laplace_1d << endl;
373  trace.info() << "dual_laplace_2d=" << dual_laplace_2d << endl;
374  trace.info() << "dual_laplace_3d=" << dual_laplace_3d << endl;
375  trace.info() << "dual_laplace_container=" << endl << MatrixXd(dual_laplace_1d.myContainer) << endl;
376  FATAL_ERROR( equal(
377  reorder_1cell_1d * calculus_1d_manual.hodge<0, DUAL>() * reorder_0cellp_1d.transpose(),
378  reorder_1cell_2d * calculus_2d_manual.hodge<0, DUAL>() * reorder_0cellp_2d.transpose()
379  ) );
380  FATAL_ERROR( equal(
381  reorder_1cell_1d * calculus_1d_manual.hodge<0, DUAL>() * reorder_0cellp_1d.transpose(),
382  reorder_1cell_3d * calculus_3d_manual.hodge<0, DUAL>() * reorder_0cellp_3d.transpose()
383  ) );
384  FATAL_ERROR( equal(
385  reorder_0cell_1d * calculus_1d_manual.hodge<1, DUAL>() * reorder_1cellp_1d.transpose(),
386  reorder_0cell_2d * calculus_2d_manual.hodge<1, DUAL>() * reorder_1cellp_2d.transpose()
387  ) );
388  FATAL_ERROR( equal(
389  reorder_0cell_1d * calculus_1d_manual.hodge<1, DUAL>() * reorder_1cellp_1d.transpose(),
390  reorder_0cell_3d * calculus_3d_manual.hodge<1, DUAL>() * reorder_1cellp_3d.transpose()
391  ) );
392  FATAL_ERROR( equal(
393  reorder_1cellp_1d * calculus_1d_manual.derivative<0, DUAL>() * reorder_0cellp_1d.transpose(),
394  reorder_1cellp_2d * calculus_2d_manual.derivative<0, DUAL>() * reorder_0cellp_2d.transpose()
395  ) );
396  FATAL_ERROR( equal(
397  reorder_1cellp_1d * calculus_1d_manual.derivative<0, DUAL>() * reorder_0cellp_1d.transpose(),
398  reorder_1cellp_3d * calculus_3d_manual.derivative<0, DUAL>() * reorder_0cellp_3d.transpose()
399  ) );
400  FATAL_ERROR( equal(dual_laplace_1d, dual_laplace_2d) );
401  FATAL_ERROR( equal(dual_laplace_1d, dual_laplace_3d) );
402  }
403 
404  { // checking dual laplace factory calculus vs manual calculus
405  const Calculus1D::DualIdentity0 reorder_0cellp_1d_factory = calculus_1d_factory.reorder<0, DUAL>(cells_1d_manual.begin(), cells_1d_manual.end());
406  const Calculus2D::DualIdentity0 reorder_0cellp_2d_factory = calculus_2d_factory.reorder<0, DUAL>(cells_2d_manual.begin(), cells_2d_manual.end());
407  const Calculus3D::DualIdentity0 reorder_0cellp_3d_factory = calculus_3d_factory.reorder<0, DUAL>(cells_3d_manual.begin(), cells_3d_manual.end());
408  FATAL_ERROR( equal(
409  reorder_0cellp_1d * calculus_1d_manual.laplace<DUAL>() * reorder_0cellp_1d.transpose(),
410  reorder_0cellp_1d_factory * calculus_1d_factory.laplace<DUAL>() * reorder_0cellp_1d_factory.transpose()
411  ) );
412  FATAL_ERROR( equal(
413  reorder_0cellp_2d * calculus_2d_manual.laplace<DUAL>() * reorder_0cellp_2d.transpose(),
414  reorder_0cellp_2d_factory * calculus_2d_factory.laplace<DUAL>() * reorder_0cellp_2d_factory.transpose()
415  ) );
416  FATAL_ERROR( equal(
417  reorder_0cellp_3d * calculus_3d_manual.laplace<DUAL>() * reorder_0cellp_3d.transpose(),
418  reorder_0cellp_3d_factory * calculus_3d_factory.laplace<DUAL>() * reorder_0cellp_3d_factory.transpose()
419  ) );
420  }
421 
422  trace.endBlock();
423 
424  trace.beginBlock("checking border");
425 
426  { // 2d ambient border
427  const Calculus2D calculus_2d_factory_no_border = CalculusFactory::createFromNSCells<1>(ncells_2d_factory.begin(), ncells_2d_factory.end(), false);
428  trace.info() << "calculus_2d_factory_no_border=" << calculus_2d_factory_no_border << endl;
429  trace.info() << "calculus_2d_factory=" << calculus_2d_factory << endl;
430  FATAL_ERROR( calculus_2d_factory.containsCell(calculus_2d_factory.myKSpace.uCell(Z2i::Point(6,0))) );
431  FATAL_ERROR( !calculus_2d_factory_no_border.containsCell(calculus_2d_factory_no_border.myKSpace.uCell(Z2i::Point(6,0))) );
432  FATAL_ERROR( calculus_2d_factory.containsCell(calculus_2d_factory.myKSpace.uCell(Z2i::Point(0,0))) );
433  FATAL_ERROR( !calculus_2d_factory_no_border.containsCell(calculus_2d_factory_no_border.myKSpace.uCell(Z2i::Point(0,0))) );
434  }
435 
436  { // 3d ambient border
437  const Calculus3D calculus_3d_factory_no_border = CalculusFactory::createFromNSCells<1>(ncells_3d_factory.begin(), ncells_3d_factory.end(), false);
438  trace.info() << "calculus_3d_factory_no_border=" << calculus_3d_factory_no_border << endl;
439  trace.info() << "calculus_3d_factory=" << calculus_3d_factory << endl;
440  FATAL_ERROR( calculus_3d_factory.containsCell(calculus_3d_factory.myKSpace.uCell(Z3i::Point(2,2,-2))) );
441  FATAL_ERROR( !calculus_3d_factory_no_border.containsCell(calculus_3d_factory_no_border.myKSpace.uCell(Z3i::Point(2,2,-2))) );
442  FATAL_ERROR( calculus_3d_factory.containsCell(calculus_3d_factory.myKSpace.uCell(Z3i::Point(0,0,0))) );
443  FATAL_ERROR( !calculus_3d_factory_no_border.containsCell(calculus_3d_factory_no_border.myKSpace.uCell(Z3i::Point(0,0,0))) );
444  }
445 
446  trace.endBlock();
447 
448  trace.beginBlock("checking area");
449 
450  {
451  const Calculus2D::Scalar area_th = calculus_2d_factory.kFormLength(0, DUAL);
452  const Calculus2D::Scalar area_primal = (
453  calculus_2d_factory.hodge<0, PRIMAL>() *
454  Calculus2D::PrimalForm0::ones(calculus_2d_factory)
455  ).myContainer.array().sum();
456  const Calculus2D::Scalar area_dual = (
457  calculus_2d_factory.hodge<0, DUAL>() *
458  Calculus2D::DualForm0::ones(calculus_2d_factory)
459  ).myContainer.array().sum();
460  trace.info() << "area_2d_th=" << area_th << endl;
461  trace.info() << "area_2d_primal=" << area_primal << endl;
462  trace.info() << "area_2d_dual=" << area_dual << endl;
463  FATAL_ERROR( area_th == area_primal );
464  FATAL_ERROR( area_th == area_dual );
465  }
466 
467  {
468  const Calculus3D::Scalar area_th = calculus_3d_factory.kFormLength(0, DUAL);
469  const Calculus3D::Scalar area_primal = (
470  calculus_3d_factory.hodge<0, PRIMAL>() *
471  Calculus3D::PrimalForm0::ones(calculus_3d_factory)
472  ).myContainer.array().sum();
473  const Calculus3D::Scalar area_dual = (
474  calculus_3d_factory.hodge<0, DUAL>() *
475  Calculus3D::DualForm0::ones(calculus_3d_factory)
476  ).myContainer.array().sum();
477  trace.info() << "area_3d_th=" << area_th << endl;
478  trace.info() << "area_3d_primal=" << area_primal << endl;
479  trace.info() << "area_3d_dual=" << area_dual << endl;
480  FATAL_ERROR( area_th == area_primal );
481  FATAL_ERROR( area_th == area_dual );
482  }
483 
484  trace.endBlock();
485 
486  trace.endBlock();
487  }
488 
489  /*{
490  trace.beginBlock("2d manifold embedding");
491 
493  typedef DiscreteExteriorCalculus<2, 2, EigenLinearAlgebraBackend> Calculus2D;
494  typedef DiscreteExteriorCalculus<2, 3, EigenLinearAlgebraBackend> Calculus3D;
496 
497  Calculus2D calculus_2d_manual;
498  {
499  for (int xx=0; xx<=4; xx++)
500  for (int yy=0; yy<=4; yy++)
501  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(xx,yy)) );
502 
503  for (int xx=0; xx<=4; xx++)
504  for (int yy=0; yy<=4; yy++)
505  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(xx,yy+4)) );
506 
507  for (int xx=0; xx<=4; xx++)
508  for (int yy=0; yy<=4; yy++)
509  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(xx,yy+8)) );
510 
511  for (int xx=0; xx<=4; xx++)
512  for (int yy=0; yy<=4; yy++)
513  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(xx,yy+12)) );
514 
515  for (int xx=0; xx<=4; xx++)
516  for (int yy=0; yy<=4; yy++)
517  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(xx,yy+16)) );
518 
519  for (int xx=0; xx<=4; xx++)
520  for (int yy=0; yy<=4; yy++)
521  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(xx+4,yy)) );
522 
523  for (int xx=0; xx<=4; xx++)
524  for (int yy=0; yy<=4; yy++)
525  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(xx+8,yy)) );
526 
527  for (int xx=0; xx<=4; xx++)
528  for (int yy=0; yy<=4; yy++)
529  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(xx+12,yy)) );
530 
531  for (int xx=0; xx<=4; xx++)
532  for (int yy=0; yy<=4; yy++)
533  calculus_2d_manual.insertSCell( calculus_2d_manual.myKSpace.sCell(Z2i::Point(xx+16,yy)) );
534 
535  calculus_2d_manual.updateIndexes();
536  }
537  trace.info() << "calculus_2d_manual=" << calculus_2d_manual << endl;
538 
539  {
540  Board2D board;
541  board << Z2i::Domain(Z2i::Point(-1,-1), Z2i::Point(10,10));
542  board << calculus_2d_manual;
543  board.saveSVG("embedding_2d_calculus_2d.svg");
544  }
545 
547  typedef std::list<Calculus2D::SCell> SCells2D;
548  SCells2D ncells_2d_factory;
550  for (int kk=0; kk<calculus_2d_manual.kFormLength(2, PRIMAL); kk++) ncells_2d_factory.push_back( calculus_2d_manual.getSCell(2, PRIMAL, kk) );
552  const Calculus2D calculus_2d_factory_weighed = CalculusFactory::createFromNSCells<2>(ncells_2d_factory.begin(), ncells_2d_factory.end(), true);
554  Calculus2D calculus_2d_factory = calculus_2d_factory_weighed;
555  calculus_2d_factory.resetSizes();
556  trace.info() << "calculus_2d_factory=" << calculus_2d_factory << endl;
557 
558  SCells2D cells_2d_manual;
559  {
560  std::set<Calculus2D::SCell> already_inserted;
561 
562  for (int xx=0; xx<=4; xx++)
563  for (int yy=0; yy<=4; yy++)
564  {
565  const Calculus2D::SCell cell = calculus_2d_manual.myKSpace.sCell(Z2i::Point(xx,yy));
566  if (already_inserted.find(cell) != already_inserted.end()) continue;
567  already_inserted.insert(cell);
568  cells_2d_manual.push_back(cell);
569  }
570 
571  for (int xx=0; xx<=4; xx++)
572  for (int yy=0; yy<=4; yy++)
573  {
574  const Calculus2D::SCell cell = calculus_2d_manual.myKSpace.sCell(Z2i::Point(xx,yy+4));
575  if (already_inserted.find(cell) != already_inserted.end()) continue;
576  already_inserted.insert(cell);
577  cells_2d_manual.push_back(cell);
578  }
579 
580  for (int xx=0; xx<=4; xx++)
581  for (int yy=0; yy<=4; yy++)
582  {
583  const Calculus2D::SCell cell = calculus_2d_manual.myKSpace.sCell(Z2i::Point(xx,yy+8));
584  if (already_inserted.find(cell) != already_inserted.end()) continue;
585  already_inserted.insert(cell);
586  cells_2d_manual.push_back(cell);
587  }
588 
589  for (int xx=0; xx<=4; xx++)
590  for (int yy=0; yy<=4; yy++)
591  {
592  const Calculus2D::SCell cell = calculus_2d_manual.myKSpace.sCell(Z2i::Point(xx,yy+12));
593  if (already_inserted.find(cell) != already_inserted.end()) continue;
594  already_inserted.insert(cell);
595  cells_2d_manual.push_back(cell);
596  }
597 
598  for (int xx=0; xx<=4; xx++)
599  for (int yy=0; yy<=4; yy++)
600  {
601  const Calculus2D::SCell cell = calculus_2d_manual.myKSpace.sCell(Z2i::Point(xx,yy+16));
602  if (already_inserted.find(cell) != already_inserted.end()) continue;
603  already_inserted.insert(cell);
604  cells_2d_manual.push_back(cell);
605  }
606 
607  for (int xx=0; xx<=4; xx++)
608  for (int yy=0; yy<=4; yy++)
609  {
610  const Calculus2D::SCell cell = calculus_2d_manual.myKSpace.sCell(Z2i::Point(xx+4,yy));
611  if (already_inserted.find(cell) != already_inserted.end()) continue;
612  already_inserted.insert(cell);
613  cells_2d_manual.push_back(cell);
614  }
615 
616  for (int xx=0; xx<=4; xx++)
617  for (int yy=0; yy<=4; yy++)
618  {
619  const Calculus2D::SCell cell = calculus_2d_manual.myKSpace.sCell(Z2i::Point(xx+8,yy));
620  if (already_inserted.find(cell) != already_inserted.end()) continue;
621  already_inserted.insert(cell);
622  cells_2d_manual.push_back(cell);
623  }
624 
625  for (int xx=0; xx<=4; xx++)
626  for (int yy=0; yy<=4; yy++)
627  {
628  const Calculus2D::SCell cell = calculus_2d_manual.myKSpace.sCell(Z2i::Point(xx+12,yy));
629  if (already_inserted.find(cell) != already_inserted.end()) continue;
630  already_inserted.insert(cell);
631  cells_2d_manual.push_back(cell);
632  }
633 
634  for (int xx=0; xx<=4; xx++)
635  for (int yy=0; yy<=4; yy++)
636  {
637  const Calculus2D::SCell cell = calculus_2d_manual.myKSpace.sCell(Z2i::Point(xx+16,yy));
638  if (already_inserted.find(cell) != already_inserted.end()) continue;
639  already_inserted.insert(cell);
640  cells_2d_manual.push_back(cell);
641  }
642  }
643 
644  Calculus3D calculus_3d_manual;
645  {
646  for (int xx=0; xx<=4; xx++)
647  for (int yy=0; yy<=4; yy++)
648  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(xx,yy,0)) );
649 
650  for (int xx=0; xx<=4; xx++)
651  for (int yy=0; yy<=4; yy++)
652  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(xx,4,yy)) );
653 
654  for (int xx=0; xx<=4; xx++)
655  for (int yy=0; yy<=4; yy++)
656  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(xx,4-yy,4),
657  xx%2 != 0 && yy%2 != 0 ? Calculus3D::KSpace::NEG : // surfels
658  yy%2 != 0 ? Calculus3D::KSpace::NEG : // y-edge
659  Calculus3D::KSpace::POS) );
660 
661  for (int xx=0; xx<=4; xx++)
662  for (int yy=0; yy<=4; yy++)
663  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(xx,-yy,4),
664  xx%2 != 0 && yy%2 != 0 ? Calculus3D::KSpace::NEG : // surfels
665  yy%2 != 0 ? Calculus3D::KSpace::NEG : // y-edge
666  Calculus3D::KSpace::POS) );
667 
668  for (int xx=0; xx<=4; xx++)
669  for (int yy=0; yy<=4; yy++)
670  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(xx,-4,4-yy),
671  xx%2 != 0 && yy%2 != 0 ? Calculus3D::KSpace::NEG : // surfels
672  yy%2 != 0 ? Calculus3D::KSpace::NEG : // y-edge
673  Calculus3D::KSpace::POS) );
674 
675  for (int xx=0; xx<=4; xx++)
676  for (int yy=0; yy<=4; yy++)
677  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(4,yy,-xx),
678  xx%2 != 0 && yy%2 != 0 ? Calculus3D::KSpace::POS : // surfels
679  xx%2 != 0 ? Calculus3D::KSpace::NEG : // x-edge
680  Calculus3D::KSpace::POS) );
681 
682  for (int xx=0; xx<=4; xx++)
683  for (int yy=0; yy<=4; yy++)
684  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(4-xx,yy,-4),
685  xx%2 != 0 && yy%2 != 0 ? Calculus3D::KSpace::NEG : // surfels
686  xx%2 != 0 ? Calculus3D::KSpace::NEG : // x-edge
687  Calculus3D::KSpace::POS) );
688 
689  for (int xx=0; xx<=4; xx++)
690  for (int yy=0; yy<=4; yy++)
691  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(-xx,yy,-4),
692  xx%2 != 0 && yy%2 != 0 ? Calculus3D::KSpace::NEG : // surfels
693  xx%2 != 0 ? Calculus3D::KSpace::NEG : // x-edge
694  Calculus3D::KSpace::POS) );
695 
696  for (int xx=0; xx<=4; xx++)
697  for (int yy=0; yy<=4; yy++)
698  calculus_3d_manual.insertSCell( calculus_3d_manual.myKSpace.sCell(Z3i::Point(-4,yy,-4+xx),
699  xx%2 != 0 && yy%2 != 0 ? Calculus3D::KSpace::NEG : // surfels
700  xx%2 != 0 ? Calculus3D::KSpace::POS : // x-edge
701  Calculus3D::KSpace::POS) );
702 
703  calculus_3d_manual.updateIndexes();
704  }
705  trace.info() << "calculus_3d_manual=" << calculus_3d_manual << endl;
706 
707 #if !defined(NOVIEWER)
708  Display3DFactory<Calculus3D::KSpace::Space, Calculus3D::KSpace>::draw(viewer2, calculus_3d_manual);
709  viewer2 << Viewer::updateDisplay;
710 #endif
711 
713  typedef std::list<Calculus3D::SCell> SCells3D;
714  SCells3D ncells_3d_factory;
716  for (int kk=0; kk<calculus_3d_manual.kFormLength(2, PRIMAL); kk++) ncells_3d_factory.push_back( calculus_3d_manual.getSCell(2, PRIMAL, kk) );
718  const Calculus3D calculus_3d_factory_weighed = CalculusFactory::createFromNSCells<2>(ncells_3d_factory.begin(), ncells_3d_factory.end(), true);
720  Calculus3D calculus_3d_factory = calculus_3d_factory_weighed;
721  calculus_3d_factory.resetSizes();
722  trace.info() << "calculus_3d_factory=" << calculus_3d_factory << endl;
723 
724  SCells3D cells_3d_manual;
725  {
726  std::set<Calculus3D::SCell> already_inserted;
727 
728  for (int xx=0; xx<=4; xx++)
729  for (int yy=0; yy<=4; yy++)
730  {
731  const Calculus3D::SCell cell = calculus_3d_manual.myKSpace.sCell(Z3i::Point(xx,yy,0));
732  if (already_inserted.find(cell) != already_inserted.end()) continue;
733  already_inserted.insert(cell);
734  cells_3d_manual.push_back(cell);
735  }
736 
737  for (int xx=0; xx<=4; xx++)
738  for (int yy=0; yy<=4; yy++)
739  {
740  const Calculus3D::SCell cell = calculus_3d_manual.myKSpace.sCell(Z3i::Point(xx,4,yy));
741  if (already_inserted.find(cell) != already_inserted.end()) continue;
742  already_inserted.insert(cell);
743  cells_3d_manual.push_back(cell);
744  }
745 
746  for (int xx=0; xx<=4; xx++)
747  for (int yy=0; yy<=4; yy++)
748  {
749  const Calculus3D::SCell cell = calculus_3d_manual.myKSpace.sCell(Z3i::Point(xx,4-yy,4),
750  xx%2 != 0 && yy%2 != 0 ? Calculus3D::KSpace::NEG : // surfels
751  yy%2 != 0 ? Calculus3D::KSpace::NEG : // y-edge
752  Calculus3D::KSpace::POS);
753  if (already_inserted.find(cell) != already_inserted.end()) continue;
754  already_inserted.insert(cell);
755  cells_3d_manual.push_back(cell);
756  }
757 
758  for (int xx=0; xx<=4; xx++)
759  for (int yy=0; yy<=4; yy++)
760  {
761  const Calculus3D::SCell cell = calculus_3d_manual.myKSpace.sCell(Z3i::Point(xx,-yy,4),
762  xx%2 != 0 && yy%2 != 0 ? Calculus3D::KSpace::NEG : // surfels
763  yy%2 != 0 ? Calculus3D::KSpace::NEG : // y-edge
764  Calculus3D::KSpace::POS);
765  if (already_inserted.find(cell) != already_inserted.end()) continue;
766  already_inserted.insert(cell);
767  cells_3d_manual.push_back(cell);
768  }
769 
770  for (int xx=0; xx<=4; xx++)
771  for (int yy=0; yy<=4; yy++)
772  {
773  const Calculus3D::SCell cell = calculus_3d_manual.myKSpace.sCell(Z3i::Point(xx,-4,4-yy),
774  xx%2 != 0 && yy%2 != 0 ? Calculus3D::KSpace::NEG : // surfels
775  yy%2 != 0 ? Calculus3D::KSpace::NEG : // y-edge
776  Calculus3D::KSpace::POS);
777  if (already_inserted.find(cell) != already_inserted.end()) continue;
778  already_inserted.insert(cell);
779  cells_3d_manual.push_back(cell);
780  }
781 
782  for (int xx=0; xx<=4; xx++)
783  for (int yy=0; yy<=4; yy++)
784  {
785  const Calculus3D::SCell cell = calculus_3d_manual.myKSpace.sCell(Z3i::Point(4,yy,-xx),
786  xx%2 != 0 && yy%2 != 0 ? Calculus3D::KSpace::POS : // surfels
787  xx%2 != 0 ? Calculus3D::KSpace::NEG : // x-edge
788  Calculus3D::KSpace::POS);
789  if (already_inserted.find(cell) != already_inserted.end()) continue;
790  already_inserted.insert(cell);
791  cells_3d_manual.push_back(cell);
792  }
793 
794  for (int xx=0; xx<=4; xx++)
795  for (int yy=0; yy<=4; yy++)
796  {
797  const Calculus3D::SCell cell = calculus_3d_manual.myKSpace.sCell(Z3i::Point(4-xx,yy,-4),
798  xx%2 != 0 && yy%2 != 0 ? Calculus3D::KSpace::NEG : // surfels
799  xx%2 != 0 ? Calculus3D::KSpace::NEG : // x-edge
800  Calculus3D::KSpace::POS);
801  if (already_inserted.find(cell) != already_inserted.end()) continue;
802  already_inserted.insert(cell);
803  cells_3d_manual.push_back(cell);
804  }
805 
806  for (int xx=0; xx<=4; xx++)
807  for (int yy=0; yy<=4; yy++)
808  {
809  const Calculus3D::SCell cell = calculus_3d_manual.myKSpace.sCell(Z3i::Point(-xx,yy,-4),
810  xx%2 != 0 && yy%2 != 0 ? Calculus3D::KSpace::NEG : // surfels
811  xx%2 != 0 ? Calculus3D::KSpace::NEG : // x-edge
812  Calculus3D::KSpace::POS);
813  if (already_inserted.find(cell) != already_inserted.end()) continue;
814  already_inserted.insert(cell);
815  cells_3d_manual.push_back(cell);
816  }
817 
818  for (int xx=0; xx<=4; xx++)
819  for (int yy=0; yy<=4; yy++)
820  {
821  const Calculus3D::SCell cell = calculus_3d_manual.myKSpace.sCell(Z3i::Point(-4,yy,-4+xx),
822  xx%2 != 0 && yy%2 != 0 ? Calculus3D::KSpace::NEG : // surfels
823  xx%2 != 0 ? Calculus3D::KSpace::POS : // x-edge
824  Calculus3D::KSpace::POS);
825  if (already_inserted.find(cell) != already_inserted.end()) continue;
826  already_inserted.insert(cell);
827  cells_3d_manual.push_back(cell);
828  }
829  }
830 
831  trace.beginBlock("checking operators");
832  trace.info() << calculus_3d_manual << endl;
833  trace.info() << cells_3d_manual.size() << endl;
834 
835  const Calculus2D::PrimalIdentity0 reorder_0cell_2d = calculus_2d_manual.reorder<0, PRIMAL>(cells_2d_manual.begin(), cells_2d_manual.end());
836  const Calculus3D::PrimalIdentity0 reorder_0cell_3d = calculus_3d_manual.reorder<0, PRIMAL>(cells_3d_manual.begin(), cells_3d_manual.end());
837  const Calculus2D::PrimalIdentity1 reorder_1cell_2d = calculus_2d_manual.reorder<1, PRIMAL>(cells_2d_manual.begin(), cells_2d_manual.end());
838  const Calculus3D::PrimalIdentity1 reorder_1cell_3d = calculus_3d_manual.reorder<1, PRIMAL>(cells_3d_manual.begin(), cells_3d_manual.end());
839  const Calculus2D::PrimalIdentity2 reorder_2cell_2d = calculus_2d_manual.reorder<2, PRIMAL>(cells_2d_manual.begin(), cells_2d_manual.end());
840  const Calculus3D::PrimalIdentity2 reorder_2cell_3d = calculus_3d_manual.reorder<2, PRIMAL>(cells_3d_manual.begin(), cells_3d_manual.end());
841  const Calculus2D::DualIdentity0 reorder_0cellp_2d = calculus_2d_manual.reorder<0, DUAL>(cells_2d_manual.begin(), cells_2d_manual.end());
842  const Calculus3D::DualIdentity0 reorder_0cellp_3d = calculus_3d_manual.reorder<0, DUAL>(cells_3d_manual.begin(), cells_3d_manual.end());
843  const Calculus2D::DualIdentity1 reorder_1cellp_2d = calculus_2d_manual.reorder<1, DUAL>(cells_2d_manual.begin(), cells_2d_manual.end());
844  const Calculus3D::DualIdentity1 reorder_1cellp_3d = calculus_3d_manual.reorder<1, DUAL>(cells_3d_manual.begin(), cells_3d_manual.end());
845  const Calculus2D::DualIdentity2 reorder_2cellp_2d = calculus_2d_manual.reorder<2, DUAL>(cells_2d_manual.begin(), cells_2d_manual.end());
846  const Calculus3D::DualIdentity2 reorder_2cellp_3d = calculus_3d_manual.reorder<2, DUAL>(cells_3d_manual.begin(), cells_3d_manual.end());
847 
848  { // check primal operators
849  const Calculus2D::PrimalIdentity0 primal_laplace_2d = reorder_0cell_2d * calculus_2d_manual.laplace<PRIMAL>() * reorder_0cell_2d.transpose();
850  const Calculus3D::PrimalIdentity0 primal_laplace_3d = reorder_0cell_3d * calculus_3d_manual.laplace<PRIMAL>() * reorder_0cell_3d.transpose();
851  trace.info() << "primal_laplace_2d=" << primal_laplace_2d << endl;
852  trace.info() << "primal_laplace_3d=" << primal_laplace_3d << endl;
853  trace.info() << "primal_laplace_container=" << endl << MatrixXd(primal_laplace_2d.myContainer) << endl;
854  FATAL_ERROR( equal(
855  reorder_2cellp_2d * calculus_2d_manual.hodge<0, PRIMAL>() * reorder_0cell_2d.transpose(),
856  reorder_2cellp_3d * calculus_3d_manual.hodge<0, PRIMAL>() * reorder_0cell_3d.transpose()
857  ) );
858  FATAL_ERROR( equal(
859  reorder_1cellp_2d * calculus_2d_manual.hodge<1, PRIMAL>() * reorder_1cell_2d.transpose(),
860  reorder_1cellp_3d * calculus_3d_manual.hodge<1, PRIMAL>() * reorder_1cell_3d.transpose()
861  ) );
862  FATAL_ERROR( equal(
863  reorder_0cellp_2d * calculus_2d_manual.hodge<2, PRIMAL>() * reorder_2cell_2d.transpose(),
864  reorder_0cellp_3d * calculus_3d_manual.hodge<2, PRIMAL>() * reorder_2cell_3d.transpose()
865  ) );
866  FATAL_ERROR( equal(
867  reorder_1cell_2d * calculus_2d_manual.derivative<0, PRIMAL>() * reorder_0cell_2d.transpose(),
868  reorder_1cell_3d * calculus_3d_manual.derivative<0, PRIMAL>() * reorder_0cell_3d.transpose()
869  ) );
870  FATAL_ERROR( equal(
871  reorder_2cell_2d * calculus_2d_manual.derivative<1, PRIMAL>() * reorder_1cell_2d.transpose(),
872  reorder_2cell_3d * calculus_3d_manual.derivative<1, PRIMAL>() * reorder_1cell_3d.transpose()
873  ) );
874  FATAL_ERROR( equal(primal_laplace_2d, primal_laplace_3d) );
875  }
876 
877  { // check dual operators
878  const Calculus2D::DualIdentity0 dual_laplace_2d = reorder_0cellp_2d * calculus_2d_manual.laplace<DUAL>() * reorder_0cellp_2d.transpose();
879  const Calculus3D::DualIdentity0 dual_laplace_3d = reorder_0cellp_3d * calculus_3d_manual.laplace<DUAL>() * reorder_0cellp_3d.transpose();
880  trace.info() << "dual_laplace_2d=" << dual_laplace_2d << endl;
881  trace.info() << "dual_laplace_3d=" << dual_laplace_3d << endl;
882  trace.info() << "dual_laplace_container=" << endl << MatrixXd(dual_laplace_2d.myContainer) << endl;
883  FATAL_ERROR( equal(
884  reorder_2cell_2d * calculus_2d_manual.hodge<0, DUAL>() * reorder_0cellp_2d.transpose(),
885  reorder_2cell_3d * calculus_3d_manual.hodge<0, DUAL>() * reorder_0cellp_3d.transpose()
886  ) );
887  FATAL_ERROR( equal(
888  reorder_1cell_2d * calculus_2d_manual.hodge<1, DUAL>() * reorder_1cellp_2d.transpose(),
889  reorder_1cell_3d * calculus_3d_manual.hodge<1, DUAL>() * reorder_1cellp_3d.transpose()
890  ) );
891  FATAL_ERROR( equal(
892  reorder_0cell_2d * calculus_2d_manual.hodge<2, DUAL>() * reorder_2cellp_2d.transpose(),
893  reorder_0cell_3d * calculus_3d_manual.hodge<2, DUAL>() * reorder_2cellp_3d.transpose()
894  ) );
895  FATAL_ERROR( equal(
896  reorder_1cellp_2d * calculus_2d_manual.derivative<0, DUAL>() * reorder_0cellp_2d.transpose(),
897  reorder_1cellp_3d * calculus_3d_manual.derivative<0, DUAL>() * reorder_0cellp_3d.transpose()
898  ) );
899  FATAL_ERROR( equal(
900  reorder_2cellp_2d * calculus_2d_manual.derivative<1, DUAL>() * reorder_1cellp_2d.transpose(),
901  reorder_2cellp_3d * calculus_3d_manual.derivative<1, DUAL>() * reorder_1cellp_3d.transpose()
902  ) );
903  FATAL_ERROR( equal(dual_laplace_2d, dual_laplace_3d) );
904  }
905 
906  { // checking dual laplace factory calculus vs manual calculus
907  const Calculus2D::DualIdentity0 reorder_0cellp_2d_factory = calculus_2d_factory.reorder<0, DUAL>(cells_2d_manual.begin(), cells_2d_manual.end());
908  const Calculus3D::DualIdentity0 reorder_0cellp_3d_factory = calculus_3d_factory.reorder<0, DUAL>(cells_3d_manual.begin(), cells_3d_manual.end());
909  FATAL_ERROR( equal(
910  reorder_0cellp_2d * calculus_2d_manual.laplace<DUAL>() * reorder_0cellp_2d.transpose(),
911  reorder_0cellp_2d_factory * calculus_2d_factory.laplace<DUAL>() * reorder_0cellp_2d_factory.transpose()
912  ) );
913  FATAL_ERROR( equal(
914  reorder_0cellp_3d * calculus_3d_manual.laplace<DUAL>() * reorder_0cellp_3d.transpose(),
915  reorder_0cellp_3d_factory * calculus_3d_factory.laplace<DUAL>() * reorder_0cellp_3d_factory.transpose()
916  ) );
917  }
918 
919  trace.endBlock();
920 
921  trace.beginBlock("checking border");
922 
923  { // 2d ambient border
925  const Calculus2D calculus_2d_factory_no_border = CalculusFactory::createFromNSCells<2>(ncells_2d_factory.begin(), ncells_2d_factory.end(), false);
927  trace.info() << "calculus_2d_factory_no_border=" << calculus_2d_factory_no_border << endl;
928  trace.info() << "calculus_2d_factory=" << calculus_2d_factory << endl;
929 
930  Board2D board;
931  board << Z2i::Domain(Z2i::Point(-1,-1), Z2i::Point(10,10));
932  board << calculus_2d_factory_no_border;
933  board.saveSVG("embedding_2d_calculus_2d_no_border.svg");
934  }
935 
936  { // 3d ambient border
938  const Calculus3D calculus_3d_factory_no_border = CalculusFactory::createFromNSCells<2>(ncells_3d_factory.begin(), ncells_3d_factory.end(), false);
940  trace.info() << "calculus_3d_factory_no_border=" << calculus_3d_factory_no_border << endl;
941  trace.info() << "calculus_3d_factory=" << calculus_3d_factory << endl;
942 
943 #if !defined(NOVIEWER)
944  Display3DFactory<Calculus3D::KSpace::Space, Calculus3D::KSpace>::draw(viewer3, calculus_3d_factory_no_border);
945  viewer3 << Viewer::updateDisplay;
946 #endif
947  }
948 
949  trace.endBlock();
950 
951  trace.beginBlock("checking area");
952 
953  {
954  const Calculus2D::Scalar area_th = calculus_2d_factory_weighed.kFormLength(0, DUAL);
955  const Calculus2D::Scalar area_primal = (
956  calculus_2d_factory_weighed.hodge<0, PRIMAL>() *
957  Calculus2D::PrimalForm0::ones(calculus_2d_factory_weighed)
958  ).myContainer.array().sum();
959  const Calculus2D::Scalar area_dual = (
960  calculus_2d_factory_weighed.hodge<0, DUAL>() *
961  Calculus2D::DualForm0::ones(calculus_2d_factory_weighed)
962  ).myContainer.array().sum();
963  trace.info() << "area_2d_th=" << area_th << endl;
964  trace.info() << "area_2d_primal=" << area_primal << endl;
965  trace.info() << "area_2d_dual=" << area_dual << endl;
966  FATAL_ERROR( area_th == area_primal );
967  FATAL_ERROR( area_th == area_dual );
968  }
969 
970  {
971  const Calculus3D::Scalar area_th = calculus_3d_factory_weighed.kFormLength(0, DUAL);
972  const Calculus3D::Scalar area_primal = (
973  calculus_3d_factory_weighed.hodge<0, PRIMAL>() *
974  Calculus3D::PrimalForm0::ones(calculus_3d_factory_weighed)
975  ).myContainer.array().sum();
976  const Calculus3D::Scalar area_dual = (
977  calculus_3d_factory_weighed.hodge<0, DUAL>() *
978  Calculus3D::DualForm0::ones(calculus_3d_factory_weighed)
979  ).myContainer.array().sum();
980  trace.info() << "area_3d_th=" << area_th << endl;
981  trace.info() << "area_3d_primal=" << area_primal << endl;
982  trace.info() << "area_3d_dual=" << area_dual << endl;
983  FATAL_ERROR( area_th == area_primal );
984  FATAL_ERROR( area_th == area_dual );
985  }
986 
987  trace.endBlock();
988 
989  trace.endBlock();
990  }*/
991 
992 #if !defined(NOVIEWER)
993  return app.exec();
994 #else
995  return 0;
996 #endif
997 }
998 
void beginBlock(const std::string &keyword="")
Aim: DiscreteExteriorCalculus represents a calculus in the dec package. This is the main structure in...
static void drawDECSignedKhalimskyCell(DGtal::Board2D &board, const DGtal::SignedKhalimskyCell< dim, TInteger > &cell)
Trace trace
Definition: Common.h:130
HyperRectDomain< Space > Domain
Definition: StdDefs.h:99
SCell sCell(const SPreCell &c) const
STL namespace.
double endBlock()
static void draw(Display3D< Space, KSpace > &display, const DGtal::DiscreteExteriorCalculus< dimEmbedded, dimAmbient, TLinearAlgebraBackend, TInteger > &calculus)
Aim: This class provides static members to create DEC structures from various other DGtal structures...
DGtal is the top-level namespace which contains all DGtal functions and types.
Board & setFillColor(const DGtal::Color &color)
Definition: Board.cpp:322
void saveSVG(const char *filename, PageSize size=Board::BoundingBox, double margin=10.0) const
Definition: Board.cpp:1012
std::ostream & info()
Represents a signed cell in a cellular grid space by its Khalimsky coordinates and a boolean value...
Structure representing an RGB triple with alpha component.
Definition: Color.h:66
Aim: This class is a model of CCellularGridSpaceND. It represents the cubical grid as a cell complex...
Aim: This class specializes a 'Board' class so as to display DGtal objects more naturally (with <<)...
Definition: Board2D.h:70