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"
9 #include "DGtal/io/viewers/Viewer3D.h"
11 #include "DGtal/io/boards/Board2D.h"
12 #include "DGtal/base/Common.h"
13 #include "DGtal/helpers/StdDefs.h"
15 using namespace DGtal;
17 using namespace Eigen;
19 template <
typename OperatorAA,
typename OperatorBB>
21 equal(
const OperatorAA& aa,
const OperatorBB& bb)
23 return MatrixXd(aa.myContainer) == MatrixXd(bb.myContainer);
26 int main(
int ,
char** )
28 #if !defined(NOVIEWER)
35 #if !defined(NOVIEWER)
36 QApplication app(argc, argv);
41 #if !defined(NOVIEWER)
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) );
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) );
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) );
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++)
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);
88 calculus_1d_manual.updateIndexes();
89 trace.
info() <<
"calculus_1d_manual=" << calculus_1d_manual << endl;
97 const Z2i::Point point(calculus_1d_manual.myKSpace.uKCoord(ii->first, 0), 0);
101 board.
saveSVG(
"embedding_1d_calculus_1d.svg");
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;
109 Calculus2D calculus_2d_manual;
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();
144 trace.
info() <<
"calculus_2d_manual=" << calculus_2d_manual << endl;
149 board << calculus_2d_manual;
150 board.
saveSVG(
"embedding_1d_calculus_2d.svg");
154 typedef std::list<Calculus2D::SCell> SCells2D;
155 SCells2D ncells_2d_factory;
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) );
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;
180 SCells2D cells_2d_manual;
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)) );
215 Calculus3D calculus_3d_manual;
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();
250 trace.
info() <<
"calculus_3d_manual=" << calculus_3d_manual << endl;
252 #if !defined(NOVIEWER)
258 typedef std::vector<Calculus3D::SCell> SCells3D;
259 SCells3D ncells_3d_factory;
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) );
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;
283 SCells3D cells_3d_manual;
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)) );
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());
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;
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()
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()
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()
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()
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()
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()
364 FATAL_ERROR( equal(primal_laplace_1d, primal_laplace_2d) );
365 FATAL_ERROR( equal(primal_laplace_1d, primal_laplace_3d) );
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;
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()
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()
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()
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()
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()
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()
400 FATAL_ERROR( equal(dual_laplace_1d, dual_laplace_2d) );
401 FATAL_ERROR( equal(dual_laplace_1d, dual_laplace_3d) );
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());
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()
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()
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()
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))) );
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))) );
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 );
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 );
992 #if !defined(NOVIEWER)
Aim: This class specializes a 'Board' class so as to display DGtal objects more naturally (with <<)....
Structure representing an RGB triple with alpha component.
Aim: This class provides static members to create DEC structures from various other DGtal structures.
Aim: DiscreteExteriorCalculus represents a calculus in the dec package. This is the main structure in...
Aim: This class is a model of CCellularGridSpaceND. It represents the cubical grid as a cell complex,...
SCell sCell(const SPreCell &c) const
From a signed cell, returns a signed cell lying into this Khalismky space.
void beginBlock(const std::string &keyword="")
Board & setFillColor(const DGtal::Color &color)
void saveSVG(const char *filename, PageSize size=Board::BoundingBox, double margin=10.0) const
MyDigitalSurface::ConstIterator ConstIterator
HyperRectDomain< Space > Domain
DGtal is the top-level namespace which contains all DGtal functions and types.
static void drawDECSignedKhalimskyCell(DGtal::Board2D &board, const DGtal::SignedKhalimskyCell< dim, TInteger > &cell)
static void draw(Display3D< Space, KSpace > &display, const DGtal::DiscreteExteriorCalculus< dimEmbedded, dimAmbient, TLinearAlgebraBackend, TInteger > &calculus)
Represents a signed cell in a cellular grid space by its Khalimsky coordinates and a boolean value.
int main(int argc, char **argv)