The Original CHomP Software
cubmaps.h
Go to the documentation of this file.
1
3
15
16// Copyright (C) 1997-2020 by Pawel Pilarczyk.
17//
18// This file is part of my research software package. This is free software:
19// you can redistribute it and/or modify it under the terms of the GNU
20// General Public License as published by the Free Software Foundation,
21// either version 3 of the License, or (at your option) any later version.
22//
23// This software is distributed in the hope that it will be useful,
24// but WITHOUT ANY WARRANTY; without even the implied warranty of
25// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
26// GNU General Public License for more details.
27//
28// You should have received a copy of the GNU General Public License
29// along with this software; see the file "license.txt". If not,
30// please, see <https://www.gnu.org/licenses/>.
31
32// Started in January 2002. Last revision: October 3, 2008.
33
34
35#ifndef _CHOMP_CUBES_CUBMAPS_H_
36#define _CHOMP_CUBES_CUBMAPS_H_
37
38#include "chomp/system/config.h"
39#include "chomp/system/textfile.h"
40#include "chomp/cubes/pointset.h"
41#include "chomp/struct/integer.h"
42#include "chomp/struct/hashsets.h"
43#include "chomp/homology/gcomplex.h"
44#include "chomp/cubes/pointbas.h"
45#include "chomp/cubes/cube.h"
46#include "chomp/cubes/cell.h"
47
48#include <iostream>
49#include <fstream>
50#include <cstdlib>
51
52namespace chomp {
53namespace homology {
54
55// --------------------------------------------------
56// ------------------- generators -------------------
57// --------------------------------------------------
58
61// Note: The main loop(s) of this function are copied from the function
62// "createprojection".
63template <class euclidom, class tCell>
64std::ostream &writegenerators (std::ostream &out, const chain<euclidom> *hom,
65 const chaincomplex<euclidom> &c, const gcomplex<tCell,euclidom> &g,
66 const int *level, int xdim, int format = 0)
67{
68 typedef typename tCell::CoordType coordType;
69 bool firstlist = true;
70 for (int d = 0; d <= c. dim (); ++ d)
71 {
72 if ((level && !level [d]) || hom [d]. empty ())
73 continue;
74 if (firstlist)
75 firstlist = false;
76 else if (format == 0)
77 out << '\n';
78 if (format == 0)
79 {
80 if (hom [d]. size () == 1)
81 {
82 out << "The generator of H_" << d <<
83 " follows:" << '\n';
84 }
85 else
86 {
87 out << "The " << hom [d]. size () <<
88 " generators of H_" << d <<
89 " follow:" << '\n';
90 }
91 }
92 const hashedset<tCell> &cset = g [d];
93 for (int_t i = 0; i < hom [d]. size (); ++ i)
94 {
95 if (format == 0)
96 {
97 if (hom [d]. size () > 1)
98 out << "generator " <<
99 (i + 1) << '\n';
100 }
101 const chain<euclidom> &lst =
102 c. gethomgen (d, hom [d]. num (i));
103 for (int_t j = 0; j < lst. size (); ++ j)
104 {
105 coordType left [tCell::MaxDim];
106 cset [lst. num (j)]. leftcoord (left);
107 coordType right [tCell::MaxDim];
108 cset [lst. num (j)]. rightcoord (right);
109 int projdim = 0;
110 for (int k = 0; k < xdim; ++ k)
111 {
112 if (left [k] != right [k])
113 ++ projdim;
114 }
115 if (projdim != d)
116 continue;
117 if (format == 0)
118 {
119 out << lst. coef (j) << " * " <<
120 tCell (left, right, xdim) <<
121 '\n';
122 }
123 else
124 {
125 for (int k = 0; k < xdim; ++ k)
126 out << left [k] << " ";
127 for (int k = 0; k < xdim; ++ k)
128 out << right [k] << " ";
129 out << lst. coef (j) << " " <<
130 (i + 1) << " " << d << '\n';
131 }
132 }
133 }
134 }
135 return out;
136} /* writegenerators */
137
138
139// --------------------------------------------------
140// ------------------- map images -------------------
141// --------------------------------------------------
142
146template <class euclidom, class tCell, class tCube>
147int_t createimages (mvcellmap<tCell,euclidom,tCube> &m,
148 const mvmap<tCube,tCube> &f1, const mvmap<tCube,tCube> &f2,
149 const hashedset<tCube> &dom1, const hashedset<tCube> &dom2)
150{
151 typedef typename tCell::CoordType coordType;
152 int_t countimages = 0;
153 coordType leftbound [tCell::MaxDim];
154 coordType rightbound [tCell::MaxDim];
155 coordType left [tCell::MaxDim];
156 coordType right [tCell::MaxDim];
157 for (int d = 0; d <= m. dim (); ++ d)
158 {
159 const hashedset<tCell> &cset = m. get (d);
160 if (cset. empty ())
161 continue;
162 const int spacedim = cset [0]. spacedim ();
163 for (int_t i = 0; i < cset. size (); ++ i)
164 {
165 cset [i]. leftcoord (left);
166 cset [i]. rightcoord (right);
167 for (int_t j = 0; j < spacedim; ++ j)
168 {
169 // compensate for space wrapping
170 if (right [j] < left [j])
171 right [j] = left [j] + 1;
172 // compute the bounds
173 leftbound [j] = static_cast<coordType>
174 (left [j] - (left [j] == right [j]));
175 rightbound [j] = static_cast<coordType>
176 (right [j] +
177 (left [j] == right [j]));
178 }
179 tRectangle<coordType> r (leftbound, rightbound,
180 spacedim);
181 const coordType *c;
182 while ((c = r. get ()) != NULL)
183 {
184 if (!tCube::PointBase::check (c, spacedim))
185 continue;
186 tCube q (c, spacedim);
187 if (dom1. check (q))
188 m. add (d, i, f1 (q));
189 if (dom2. check (q))
190 m. add (d, i, f2 (q));
191 }
192 countimages += m (cset [i]). size ();
193 }
194 }
195 return countimages;
196} /* createimages */
197
199template <class euclidom, class tCell, class tCube>
200inline int_t createimages (mvcellmap<tCell,euclidom,tCube> &m,
201 const mvmap<tCube,tCube> &f, const hashedset<tCube> &dom)
202{
203 mvmap<tCube,tCube> emptymap;
204 hashedset<tCube> emptyset;
205 return createimages (m, f, emptymap, dom, emptyset);
206} /* createimages */
207
210template <class euclidom, class tCell, class tCube>
211int_t createimages (mvcellmap<tCell,euclidom,tCube> &m,
212 const mvmap<tCube,tCube> &f1, const mvmap<tCube,tCube> &f2)
213{
214 const hashedset<tCube> &dom1 = f1. getdomain ();
215 const hashedset<tCube> &dom2 = f2. getdomain ();
216 return createimages (m, f1, f2, dom1, dom2);
217} /* createimages */
218
220template <class euclidom, class tCell, class tCube>
221inline int_t createimages (mvcellmap<tCell,euclidom,tCube> &m,
222 const mvmap<tCube,tCube> &f)
223{
224 mvmap<tCube,tCube> emptymap;
225 return createimages (m, f, emptymap);
226} /* createimages */
227
228
229// --------------------------------------------------
230// ------------------ projections -------------------
231// --------------------------------------------------
232
236template <class euclidom, class tCell>
237void createprojection (const gcomplex<tCell,euclidom> &Fcompl,
238 const gcomplex<tCell,euclidom> &Ycompl,
239 chainmap<euclidom> &cmap,
240 int offset, int outdim, int discarddim, int *level = NULL)
241{
242 typedef typename tCell::CoordType coordType;
243 // go through the list of all the dimensions which are of concern
244 for (int d = 0; d <= Ycompl. dim (); ++ d)
245 {
246 if ((!level || level [d]) && (Fcompl. dim () >= d))
247 {
248 // take sets of cells of this dimension
249 const hashedset<tCell> &Fset = Fcompl [d];
250 if (Fset. empty ())
251 continue;
252 const hashedset<tCell> &Yset = Ycompl [d];
253 if (Yset. empty ())
254 continue;
255
256 // go through the list of cells in Fcompl of dim. d
257 for (int_t i = 0; i < Fset. size (); ++ i)
258 {
259 // get this cell and its coordinates
260 const tCell &Fcell = Fset [i];
261 coordType left [tCell::MaxDim];
262 Fcell. leftcoord (left);
263 coordType right [tCell::MaxDim];
264 Fcell. rightcoord (right);
265
266 // check if this cell has no width in the
267 // directions that are discarded
268 register int j;
269 for (j = 0; j < offset; ++ j)
270 if (left [j] != right [j])
271 {
272 j = offset + 33;
273 break;
274 }
275 if (j > offset)
276 continue;
277 for (j = 0; j < discarddim; ++ j)
278 if (left [offset + outdim + j] !=
279 right [offset + outdim + j])
280 {
281 j = discarddim + 33;
282 break;
283 }
284 if (j > discarddim)
285 continue;
286
287 // create the projected cell
288 if (!(tCell::PointBase::check
289 (left + offset, outdim)))
290 continue;
291 if (!(tCell::PointBase::check
292 (right + offset, outdim)))
293 continue;
294 // tCell projected (left + offset,
295 // right + offset, outdim);
296 tCell projected (Fcell, offset, outdim);
297
298 // find its number in Y
299 int_t nr = Yset. getnumber (projected);
300
301 // if not found, discard it
302 if (nr < 0)
303 continue;
304
305 // add the pair to the projection map
306 euclidom e;
307 e = 1;
308 cmap. add (d, nr, i, e);
309 }
310 }
311 }
312 return;
313} /* createprojection */
314
319template <class euclidom, class tCell>
320void project (const gcomplex<tCell,euclidom> &c,
321 gcomplex<tCell,euclidom> &img,
322 const gcomplex<tCell,euclidom> &only,
323 int offset, int outdim, int discarddim, const int *level,
324 bool watchforimages)
325{
326 typedef typename tCell::CoordType coordType;
327
328 // go through the list of all the dimensions which are of concern
329 for (int d = 0; d <= c. dim (); ++ d)
330 {
331 if ((!level || level [d]) && (c. dim () >= d))
332 {
333 // take sets of cells of this dimension
334 const hashedset<tCell> &cset = c [d];
335 if (cset. empty ())
336 continue;
337
338 // go through the list of cells in c of dimension d
339 for (int_t i = 0; i < cset. size (); ++ i)
340 {
341 // get this cell and its coordinates
342 const tCell &thecell = cset [i];
343 coordType left [tCell::MaxDim];
344 thecell. leftcoord (left);
345 coordType right [tCell::MaxDim];
346 thecell. rightcoord (right);
347
348 // check if this cell has no width in the
349 // directions that are discarded
350 int j;
351 for (j = 0; j < offset; ++ j)
352 if (left [j] != right [j])
353 {
354 j = offset + 33;
355 break;
356 }
357 if (j > offset)
358 continue;
359 for (j = 0; j < discarddim; ++ j)
360 if (left [offset + outdim + j] !=
361 right [offset + outdim + j])
362 {
363 j = discarddim + 33;
364 break;
365 }
366 if (j > discarddim)
367 continue;
368
369 // add the projected cell to the complex
370 if (!(tCell::PointBase::check
371 (left + offset, outdim))
372 || !(tCell::PointBase::check
373 (right + offset, outdim)))
374 {
375 if (watchforimages)
376 throw "Inclusion undefined!";
377 else
378 continue;
379 }
380 tCell projected (left + offset,
381 right + offset, outdim);
382 if ((only. dim () >= 0) &&
383 !only. check (projected))
384 {
385 if (watchforimages)
386 throw "Inclusion undefined.";
387 else
388 continue;
389 }
390 img. add (projected);
391 }
392 }
393 }
394 return;
395} /* project */
396
397
398// --------------------------------------------------
399// ------------------- read cubes -------------------
400// --------------------------------------------------
401
405template <class tCube>
406std::istream &readdomain (std::istream &in, hashedset<tCube> &dom)
407{
408 ignorecomments (in);
409 if ((in. peek () != 'S') && (in. peek () != 's'))
410 {
411 mvmap<tCube,tCube> Fdummy;
412 return readdomain (in, dom, Fdummy);
413 }
414
415 while (in. peek () != EOF)
416 {
417 if (in. peek () == '[')
418 {
419 in. get ();
420 tCube q;
421 in >> q;
422 if (!in)
423 throw "Can't read the file.";
424 dom. add (q);
425 }
426 ignoreline (in);
427 ignorecomments (in);
428 }
429
430 return in;
431} /* readdomain */
432
436template <class tCube>
437std::istream &readimage (std::istream &in, hashedset<tCube> &img)
438{
439 ignorecomments (in);
440 if ((in. peek () != 'S') && (in. peek () != 's'))
441 {
442 mvmap<tCube,tCube> Fdummy;
443 return readimage (in, img, Fdummy);
444 }
445
446 typename tCube::CoordType left [tCube::MaxDim];
447 typename tCube::CoordType right [tCube::MaxDim];
448 while (in. peek () != EOF)
449 {
450 if (in. peek () == '[')
451 {
452 in. get ();
453 tCube dummy;
454 in >> dummy;
455 in >> dummy;
456 ignorecomments (in);
457 in. get ();
458 ignorecomments (in);
459 in. get ();
460 tCube q1, q2;
461 in >> q1 >> q2;
462 if (!in)
463 throw "Can't read the file.";
464 ignorecomments (in);
465 in. get ();
466 ignorecomments (in);
467 int dim = q1. dim ();
468 q1. coord (left);
469 q2. coord (right);
470 tRectangle<typename tCube::CoordType> r
471 (left, right, dim);
472 const typename tCube::CoordType *c;
473 while ((c = r. get ()) != NULL)
474 img. add (tCube (c, dim));
475 }
476 else
477 ignoreline (in);
478 ignorecomments (in);
479 }
480
481 return in;
482} /* readimage */
483
487template <class tCube>
488std::istream &readimage (std::istream &in, const hashedset<tCube> &dom,
489 hashedset<tCube> &img)
490{
491 ignorecomments (in);
492 // the general format: [x,y,z] -> {[a,b,c], [d,e,f]}
493 if ((in. peek () != 'S') && (in. peek () != 's'))
494 {
495 while (in. peek () != EOF)
496 {
497 if ((closingparenthesis (in. peek ()) == EOF) &&
498 !std::isdigit (in. peek ()))
499 {
500 ignoreline (in);
501 ignorecomments (in);
502 continue;
503 }
504 tCube q;
505 in >> q;
506 bool ignore = !dom. check (q);
507 ignorecomments (in);
508 while (in. peek () == '-')
509 in. get ();
510 in. get (); // '>'
511 ignorecomments (in);
512 int opening = in. get ();
513 int closing = closingparenthesis (opening);
514 if (closing == EOF)
515 throw "An opening brace '{' expected.";
516 while ((in. peek () != EOF) &&
517 (in. peek () != closing))
518 {
519 if (ignore)
520 {
521 if (in. get () == opening)
522 {
523 while ((in. peek () != EOF) &&
524 (in. peek () != closing))
525 {
526 in. get ();
527 }
528 in. get (); // '}'
529 }
530 }
531 else
532 {
533 in >> q;
534 img. add (q);
535 ignorecomments (in);
536 if (in. peek () == ',')
537 {
538 in. get ();
539 ignorecomments (in);
540 }
541 }
542 }
543 // in. get (); // '}'
544 ignoreline (in);
545 ignorecomments (in);
546 }
547 }
548
549 typename tCube::CoordType left [tCube::MaxDim];
550 typename tCube::CoordType right [tCube::MaxDim];
551 while (in. peek () != EOF)
552 {
553 if (in. peek () == '[')
554 {
555 in. get ();
556 tCube domcube1, domcube2;
557 in >> domcube1;
558 if (!dom. check (domcube1))
559 {
560 ignoreline (in);
561 ignorecomments (in);
562 continue;
563 }
564 in >> domcube2;
565 ignorecomments (in);
566 in. get ();
567 ignorecomments (in);
568 in. get ();
569 tCube q1, q2;
570 in >> q1 >> q2;
571 if (!in)
572 throw "Can't read the file.";
573 ignorecomments (in);
574 in. get ();
575 ignorecomments (in);
576 if (dom. check (domcube1))
577 {
578 int dim = q1. dim ();
579 q1. coord (left);
580 q2. coord (right);
581 tRectangle<typename tCube::CoordType> r
582 (left, right, dim);
583 const typename tCube::CoordType *c;
584 while ((c = r. get ()) != NULL)
585 img. add (tCube (c, dim));
586 }
587 }
588 else
589 ignoreline (in);
590 ignorecomments (in);
591 }
592
593 return in;
594} /* readimage */
595
597template <class tCube>
598std::istream &readselective (std::istream &in, const hashedset<tCube> &dom1,
599 const hashedset<tCube> &dom2, mvmap<tCube,tCube> &m)
600{
601 if (dom1. empty () && dom2. empty ())
602 {
603 sout << "Warning: The domain of the map is empty.\n";
604 return in;
605 }
606
607 ignorecomments (in);
608 if ((in. peek () != 'S') && (in. peek () != 's'))
609 return readselective (in, m, dom1, dom2);
610
611 typename tCube::CoordType left [tCube::MaxDim];
612 typename tCube::CoordType right [tCube::MaxDim];
613 while (in. peek () != EOF)
614 {
615 if (in. peek () == '[')
616 {
617 in. get ();
618 tCube domcube;
619 in >> domcube;
620 int dim = domcube. dim ();
621 if (dom1. check (domcube) || dom2. check (domcube))
622 {
623 tCube q1, q2;
624 in >> q1; // (ignored)
625 ignorecomments (in);
626 in. get ();
627 ignorecomments (in);
628 in. get ();
629 in >> q1 >> q2;
630 if (!in)
631 throw "Can't read the file.";
632 hashedset<tCube> &img = m [domcube];
633 q1. coord (left);
634 q2. coord (right);
635 tRectangle<typename tCube::CoordType> r
636 (left, right, dim);
637 const typename tCube::CoordType *c;
638 while ((c = r. get ()) != NULL)
639 img. add (tCube (c, dim));
640 }
641 }
642 ignoreline (in);
643 ignorecomments (in);
644 }
645
646 return in;
647} /* readselective */
648
652template <class tCube>
653inline std::istream &readselective (std::istream &in,
654 const hashedset<tCube> &dom, mvmap<tCube,tCube> &m)
655{
656 hashedset<tCube> empty;
657 return readselective (in, dom, empty, m);
658} /* readselective */
659
660
661} // namespace homology
662} // namespace chomp
663
664#endif // _CHOMP_CUBES_CUBMAPS_H_
665
667
cell.h
This file includes header files with various definitions of cubical cells and defines the standard ty...
chomp::homology::chain
This class defines objects which represent chains as finite sequences of elements identified by integ...
Definition: chains.h:94
chomp::homology::chaincomplex
This is an implementation of a chain complex over an arbitrary ring.
Definition: chains.h:2662
chomp::homology::chainmap
This class defines a chain map between two chain complexes.
Definition: chains.h:3243
chomp::homology::gcomplex
The class that defines a geometric complex - a set of cells (cubes, simplices, etc).
Definition: gcomplex.h:85
chomp::homology::mvcellmap
This class represents a multivalued map whose domain is a geometric complex.
Definition: gcomplex.h:1530
chomp::homology::tRectangle
This class can be used for iterating a rectangular set of points, given its left and right bound.
Definition: pointset.h:1627
config.h
This file contains some precompiler definitions which indicate the operating system and/or compiler u...
cube.h
This file includes header files with various definitions of full cubes and defines the standard types...
gcomplex.h
This file contains a definition of a general geometric complex which represents a polyhedron.
int_t
int int_t
Index type for indexing arrays, counting cubes, etc.
Definition: config.h:115
hashsets.h
This file contains the definition of the container "hashedset" which can be used to represent a set o...
integer.h
This file defines a class "integer" which represents the ring of integers or the field of integers mo...
chomp::homology::writegenerators
std::ostream & writegenerators(std::ostream &out, const chain< euclidom > *hom, const chaincomplex< euclidom > &c, const gcomplex< tCell, euclidom > &g, const int *level, int xdim, int format=0)
Writes projected homology generators of a cubical complex to a file.
Definition: cubmaps.h:64
chomp::homology::ignoreline
void ignoreline(std::istream &in)
Ignores the input characters until the end of a line, including this end of the line.
Definition: textfile.h:376
chomp::homology::sout
outputstream sout
A replacement for standard output stream, with optional logging and other features provided by the cl...
chomp::homology::closingparenthesis
int closingparenthesis(int ch)
Returns the matching closing parenthesis for the given opening one or EOF if none.
Definition: textfile.h:411
chomp::homology::project
void project(const gcomplex< tCell, euclidom > &c, gcomplex< tCell, euclidom > &img, const gcomplex< tCell, euclidom > &only, int offset, int outdim, int discarddim, const int *level, bool watchforimages)
Creates the image of the projection from the set of cubical cells in the given geometric complex to t...
Definition: cubmaps.h:320
chomp::homology::createprojection
void createprojection(const gcomplex< tCell, euclidom > &Fcompl, const gcomplex< tCell, euclidom > &Ycompl, chainmap< euclidom > &cmap, int offset, int outdim, int discarddim, int *level=NULL)
Creates the chain map of the projection from a cell complex of the graph of a map to a cell complex o...
Definition: cubmaps.h:237
chomp::homology::readimage
std::istream & readimage(std::istream &in, hashedset< tCube > &img)
Reads the image of a multivalued cubical map.
Definition: cubmaps.h:437
chomp::homology::createimages
int_t createimages(mvcellmap< tCell, euclidom, tCube > &m, const mvmap< tCube, tCube > &f1, const mvmap< tCube, tCube > &f2, const hashedset< tCube > &dom1, const hashedset< tCube > &dom2)
Creates images of cells in 'm' as unions of cubes determined by f1 and f2.
Definition: cubmaps.h:147
chomp::homology::readselective
std::istream & readselective(std::istream &in, const hashedset< tCube > &dom1, const hashedset< tCube > &dom2, mvmap< tCube, tCube > &m)
Reads the restriction of a multivalued map to the given pair of sets.
Definition: cubmaps.h:598
chomp::homology::ignorecomments
void ignorecomments(std::istream &in)
Ignores white characters (spaces, tabulators, CRs and LFs), as well as comments from the input text f...
Definition: textfile.h:385
chomp::homology::readdomain
std::istream & readdomain(std::istream &in, hashedset< tCube > &dom)
Reads the domain of a multivalued cubical map.
Definition: cubmaps.h:406
chomp
This namespace contains the entire CHomP library interface.
Definition: bitmaps.h:51
pointbas.h
This file contains the definition of a point base class which is used for indexing all the points (n-...
pointset.h
This file contains the definition of a set of n-dimensional points with integer coordinates and sever...
textfile.h
This file contains some useful functions related to the text input/output procedures.
Generated by doxygen 1.9.4