diff --git a/docs/src/PolyhedralGeometry/intro.md b/docs/src/PolyhedralGeometry/intro.md
index f11eeb31b51e..e9e383063254 100644
--- a/docs/src/PolyhedralGeometry/intro.md
+++ b/docs/src/PolyhedralGeometry/intro.md
@@ -159,7 +159,13 @@ Some methods will require input or return output in form of an `IncidenceMatrix`
 IncidenceMatrix
 ```
 
-From the example it can be seen that this type supports `julia`'s matrix functionality. There are also functions to retrieve specific rows or columns as a `Set` over the non-zero indices.
+The unique nature of the `IncidenceMatrix` allows for different ways of construction:
+
+```@docs
+incidence_matrix
+```
+
+From the examples it can be seen that this type supports `julia`'s matrix functionality. There are also functions to retrieve specific rows or columns as a `Set` over the non-zero indices.
 
 ```@docs
 row(i::IncidenceMatrix, n::Int)
diff --git a/experimental/Schemes/src/ToricBlowups/constructors.jl b/experimental/Schemes/src/ToricBlowups/constructors.jl
index c9a192240382..0263c5d0028c 100644
--- a/experimental/Schemes/src/ToricBlowups/constructors.jl
+++ b/experimental/Schemes/src/ToricBlowups/constructors.jl
@@ -137,7 +137,7 @@ julia> rs = [1 1; -1 1]
   1  1
  -1  1
 
-julia> max_cones = IncidenceMatrix([[1, 2]])
+julia> max_cones = incidence_matrix([[1, 2]])
 1×2 IncidenceMatrix
 [1, 2]
 
diff --git a/src/AlgebraicGeometry/ToricVarieties/AlgebraicCycles/special_attributes.jl b/src/AlgebraicGeometry/ToricVarieties/AlgebraicCycles/special_attributes.jl
index 3f03add700af..c61df9438095 100644
--- a/src/AlgebraicGeometry/ToricVarieties/AlgebraicCycles/special_attributes.jl
+++ b/src/AlgebraicGeometry/ToricVarieties/AlgebraicCycles/special_attributes.jl
@@ -27,7 +27,7 @@ true
 julia> ngens(chow_ring(p2))
 3
 
-julia> v = normal_toric_variety(IncidenceMatrix([[1], [2], [3]]), [[1, 0], [0, 1], [-1, -1]])
+julia> v = normal_toric_variety(incidence_matrix([[1], [2], [3]]), [[1, 0], [0, 1], [-1, -1]])
 Normal toric variety
 
 julia> is_complete(v)
diff --git a/src/AlgebraicGeometry/ToricVarieties/CohomologyClasses/methods.jl b/src/AlgebraicGeometry/ToricVarieties/CohomologyClasses/methods.jl
index 23eca4350f2f..21478f09306b 100644
--- a/src/AlgebraicGeometry/ToricVarieties/CohomologyClasses/methods.jl
+++ b/src/AlgebraicGeometry/ToricVarieties/CohomologyClasses/methods.jl
@@ -63,7 +63,7 @@ julia> m = 2;
 
 julia> ray_generators = [e1, -e1, e2, e3, - e2 - e3 - m * e1];
 
-julia> max_cones = IncidenceMatrix([[1,3,4], [1,3,5], [1,4,5], [2,3,4], [2,3,5], [2,4,5]]);
+julia> max_cones = incidence_matrix([[1,3,4], [1,3,5], [1,4,5], [2,3,4], [2,3,5], [2,4,5]]);
 
 julia> X = normal_toric_variety(max_cones, ray_generators; non_redundant = true)
 Normal toric variety
diff --git a/src/AlgebraicGeometry/ToricVarieties/NormalToricVarieties/constructors.jl b/src/AlgebraicGeometry/ToricVarieties/NormalToricVarieties/constructors.jl
index ecce8f513743..507d2c082daf 100644
--- a/src/AlgebraicGeometry/ToricVarieties/NormalToricVarieties/constructors.jl
+++ b/src/AlgebraicGeometry/ToricVarieties/NormalToricVarieties/constructors.jl
@@ -82,7 +82,7 @@ julia> ray_generators = [[1,0], [0, 1], [-1, 5], [0, -1]]
  [-1, 5]
  [0, -1]
 
-julia> max_cones = IncidenceMatrix([[1, 2], [2, 3], [3, 4], [4, 1]])
+julia> max_cones = incidence_matrix([[1, 2], [2, 3], [3, 4], [4, 1]])
 4×4 IncidenceMatrix
 [1, 2]
 [2, 3]
diff --git a/src/PolyhedralGeometry/PolyhedralComplex/constructors.jl b/src/PolyhedralGeometry/PolyhedralComplex/constructors.jl
index 6b5217e01143..db9e9c56125d 100644
--- a/src/PolyhedralGeometry/PolyhedralComplex/constructors.jl
+++ b/src/PolyhedralGeometry/PolyhedralComplex/constructors.jl
@@ -40,7 +40,7 @@ points and the rows represent the polyhedra.
 
 # Examples
 ```jldoctest
-julia> IM = IncidenceMatrix([[1,2,3],[1,3,4]])
+julia> IM = incidence_matrix([[1,2,3],[1,3,4]])
 2×4 IncidenceMatrix
 [1, 2, 3]
 [1, 3, 4]
@@ -61,7 +61,7 @@ Polyhedral complex with rays and lineality:
 ```jldoctest
 julia> VR = [0 0 0; 1 0 0; 0 1 0; -1 0 0];
 
-julia> IM = IncidenceMatrix([[1,2,3],[1,3,4]]);
+julia> IM = incidence_matrix([[1,2,3],[1,3,4]]);
 
 julia> far_vertices = [2,3,4];
 
diff --git a/src/PolyhedralGeometry/PolyhedralComplex/properties.jl b/src/PolyhedralGeometry/PolyhedralComplex/properties.jl
index f769ef684cdb..670ae81cc9e6 100644
--- a/src/PolyhedralGeometry/PolyhedralComplex/properties.jl
+++ b/src/PolyhedralGeometry/PolyhedralComplex/properties.jl
@@ -5,7 +5,7 @@ Return the ambient dimension of `PC`.
 
 # Examples
 ```jldoctest
-julia> IM = IncidenceMatrix([[1,2,3],[1,3,4]])
+julia> IM = incidence_matrix([[1,2,3],[1,3,4]])
 2×4 IncidenceMatrix
 [1, 2, 3]
 [1, 3, 4]
@@ -41,7 +41,7 @@ Optional arguments for `as` include
 
 # Examples
 ```jldoctest
-julia> IM = IncidenceMatrix([[1,2,3],[1,3,4]]);
+julia> IM = incidence_matrix([[1,2,3],[1,3,4]]);
 
 julia> V = [0 0; 1 0; 1 1; 0 1];
 
@@ -63,7 +63,7 @@ julia> matrix(QQ, vertices(PointVector, PC))
 ```
 The following complex has no vertices:
 ```
-julia> IM = IncidenceMatrix([[1,2,3],[1,3,4]]);
+julia> IM = incidence_matrix([[1,2,3],[1,3,4]]);
 
 julia> VR = [0 0 0; 1 0 0; 0 1 0; -1 0 0];
 
@@ -139,7 +139,7 @@ function is mainly a helper function for [`maximal_polyhedra`](@ref maximal_poly
 
 # Examples
 ```jldoctest
-julia> IM = IncidenceMatrix([[1,2,3],[1,3,4]]);
+julia> IM = incidence_matrix([[1,2,3],[1,3,4]]);
 
 julia> VR = [0 0 0; 1 0 0; 0 1 0; -1 0 0];
 
@@ -193,7 +193,7 @@ See also [`rays`](@ref rays(PC::PolyhedralComplex{T}) where {T<:scalar_types}) a
 ```jldoctest
 julia> VR = [0 0 0; 1 0 0; 0 1 0; -1 0 0];
 
-julia> IM = IncidenceMatrix([[1,2,3],[1,3,4]]);
+julia> IM = incidence_matrix([[1,2,3],[1,3,4]]);
 
 julia> far_vertices = [2,3,4];
 
@@ -259,7 +259,7 @@ See also [`vertices`](@ref vertices(as::Type{PointVector{T}}, PC::PolyhedralComp
 ```jldoctest
 julia> VR = [0 0 0; 1 0 0; 0 1 0; -1 0 0];
 
-julia> IM = IncidenceMatrix([[1,2,3],[1,3,4]]);
+julia> IM = incidence_matrix([[1,2,3],[1,3,4]]);
 
 julia> far_vertices = [2,3,4];
 
@@ -317,7 +317,7 @@ Optional arguments for `as` include
 
 # Examples
 ```jldoctest
-julia> IM = IncidenceMatrix([[1,2,3],[1,3,4]]);
+julia> IM = incidence_matrix([[1,2,3],[1,3,4]]);
 
 julia> VR = [0 0; 1 0; 1 1; 0 1];
 
@@ -333,7 +333,7 @@ julia> matrix(QQ, rays(RayVector, PC))
 ```
 The following complex has no vertices:
 ```
-julia> IM = IncidenceMatrix([[1,2,3],[1,3,4]]);
+julia> IM = incidence_matrix([[1,2,3],[1,3,4]]);
 
 julia> VR = [0 0 0; 1 0 0; 0 1 0; -1 0 0];
 
@@ -398,7 +398,7 @@ refer to the output of [`vertices_and_rays`](@ref vertices_and_rays(PC::Polyhedr
 
 # Examples
 ```jldoctest
-julia> IM = IncidenceMatrix([[1,2,3],[1,3,4]])
+julia> IM = incidence_matrix([[1,2,3],[1,3,4]])
 2×4 IncidenceMatrix
 [1, 2, 3]
 [1, 3, 4]
@@ -435,7 +435,7 @@ Return the number of maximal polyhedra of `PC`
 
 # Examples
 ```jldoctest
-julia> IM = IncidenceMatrix([[1,2,3],[1,3,4]])
+julia> IM = incidence_matrix([[1,2,3],[1,3,4]])
 2×4 IncidenceMatrix
 [1, 2, 3]
 [1, 3, 4]
@@ -464,7 +464,7 @@ Determine whether the polyhedral complex is simplicial.
 
 # Examples
 ```jldoctest
-julia> IM = IncidenceMatrix([[1,2,3],[1,3,4]]);
+julia> IM = incidence_matrix([[1,2,3],[1,3,4]]);
 
 julia> VR = [0 0; 1 0; 1 1; 0 1];
 
@@ -484,7 +484,7 @@ Determine whether the polyhedral complex is pure.
 
 # Examples
 ```jldoctest
-julia> IM = IncidenceMatrix([[1,2,3],[1,3,4]]);
+julia> IM = incidence_matrix([[1,2,3],[1,3,4]]);
 
 julia> VR = [0 0; 1 0; 1 1; 0 1];
 
@@ -504,7 +504,7 @@ Compute the dimension of the polyhedral complex.
 
 # Examples
 ```jldoctest
-julia> IM = IncidenceMatrix([[1,2,3],[1,3,4]]);
+julia> IM = incidence_matrix([[1,2,3],[1,3,4]]);
 
 julia> VR = [0 0; 1 0; 1 1; 0 1];
 
@@ -524,7 +524,7 @@ Return the polyhedra of a given dimension in the polyhedral complex `PC`.
 
 # Examples
 ```jldoctest
-julia> IM = IncidenceMatrix([[1,2,3],[1,3,4]]);
+julia> IM = incidence_matrix([[1,2,3],[1,3,4]]);
 
 julia> VR = [0 0; 1 0; 1 1; 0 1];
 
@@ -596,7 +596,7 @@ Return the lineality space of `PC`.
 ```jldoctest
 julia> VR = [0 0 0; 1 0 0; 0 1 0; -1 0 0];
 
-julia> IM = IncidenceMatrix([[1,2,3],[1,3,4]]);
+julia> IM = incidence_matrix([[1,2,3],[1,3,4]]);
 
 julia> far_vertices = [2,3,4];
 
@@ -636,7 +636,7 @@ Return the lineality dimension of `PC`.
 ```jldoctest
 julia> VR = [0 0 0; 1 0 0; 0 1 0; -1 0 0];
 
-julia> IM = IncidenceMatrix([[1,2,3],[1,3,4]]);
+julia> IM = incidence_matrix([[1,2,3],[1,3,4]]);
 
 julia> far_vertices = [2,3,4];
 
@@ -661,7 +661,7 @@ faces of $PC$ of dimension $i$.
 ```jldoctest
 julia> VR = [0 0; 1 0; -1 0; 0 1];
 
-julia> IM = IncidenceMatrix([[1,2,4],[1,3,4]]);
+julia> IM = incidence_matrix([[1,2,4],[1,3,4]]);
 
 julia> far_vertices = [2,3,4];
 
@@ -689,7 +689,7 @@ Return the number of rays of `PC`.
 ```jldoctest
 julia> VR = [0 0; 1 0; -1 0; 0 1];
 
-julia> IM = IncidenceMatrix([[1,2,4],[1,3,4]]);
+julia> IM = incidence_matrix([[1,2,4],[1,3,4]]);
 
 julia> far_vertices = [2,3,4];
 
@@ -711,7 +711,7 @@ Return the number of vertices of `PC`.
 ```jldoctest
 julia> VR = [0 0; 1 0; -1 0; 0 1];
 
-julia> IM = IncidenceMatrix([[1,2,4],[1,3,4]]);
+julia> IM = incidence_matrix([[1,2,4],[1,3,4]]);
 
 julia> far_vertices = [2,3,4];
 
@@ -733,7 +733,7 @@ Return the total number of polyhedra in the polyhedral complex `PC`.
 ```jldoctest
 julia> VR = [0 0; 1 0; -1 0; 0 1];
 
-julia> IM = IncidenceMatrix([[1,2,4],[1,3,4]]);
+julia> IM = incidence_matrix([[1,2,4],[1,3,4]]);
 
 julia> far_vertices = [2,3,4];
 
@@ -754,7 +754,7 @@ Compute the codimension of a polyhedral complex.
 ```
 julia> VR = [0 0; 1 0; -1 0; 0 1];
 
-julia> IM = IncidenceMatrix([[1,2],[1,3],[1,4]]);
+julia> IM = incidence_matrix([[1,2],[1,3],[1,4]]);
 
 julia> far_vertices = [2,3,4];
 
@@ -777,7 +777,7 @@ subset of some $\mathbb{R}^n$.
 ```jldoctest
 julia> VR = [0 0; 1 0; -1 0; 0 1];
 
-julia> IM = IncidenceMatrix([[1,2],[1,3],[1,4]]);
+julia> IM = incidence_matrix([[1,2],[1,3],[1,4]]);
 
 julia> PC = polyhedral_complex(IM, VR)
 Polyhedral complex in ambient dimension 2
diff --git a/src/PolyhedralGeometry/PolyhedralComplex/standard_constructions.jl b/src/PolyhedralGeometry/PolyhedralComplex/standard_constructions.jl
index 1eb7cf43e588..711c3a6d87c4 100644
--- a/src/PolyhedralGeometry/PolyhedralComplex/standard_constructions.jl
+++ b/src/PolyhedralGeometry/PolyhedralComplex/standard_constructions.jl
@@ -5,7 +5,7 @@ Return the common refinement of two polyhedral complexes.
 
 # Examples
 ```jldoctest
-julia> IM = IncidenceMatrix([[1,2,3]])
+julia> IM = incidence_matrix([[1,2,3]])
 1×3 IncidenceMatrix
 [1, 2, 3]
 
@@ -51,7 +51,7 @@ Return the k-skeleton of a polyhedral complex.
 
 # Examples
 ```jldoctest
-julia> IM = IncidenceMatrix([[1,2,3]])
+julia> IM = incidence_matrix([[1,2,3]])
 1×3 IncidenceMatrix
 [1, 2, 3]
 
diff --git a/src/PolyhedralGeometry/PolyhedralFan/constructors.jl b/src/PolyhedralGeometry/PolyhedralFan/constructors.jl
index f54219634937..9005ace2a38d 100644
--- a/src/PolyhedralGeometry/PolyhedralFan/constructors.jl
+++ b/src/PolyhedralGeometry/PolyhedralFan/constructors.jl
@@ -43,7 +43,7 @@ To obtain the upper half-space of the plane:
 ```jldoctest
 julia> R = [1 0; 1 1; 0 1; -1 0; 0 -1];
 
-julia> IM=IncidenceMatrix([[1,2],[2,3],[3,4],[4,5],[1,5]]);
+julia> IM = incidence_matrix([[1,2],[2,3],[3,4],[4,5],[1,5]]);
 
 julia> PF=polyhedral_fan(IM, R)
 Polyhedral fan in ambient dimension 2
@@ -55,7 +55,7 @@ julia> R = [1 0 0; 0 0 1];
 
 julia> L = [0 1 0];
 
-julia> IM = IncidenceMatrix([[1],[2]]);
+julia> IM = incidence_matrix([[1],[2]]);
 
 julia> PF=polyhedral_fan(IM, R, L)
 Polyhedral fan in ambient dimension 3
diff --git a/src/PolyhedralGeometry/PolyhedralFan/properties.jl b/src/PolyhedralGeometry/PolyhedralFan/properties.jl
index 13d9a3fbeec3..2b157f65c563 100644
--- a/src/PolyhedralGeometry/PolyhedralFan/properties.jl
+++ b/src/PolyhedralGeometry/PolyhedralFan/properties.jl
@@ -50,7 +50,7 @@ julia> matrix(QQ, rays(NF))
 ```
 The following fan has no rays:
 ```
-julia> IM = IncidenceMatrix([[1,2],[2,3]]);
+julia> IM = incidence_matrix([[1,2],[2,3]]);
 
 julia> R = [1 0 0; 0 1 0; -1 0 0];
 
@@ -287,7 +287,7 @@ Return the dimension of `PF`.
 This fan in the plane contains a 2-dimensional cone and is thus 2-dimensional
 itself.
 ```jldoctest
-julia> PF = polyhedral_fan(IncidenceMatrix([[1, 2], [3]]), [1 0; 0 1; -1 -1]);
+julia> PF = polyhedral_fan(incidence_matrix([[1, 2], [3]]), [1 0; 0 1; -1 -1]);
 
 julia> dim(PF)
 2
@@ -304,7 +304,7 @@ Return the number of maximal cones of `PF`.
 The cones given in this construction are non-redundant. Thus there are two
 maximal cones.
 ```jldoctest
-julia> PF = polyhedral_fan(IncidenceMatrix([[1, 2], [3]]), [1 0; 0 1; -1 -1]);
+julia> PF = polyhedral_fan(incidence_matrix([[1, 2], [3]]), [1 0; 0 1; -1 -1]);
 
 julia> n_maximal_cones(PF)
 2
@@ -321,7 +321,7 @@ Return the number of cones of `PF`.
 The cones given in this construction are non-redundant. There are six
 cones in this fan.
 ```jldoctest
-julia> PF = polyhedral_fan(IncidenceMatrix([[1, 2], [3]]), [1 0; 0 1; -1 -1])
+julia> PF = polyhedral_fan(incidence_matrix([[1, 2], [3]]), [1 0; 0 1; -1 -1])
 Polyhedral fan in ambient dimension 2
 
 julia> n_cones(PF)
@@ -441,7 +441,7 @@ This fan consists of two cones, one containing all the points with $y ≤ 0$ and
 one containing all the points with $y ≥ 0$. The fan's lineality is the common
 lineality of these two cones, i.e. in $x$-direction.
 ```jldoctest
-julia> PF = polyhedral_fan(IncidenceMatrix([[1, 2, 3], [3, 4, 1]]), [1 0; 0 1; -1 0; 0 -1])
+julia> PF = polyhedral_fan(incidence_matrix([[1, 2, 3], [3, 4, 1]]), [1 0; 0 1; -1 0; 0 -1])
 Polyhedral fan in ambient dimension 2
 
 julia> lineality_space(PF)
@@ -500,7 +500,7 @@ Determine whether `PF` is smooth.
 Even though the cones of this fan cover the positive orthant together, one of
 these und thus the whole fan is not smooth.
 ```jldoctest
-julia> PF = polyhedral_fan(IncidenceMatrix([[1, 2], [2, 3]]), [0 1; 2 1; 1 0]);
+julia> PF = polyhedral_fan(incidence_matrix([[1, 2], [2, 3]]), [0 1; 2 1; 1 0]);
 
 julia> is_smooth(PF)
 false
@@ -516,7 +516,7 @@ Determine whether `PF` is regular, i.e. the normal fan of a polytope.
 # Examples
 This fan is not complete and thus not regular.
 ```jldoctest
-julia> PF = polyhedral_fan(IncidenceMatrix([[1, 2], [3]]), [1 0; 0 1; -1 -1]);
+julia> PF = polyhedral_fan(incidence_matrix([[1, 2], [3]]), [1 0; 0 1; -1 -1]);
 
 julia> is_regular(PF)
 false
@@ -531,7 +531,7 @@ Determine whether `PF` is pure, i.e. all maximal cones have the same dimension.
 
 # Examples
 ```jldoctest
-julia> PF = polyhedral_fan(IncidenceMatrix([[1, 2], [3]]), [1 0; 0 1; -1 -1]);
+julia> PF = polyhedral_fan(incidence_matrix([[1, 2], [3]]), [1 0; 0 1; -1 -1]);
 
 julia> is_pure(PF)
 false
@@ -547,7 +547,7 @@ dimension.
 
 # Examples
 ```jldoctest
-julia> PF = polyhedral_fan(IncidenceMatrix([[1, 2], [3]]), [1 0; 0 1; -1 -1]);
+julia> PF = polyhedral_fan(incidence_matrix([[1, 2], [3]]), [1 0; 0 1; -1 -1]);
 
 julia> is_fulldimensional(PF)
 true
diff --git a/src/PolyhedralGeometry/SubdivisionOfPoints/constructors.jl b/src/PolyhedralGeometry/SubdivisionOfPoints/constructors.jl
index 7411b26fee52..5373f5dd966c 100644
--- a/src/PolyhedralGeometry/SubdivisionOfPoints/constructors.jl
+++ b/src/PolyhedralGeometry/SubdivisionOfPoints/constructors.jl
@@ -52,7 +52,7 @@ triangulation.
 ```jldoctest
 julia> moaepts = [4 0 0; 0 4 0; 0 0 4; 2 1 1; 1 2 1; 1 1 2];
 
-julia> moaeimnonreg0 = IncidenceMatrix([[4,5,6],[1,4,2],[2,4,5],[2,3,5],[3,5,6],[1,3,6],[1,4,6]]);
+julia> moaeimnonreg0 = incidence_matrix([[4,5,6],[1,4,2],[2,4,5],[2,3,5],[3,5,6],[1,3,6],[1,4,6]]);
 
 julia> MOAE = subdivision_of_points(moaepts, moaeimnonreg0)
 Subdivision of points in ambient dimension 3
diff --git a/src/PolyhedralGeometry/SubdivisionOfPoints/functions.jl b/src/PolyhedralGeometry/SubdivisionOfPoints/functions.jl
index 50775eeb099f..72b5b251074f 100644
--- a/src/PolyhedralGeometry/SubdivisionOfPoints/functions.jl
+++ b/src/PolyhedralGeometry/SubdivisionOfPoints/functions.jl
@@ -16,7 +16,7 @@ weights, but it will not be full-dimensional.
 ```jldoctest
 julia> moaepts = [4 0 0; 0 4 0; 0 0 4; 2 1 1; 1 2 1; 1 1 2];
 
-julia> moaeimnonreg0 = IncidenceMatrix([[4,5,6],[1,4,2],[2,4,5],[2,3,5],[3,5,6],[1,3,6],[1,4,6]]);
+julia> moaeimnonreg0 = incidence_matrix([[4,5,6],[1,4,2],[2,4,5],[2,3,5],[3,5,6],[1,3,6],[1,4,6]]);
 
 julia> MOAE = subdivision_of_points(moaepts, moaeimnonreg0)
 Subdivision of points in ambient dimension 3
diff --git a/src/PolyhedralGeometry/SubdivisionOfPoints/properties.jl b/src/PolyhedralGeometry/SubdivisionOfPoints/properties.jl
index 761ce6ac8811..f7e03198f29c 100644
--- a/src/PolyhedralGeometry/SubdivisionOfPoints/properties.jl
+++ b/src/PolyhedralGeometry/SubdivisionOfPoints/properties.jl
@@ -14,7 +14,7 @@ Display the points of the "mother of all examples" non-regular triangulation.
 ```jldoctest
 julia> moaepts = [4 0 0; 0 4 0; 0 0 4; 2 1 1; 1 2 1; 1 1 2];
 
-julia> moaeimnonreg0 = IncidenceMatrix([[4,5,6],[1,4,2],[2,4,5],[2,3,5],[3,5,6],[1,3,6],[1,4,6]]);
+julia> moaeimnonreg0 = incidence_matrix([[4,5,6],[1,4,2],[2,4,5],[2,3,5],[3,5,6],[1,3,6],[1,4,6]]);
 
 julia> MOAE = subdivision_of_points(moaepts, moaeimnonreg0);
 
@@ -62,7 +62,7 @@ julia> moaepts = [4 0 0; 0 4 0; 0 0 4; 2 1 1; 1 2 1; 1 1 2]
  1  2  1
  1  1  2
 
-julia> moaeimnonreg0 = IncidenceMatrix([[4,5,6],[1,4,2],[2,4,5],[2,3,5],[3,5,6],[1,3,6],[1,4,6]])
+julia> moaeimnonreg0 = incidence_matrix([[4,5,6],[1,4,2],[2,4,5],[2,3,5],[3,5,6],[1,3,6],[1,4,6]])
 7×6 IncidenceMatrix
 [4, 5, 6]
 [1, 2, 4]
@@ -245,7 +245,7 @@ triangulation of six points.
 ```jldoctest
 julia> moaepts = [4 0 0; 0 4 0; 0 0 4; 2 1 1; 1 2 1; 1 1 2];
 
-julia> moaeimnonreg0 = IncidenceMatrix([[4,5,6],[1,4,2],[2,4,5],[2,3,5],[3,5,6],[1,3,6],[1,4,6]]);
+julia> moaeimnonreg0 = incidence_matrix([[4,5,6],[1,4,2],[2,4,5],[2,3,5],[3,5,6],[1,3,6],[1,4,6]]);
 
 julia> MOAE = subdivision_of_points(moaepts, moaeimnonreg0);
 
diff --git a/src/PolyhedralGeometry/helpers.jl b/src/PolyhedralGeometry/helpers.jl
index 8c4c1ff3ee38..dedeebea8ba9 100644
--- a/src/PolyhedralGeometry/helpers.jl
+++ b/src/PolyhedralGeometry/helpers.jl
@@ -6,9 +6,9 @@ import Polymake: IncidenceMatrix
 A matrix with boolean entries. Each row corresponds to a fixed element of a collection of mathematical objects and the same holds for the columns and a second (possibly equal) collection. A `1` at entry `(i, j)` is interpreted as an incidence between object `i` of the first collection and object `j` of the second one.
 
 # Examples
-Note that the input and print of an `IncidenceMatrix` lists the non-zero indices for each row.
+Note that the input of this example and the print of an `IncidenceMatrix` list the non-zero indices for each row.
 ```jldoctest
-julia> IM = IncidenceMatrix([[1,2,3],[4,5,6]])
+julia> IM = incidence_matrix([[1,2,3],[4,5,6]])
 2×6 IncidenceMatrix
 [1, 2, 3]
 [4, 5, 6]
@@ -27,6 +27,109 @@ julia> IM[:, 4]
 """
 IncidenceMatrix
 
+@doc raw"""
+    incidence_matrix(r::Base.Integer, c::Base.Integer)
+
+Return an `IncidenceMatrix` of size r x c whose entries are all `false`.
+
+# Examples
+```jldoctest
+julia> IM = incidence_matrix(8, 5)
+8×5 IncidenceMatrix
+[]
+[]
+[]
+[]
+[]
+[]
+[]
+[]
+
+```
+"""
+incidence_matrix(r::Base.Integer, c::Base.Integer) = IncidenceMatrix(undef, r, c)
+
+@doc raw"""
+    incidence_matrix(mat::Union{AbstractMatrix{Bool}, IncidenceMatrix})
+
+Convert `mat` to an `IncidenceMatrix`.
+
+# Examples
+```jldoctest
+julia> IM = incidence_matrix([true false true false true false; false true false true false true])
+2×6 IncidenceMatrix
+[1, 3, 5]
+[2, 4, 6]
+
+```
+"""
+incidence_matrix(mat::Union{AbstractMatrix{Bool},IncidenceMatrix}) = IncidenceMatrix(mat)
+
+@doc raw"""
+    incidence_matrix(mat::AbstractMatrix)
+
+Convert the `0`/`1` matrix `mat` to an `IncidenceMatrix`. Entries become `true` if the initial entry is `1` and `false` if the initial entry is `0`.
+
+# Examples
+```jldoctest
+julia> IM = incidence_matrix([1 0 1 0 1 0; 0 1 0 1 0 1])
+2×6 IncidenceMatrix
+[1, 3, 5]
+[2, 4, 6]
+
+```
+"""
+function incidence_matrix(mat::AbstractMatrix)
+  m, n = size(mat)
+  for i in 1:m
+    for j in 1:n
+      iszero(mat[i, j]) || isone(mat[i, j]) ||
+        throw(
+          ArgumentError("incidence_matrix requires matrices with 0/1 or boolean entries.")
+        )
+    end
+  end
+  return IncidenceMatrix(mat)
+end
+
+@doc raw"""
+    incidence_matrix(r::Base.Integer, c::Base.Integer, incidenceRows::AbstractVector{<:AbstractVector{<:Base.Integer}})
+
+Return an `IncidenceMatrix` of size r x c. The i-th element of `incidenceRows` lists the indices of the `true` entries of the i-th row.
+
+# Examples
+```jldoctest
+julia> IM = incidence_matrix(3, 4, [[2, 3], [1]])
+3×4 IncidenceMatrix
+[2, 3]
+[1]
+[]
+
+```
+"""
+incidence_matrix(
+  r::Base.Integer,
+  c::Base.Integer,
+  incidenceRows::AbstractVector{<:AbstractVector{<:Base.Integer}},
+) = IncidenceMatrix(r, c, incidenceRows)
+
+@doc raw"""
+    incidence_matrix(incidenceRows::AbstractVector{<:AbstractVector{<:Base.Integer}})
+
+Return an `IncidenceMatrix` where the i-th element of `incidenceRows` lists the indices of the `true` entries of the i-th row. The dimensions of the result are the smallest possible row and column count that can be deduced from the input.
+
+# Examples
+```jldoctest
+julia> IM = incidence_matrix([[2, 3], [1]])
+2×3 IncidenceMatrix
+[2, 3]
+[1]
+
+```
+"""
+incidence_matrix(incidenceRows::AbstractVector{<:AbstractVector{<:Base.Integer}}) =
+  IncidenceMatrix(incidenceRows)
+
 number_of_rows(i::IncidenceMatrix) = Polymake.nrows(i)
 number_of_columns(i::IncidenceMatrix) = Polymake.ncols(i)
 
@@ -40,7 +143,7 @@ Return the indices where the `n`-th row of `i` is `true`, as a `Set{Int}`.
 
 # Examples
 ```jldoctest
-julia> IM = IncidenceMatrix([[1,2,3],[4,5,6]])
+julia> IM = incidence_matrix([[1,2,3],[4,5,6]])
 2×6 IncidenceMatrix
 [1, 2, 3]
 [4, 5, 6]
@@ -62,7 +165,7 @@ Return the indices where the `n`-th column of `i` is `true`, as a `Set{Int}`.
 
 # Examples
 ```jldoctest
-julia> IM = IncidenceMatrix([[1,2,3],[4,5,6]])
+julia> IM = incidence_matrix([[1,2,3],[4,5,6]])
 2×6 IncidenceMatrix
 [1, 2, 3]
 [4, 5, 6]
diff --git a/src/PolyhedralGeometry/iterators.jl b/src/PolyhedralGeometry/iterators.jl
index f87b07a8d9d6..d7b56e325448 100644
--- a/src/PolyhedralGeometry/iterators.jl
+++ b/src/PolyhedralGeometry/iterators.jl
@@ -324,6 +324,8 @@ function IncidenceMatrix(iter::SubObjectIterator)
   end
 end
 
+incidence_matrix(iter::SubObjectIterator) = IncidenceMatrix(iter)
+
 # primitive generators only for ray based iterators
 matrix(R::ZZRing, iter::SubObjectIterator{RayVector{QQFieldElem}}) =
   matrix(R, Polymake.common.primitive(matrix_for_polymake(iter)))
diff --git a/src/PolyhedralGeometry/triangulations.jl b/src/PolyhedralGeometry/triangulations.jl
index 774f29a27e00..084cf77e19af 100644
--- a/src/PolyhedralGeometry/triangulations.jl
+++ b/src/PolyhedralGeometry/triangulations.jl
@@ -471,7 +471,7 @@ Compute a triangulation of the square
 ```jldoctest
 julia> C = cube(2);
 
-julia> cells = IncidenceMatrix([[1,2,3],[2,3,4]]);
+julia> cells = incidence_matrix([[1,2,3],[2,3,4]]);
 
 julia> S = subdivision_of_points(C, cells)
 Subdivision of points in ambient dimension 2
diff --git a/src/TropicalGeometry/curve.jl b/src/TropicalGeometry/curve.jl
index e6f70218c390..32ef9762d5eb 100644
--- a/src/TropicalGeometry/curve.jl
+++ b/src/TropicalGeometry/curve.jl
@@ -59,7 +59,7 @@ Return the embedded tropical curve consisting of the polyhedral complex `Sigma`
 ```jldoctest
 julia> verticesAndRays = [0 0; 1 0; 0 1; -1 -1];
 
-julia> incidenceMatrix = IncidenceMatrix([[1,2],[1,3],[1,4]]);
+julia> incidenceMatrix = incidence_matrix([[1,2],[1,3],[1,4]]);
 
 julia> rayIndices = [2,3,4];
 
@@ -222,7 +222,7 @@ end
 
 # # Examples
 # ```jldoctest
-# julia> IM = IncidenceMatrix([[1,2],[1,3],[1,4],[2,3],[2,4],[3,4]]);
+# julia> IM = incidence_matrix([[1,2],[1,3],[1,4],[2,3],[2,4],[3,4]]);
 
 # julia> tc = TropicalCurve(IM)
 # ERROR: MethodError: no method matching TropicalCurve(::Polymake.LibPolymake.IncidenceMatrixAllocated{Polymake.LibPolymake.NonSymmetric})
@@ -291,7 +291,7 @@ end
 
 # # Examples
 # ```jldoctest
-# julia> IM = IncidenceMatrix([[1,2],[1,3],[1,4],[2,3],[2,4],[3,4]]);
+# julia> IM = incidence_matrix([[1,2],[1,3],[1,4],[2,3],[2,4],[3,4]]);
 
 # julia> tc = TropicalCurve(IM)
 # ERROR: MethodError: no method matching TropicalCurve(::Polymake.LibPolymake.IncidenceMatrixAllocated{Polymake.LibPolymake.NonSymmetric})
@@ -352,7 +352,7 @@ end
 
 # # Examples
 # ```jldoctest
-# julia> IM = IncidenceMatrix([[1,2],[1,3],[1,4],[2,3],[2,4],[3,4]]);
+# julia> IM = incidence_matrix([[1,2],[1,3],[1,4],[2,3],[2,4],[3,4]]);
 
 # julia> tc = TropicalCurve(IM)
 # ERROR: MethodError: no method matching TropicalCurve(::Polymake.LibPolymake.IncidenceMatrixAllocated{Polymake.LibPolymake.NonSymmetric})
@@ -411,7 +411,7 @@ end
 # ```jldoctest
 # julia> cg = complete_graph(5);
 
-# julia> IM1=IncidenceMatrix([[src(e), dst(e)] for e in edges(cg)])
+# julia> IM1 = incidence_matrix([[src(e), dst(e)] for e in edges(cg)])
 # 10×5 IncidenceMatrix
 # [1, 2]
 # [1, 3]
@@ -438,7 +438,7 @@ end
 
 # julia> cg2 = complete_graph(3);
 
-# julia> IM2=IncidenceMatrix([[src(e), dst(e)] for e in edges(cg2)])
+# julia> IM2 = incidence_matrix([[src(e), dst(e)] for e in edges(cg2)])
 # 3×3 IncidenceMatrix
 # [1, 2]
 # [1, 3]
@@ -456,7 +456,7 @@ end
 #  [1] top-level scope
 #    @ none:1
 
-# julia> IM3 = IncidenceMatrix([[1,2],[2,3],[3,4],[4,5],[1,5]])
+# julia> IM3 = incidence_matrix([[1,2],[2,3],[3,4],[4,5],[1,5]])
 # 5×5 IncidenceMatrix
 # [1, 2]
 # [2, 3]
@@ -500,7 +500,7 @@ end
 # # ```julia
 # # julia> using Revise, Plots, Oscar, Test;
 
-# # julia> IM = IncidenceMatrix([[1,2],[1,3],[1,4]]);
+# # julia> IM = incidence_matrix([[1,2],[1,3],[1,4]]);
 
 # # julia> VR = [0 0; 1 0; -1 0; 0 1];
 
diff --git a/src/TropicalGeometry/variety.jl b/src/TropicalGeometry/variety.jl
index 7c2822564274..66f481378f70 100644
--- a/src/TropicalGeometry/variety.jl
+++ b/src/TropicalGeometry/variety.jl
@@ -46,7 +46,7 @@ Return the `TropicalVariety` whose polyhedral complex is `Sigma` with multiplici
 
 # Examples
 ```jldoctest
-julia> Sigma = polyhedral_complex(IncidenceMatrix([[1],[2]]), [[0],[1]])
+julia> Sigma = polyhedral_complex(incidence_matrix([[1],[2]]), [[0],[1]])
 Polyhedral complex in ambient dimension 1
 
 julia> tropical_variety(Sigma)
diff --git a/test/AlgebraicGeometry/ToricVarieties/affine_normal_varieties.jl b/test/AlgebraicGeometry/ToricVarieties/affine_normal_varieties.jl
index b412a3037c33..a21649946a33 100644
--- a/test/AlgebraicGeometry/ToricVarieties/affine_normal_varieties.jl
+++ b/test/AlgebraicGeometry/ToricVarieties/affine_normal_varieties.jl
@@ -5,7 +5,7 @@
   antv3 = affine_normal_toric_variety(antv2)
   antv4 = affine_normal_toric_variety(Oscar.positive_hull([1 0]))
   antv5 = affine_space(NormalToricVariety, 2)
-  antv6 = normal_toric_variety(IncidenceMatrix([[1,2,3,4]]), [[1, 0, 0], [1, 0, 1], [1, 1, 1], [1, 1, 0]])
+  antv6 = normal_toric_variety(incidence_matrix([[1,2,3,4]]), [[1, 0, 0], [1, 0, 1], [1, 1, 1], [1, 1, 0]])
   
   @testset "Basic properties" begin
     @test is_smooth(antv) == false
diff --git a/test/AlgebraicGeometry/ToricVarieties/direct_products.jl b/test/AlgebraicGeometry/ToricVarieties/direct_products.jl
index a617af339b27..05ced475ffe1 100644
--- a/test/AlgebraicGeometry/ToricVarieties/direct_products.jl
+++ b/test/AlgebraicGeometry/ToricVarieties/direct_products.jl
@@ -1,6 +1,6 @@
 @testset "Direct products" begin
     
-  F5 = normal_toric_variety(IncidenceMatrix([[1, 2], [2, 3], [3, 4], [4, 1]]), [[1, 0], [0, 1], [-1, 5], [0, -1]])
+  F5 = normal_toric_variety(incidence_matrix([[1, 2], [2, 3], [3, 4], [4, 1]]), [[1, 0], [0, 1], [-1, 5], [0, -1]])
     P2 = normal_toric_variety(normal_fan(Oscar.simplex(2)))
     variety = F5 * P2
     
diff --git a/test/AlgebraicGeometry/ToricVarieties/intersection_numbers.jl b/test/AlgebraicGeometry/ToricVarieties/intersection_numbers.jl
index e83d5a469cf9..384dfd4c369a 100644
--- a/test/AlgebraicGeometry/ToricVarieties/intersection_numbers.jl
+++ b/test/AlgebraicGeometry/ToricVarieties/intersection_numbers.jl
@@ -2,9 +2,9 @@
   
   antv = affine_normal_toric_variety(Oscar.positive_hull([1 1; -1 1]))
   
-  antv2 = normal_toric_variety(IncidenceMatrix([[1,2,3,4]]), [[1, 0, 0], [1, 0, 1], [1, 1, 1], [1, 1, 0]])
+  antv2 = normal_toric_variety(incidence_matrix([[1,2,3,4]]), [[1, 0, 0], [1, 0, 1], [1, 1, 1], [1, 1, 0]])
   
-  v = normal_toric_variety(IncidenceMatrix([[1], [2], [3]]), [[1, 0], [0, 1], [-1, -1]])
+  v = normal_toric_variety(incidence_matrix([[1], [2], [3]]), [[1, 0], [0, 1], [-1, -1]])
     
   dP1 = del_pezzo_surface(NormalToricVariety, 1)
   c0 = cohomology_class(dP1, gens(cohomology_ring(dP1))[1])
diff --git a/test/AlgebraicGeometry/ToricVarieties/line_bundle_cohomologies.jl b/test/AlgebraicGeometry/ToricVarieties/line_bundle_cohomologies.jl
index ec7bc5eed471..ba7074feeb1c 100644
--- a/test/AlgebraicGeometry/ToricVarieties/line_bundle_cohomologies.jl
+++ b/test/AlgebraicGeometry/ToricVarieties/line_bundle_cohomologies.jl
@@ -5,7 +5,7 @@
   F5 = hirzebruch_surface(NormalToricVariety, 5)
   
   ray_generators = [[1, 0, 0,-2,-3], [0, 1, 0,-2,-3], [0, 0, 1,-2,-3], [-1,-1,-1,-2,-3], [0, 0, 0, 1, 0], [0, 0, 0, 0, 1], [0, 0, 0,-2,-3]]
-  max_cones = IncidenceMatrix([[1, 2, 3, 5, 6], [1, 2, 3, 5, 7], [1, 2, 3, 6, 7], [2, 3, 4, 5, 6], [2, 3, 4, 5, 7], [2, 3, 4, 6, 7], [1, 3, 4, 5, 6], [1, 3, 4, 5, 7], [1, 3, 4, 6, 7], [1, 2, 4, 5, 6], [1, 2, 4, 5, 7], [1, 2, 4, 6, 7]])
+  max_cones = incidence_matrix([[1, 2, 3, 5, 6], [1, 2, 3, 5, 7], [1, 2, 3, 6, 7], [2, 3, 4, 5, 6], [2, 3, 4, 5, 7], [2, 3, 4, 6, 7], [1, 3, 4, 5, 6], [1, 3, 4, 5, 7], [1, 3, 4, 6, 7], [1, 2, 4, 5, 6], [1, 2, 4, 5, 7], [1, 2, 4, 6, 7]])
   weierstrass_over_p3 = normal_toric_variety(max_cones, ray_generators; non_redundant = true)
   
   l = toric_line_bundle(dP3, [1, 2, 3, 4])
diff --git a/test/AlgebraicGeometry/ToricVarieties/subvarieties.jl b/test/AlgebraicGeometry/ToricVarieties/subvarieties.jl
index d54d38dc197d..f50c9a3305a1 100644
--- a/test/AlgebraicGeometry/ToricVarieties/subvarieties.jl
+++ b/test/AlgebraicGeometry/ToricVarieties/subvarieties.jl
@@ -1,6 +1,6 @@
 @testset "Closed subvarieties" begin
   
-  antv = normal_toric_variety(IncidenceMatrix([[1,2,3,4]]), [[1, 0, 0], [1, 0, 1], [1, 1, 1], [1, 1, 0]])
+  antv = normal_toric_variety(incidence_matrix([[1,2,3,4]]), [[1, 0, 0], [1, 0, 1], [1, 1, 1], [1, 1, 0]])
   
   ntv = normal_toric_variety(cube(2))
   (x1, x2, y1, y2) = gens(cox_ring(ntv));
diff --git a/test/AlgebraicGeometry/ToricVarieties/toric_divisor_classes.jl b/test/AlgebraicGeometry/ToricVarieties/toric_divisor_classes.jl
index 775e514cd954..724830921061 100644
--- a/test/AlgebraicGeometry/ToricVarieties/toric_divisor_classes.jl
+++ b/test/AlgebraicGeometry/ToricVarieties/toric_divisor_classes.jl
@@ -5,7 +5,7 @@
   P2 = projective_space(NormalToricVariety, 2)
   
   rs = [[0, 1, 0], [-1, 1, 1], [-1, 0, 0], [-1, -2, -2], [0, -1, -1], [0, 0, 1], [-1, -1, -1], [1, 0, 0]]
-  cs = IncidenceMatrix([[1, 2, 3], [1, 4, 5], [3, 6, 7], [4, 5, 6], [1, 4, 7], [5, 6, 8], [1, 6, 8], [2, 3, 6], [1, 5, 8], [1, 3, 7], [1, 2, 6], [4, 6, 7]])
+  cs = incidence_matrix([[1, 2, 3], [1, 4, 5], [3, 6, 7], [4, 5, 6], [1, 4, 7], [5, 6, 8], [1, 6, 8], [2, 3, 6], [1, 5, 8], [1, 3, 7], [1, 2, 6], [4, 6, 7]])
   test_space = normal_toric_variety(cs, rs)
 
   DC = toric_divisor_class(F5, [ZZRingElem(0), ZZRingElem(0)])
diff --git a/test/AlgebraicGeometry/ToricVarieties/toric_schemes.jl b/test/AlgebraicGeometry/ToricVarieties/toric_schemes.jl
index b723ea89f198..8eb6b42cbd9f 100644
--- a/test/AlgebraicGeometry/ToricVarieties/toric_schemes.jl
+++ b/test/AlgebraicGeometry/ToricVarieties/toric_schemes.jl
@@ -97,7 +97,7 @@
 end
 
 @testset "Lazy gluings" begin
-  f = polyhedral_fan(IncidenceMatrix([[1, 2, 3],[1, 4, 5, 6]]), [0 0 1; 1 0 1; 0 1 0; -1 0 1; -1 -1 1; 0 -1 1])
+  f = polyhedral_fan(incidence_matrix([[1, 2, 3],[1, 4, 5, 6]]), [0 0 1; 1 0 1; 0 1 0; -1 0 1; -1 -1 1; 0 -1 1])
   number_of_maximal_cones(f)
   ntv = normal_toric_variety(f)
   X = Oscar.underlying_scheme(ntv)
diff --git a/test/Combinatorics/Graph.jl b/test/Combinatorics/Graph.jl
index 0e6ef72156b8..d742171d7b67 100644
--- a/test/Combinatorics/Graph.jl
+++ b/test/Combinatorics/Graph.jl
@@ -67,7 +67,7 @@
         @test n_vertices(egplc) == 2
         @test n_edges(egplc) == 1
 
-        @test incidence_matrix(egtriangle) == IncidenceMatrix([[1,2],[1,3],[2,3]])
+        @test incidence_matrix(egtriangle) == incidence_matrix([[1,2],[1,3],[2,3]])
 
         @test is_isomorphic(dual_graph(convex_hull([0 0 0; 1 0 0], nothing, [0 1 0])), Graph{Undirected}(2))
         @test is_isomorphic(dual_graph(convex_hull([0 0 0], [0 0 1; 0 1 0; 1 0 0])), complete_graph(3))
diff --git a/test/Combinatorics/SimplicialComplexes.jl b/test/Combinatorics/SimplicialComplexes.jl
index c89c3baa1856..3ba32a31e8b6 100644
--- a/test/Combinatorics/SimplicialComplexes.jl
+++ b/test/Combinatorics/SimplicialComplexes.jl
@@ -4,7 +4,7 @@
 
     sphere = simplicial_complex([[1, 2, 3], [1, 2, 4], [1, 3, 4], [2, 3, 4]])
 
-    sphere2 = simplicial_complex(IncidenceMatrix([[1, 2, 3], [1, 2, 4], [1, 3, 4], [2, 3, 4]]))
+    sphere2 = simplicial_complex(incidence_matrix([[1, 2, 3], [1, 2, 4], [1, 3, 4], [2, 3, 4]]))
 
     not_known_to_be_a_ball = SimplicialComplex(Polymake.topaz.SimplicialComplex(FACETS=[[0,1]], BALL=nothing))
 
diff --git a/test/PolyhedralGeometry/cone.jl b/test/PolyhedralGeometry/cone.jl
index ba97410c6449..8adcfc125b94 100644
--- a/test/PolyhedralGeometry/cone.jl
+++ b/test/PolyhedralGeometry/cone.jl
@@ -57,7 +57,7 @@
       @test linear_inequality_matrix(facets(S, Cone1)) == _oscar_matrix_from_property(f, facets(S, Cone1))
       @test Oscar.linear_matrix_for_polymake(facets(S, Cone1)) == _polymake_matrix_from_property(facets(S, Cone1))
       @test _check_im_perm_rows(ray_indices(facets(S, Cone1)), [[1], [2]])
-      @test _check_im_perm_rows(IncidenceMatrix(facets(S, Cone1)), [[1], [2]])
+      @test _check_im_perm_rows(incidence_matrix(facets(S, Cone1)), [[1], [2]])
     end
     @test _check_im_perm_rows(facets(IncidenceMatrix, Cone1), [[1], [2]])
     @test facets(Halfspace, Cone1) isa SubObjectIterator{LinearHalfspace{T}}
@@ -96,13 +96,13 @@
     @test length(faces(Cone4, 1)) == 2
     @test issetequal(faces(Cone2, 2), positive_hull.(Ref(f), [[1 0 0], [0 0 1]], [[0 1 0]]))
     @test _check_im_perm_rows(ray_indices(faces(Cone2, 2)), [[1], [2]])
-    @test _check_im_perm_rows(IncidenceMatrix(faces(Cone2, 2)), [[1], [2]])
+    @test _check_im_perm_rows(incidence_matrix(faces(Cone2, 2)), [[1], [2]])
     @test _check_im_perm_rows(faces(IncidenceMatrix, Cone2, 2), [[1], [2]])
     @test issetequal(faces(Cone4, 1), positive_hull.(Ref(f), [[0 0 1], [1 0 0]]))
     @test _check_im_perm_rows(ray_indices(faces(Cone4, 1)), [[1], [2]])
-    @test _check_im_perm_rows(IncidenceMatrix(faces(Cone4, 1)), [[1], [2]])
+    @test _check_im_perm_rows(incidence_matrix(faces(Cone4, 1)), [[1], [2]])
     @test _check_im_perm_rows(faces(IncidenceMatrix, Cone4, 1), [[1], [2]])
-    @test _check_im_perm_rows(IncidenceMatrix(faces(Cone5, 1)), [[1], [2], [3], [4]])
+    @test _check_im_perm_rows(incidence_matrix(faces(Cone5, 1)), [[1], [2], [3], [4]])
     @test isnothing(faces(Cone2, 1))
 
     @test f_vector(Cone5) == [4, 4]
diff --git a/test/PolyhedralGeometry/extended.jl b/test/PolyhedralGeometry/extended.jl
index 10de0db2606e..c7fc6ee058a6 100644
--- a/test/PolyhedralGeometry/extended.jl
+++ b/test/PolyhedralGeometry/extended.jl
@@ -195,8 +195,8 @@
     @test_throws ArgumentError convex_hull(vertices(Pos_poly), collect(vertices(Pos_poly)))
     @test_throws ArgumentError positive_hull(collect(vertices(Pos_poly)))
 
-    @test_throws ArgumentError IncidenceMatrix(lineality_space(Pos_poly))
-    IM = IncidenceMatrix([[1]])
+    @test_throws ArgumentError incidence_matrix(lineality_space(Pos_poly))
+    IM = incidence_matrix([[1]])
     lincone = positive_hull([1 0 0], [0 1 0])
 
     @test positive_hull(rays_modulo_lineality(lincone)...) == lincone
diff --git a/test/PolyhedralGeometry/group.jl b/test/PolyhedralGeometry/group.jl
index a0945be71ea9..e749e1c87b1b 100644
--- a/test/PolyhedralGeometry/group.jl
+++ b/test/PolyhedralGeometry/group.jl
@@ -45,7 +45,7 @@
     M = matroid_from_nonbases([[1, 2, 3, 4], [1, 2, 5, 6], [1, 3, 5, 7], [3, 4, 6, 8]], 8)
     GM = automorphism_group(M)
     @test is_trivial(GM)
-    IM = IncidenceMatrix(bases(M))
+    IM = incidence_matrix(bases(M))
     GIM = automorphism_group(IM)
     @test length(GIM) == 2
     @test haskey(GIM, :on_cols)
diff --git a/test/PolyhedralGeometry/lineality.jl b/test/PolyhedralGeometry/lineality.jl
index 4b76e87430ec..abdf1280b44c 100644
--- a/test/PolyhedralGeometry/lineality.jl
+++ b/test/PolyhedralGeometry/lineality.jl
@@ -82,7 +82,7 @@
 
   @testset "PolyhedralComplex" begin
     VR = [0 0 0; 1 0 0; 0 1 0; -1 0 0]
-    IM = IncidenceMatrix([[1, 2, 3], [1, 3, 4]])
+    IM = incidence_matrix([[1, 2, 3], [1, 3, 4]])
     far_vertices = [2, 3, 4]
     L = [0 0 1]
     PC = polyhedral_complex(IM, VR, far_vertices, L)
diff --git a/test/PolyhedralGeometry/polyhedral_complex.jl b/test/PolyhedralGeometry/polyhedral_complex.jl
index 45aaa618e7ad..605ca05a1da6 100644
--- a/test/PolyhedralGeometry/polyhedral_complex.jl
+++ b/test/PolyhedralGeometry/polyhedral_complex.jl
@@ -1,5 +1,5 @@
 @testset "PolyhedralComplex{$T}" for (f, T) in _prepare_scalar_types()
-  I = IncidenceMatrix([[1, 2, 3], [2, 4]])
+  I = incidence_matrix([[1, 2, 3], [2, 4]])
   P = f.([0 0; 1 0; 0 1; 1 1])
   P2 = f.([0 0 0; 1 0 0; 0 1 0; 1 1 0])
   F = [4]
@@ -86,7 +86,7 @@
     @test vertex_indices(maximal_polyhedra(PCFLN)) == I[:, 1:3]
     @test ray_indices(maximal_polyhedra(PCFLN)) == I[:, 4:4]
     @test vertex_and_ray_indices(maximal_polyhedra(PCFLN)) == I
-    @test IncidenceMatrix(maximal_polyhedra(PCFLN)) == I
+    @test incidence_matrix(maximal_polyhedra(PCFLN)) == I
     @test maximal_polyhedra(IncidenceMatrix, PCFLN) == I
 
     @test polyhedral_complex(maximal_polyhedra(PCF)) isa PolyhedralComplex
@@ -109,7 +109,7 @@
   @testset "Fan conversion" begin
     F1 = normal_fan(cube(f, 2))
     F2 = normal_fan(convex_hull(f, [0 0; 1 0]))
-    IM = IncidenceMatrix([[1, 2], [2, 3], [4]])
+    IM = incidence_matrix([[1, 2], [2, 3], [4]])
     R = [0 1 0; 0 0 1; 0 -1 0; 0 -1 -1]
     F3 = polyhedral_fan(f, IM, R, [1 0 0])
     for F in [F1, F2, F3]
diff --git a/test/PolyhedralGeometry/polyhedral_fan.jl b/test/PolyhedralGeometry/polyhedral_fan.jl
index e012d5615d7b..133ba92b6a59 100644
--- a/test/PolyhedralGeometry/polyhedral_fan.jl
+++ b/test/PolyhedralGeometry/polyhedral_fan.jl
@@ -10,8 +10,8 @@
   @test polyhedral_fan([Cone4, Cone5]) isa PolyhedralFan{T}
   F0 = polyhedral_fan([Cone4, Cone5])
   I3 = [1 0 0; 0 1 0; 0 0 1]
-  incidence1 = IncidenceMatrix([[1, 2], [2, 3]])
-  incidence2 = IncidenceMatrix([[1, 2]])
+  incidence1 = incidence_matrix([[1, 2], [2, 3]])
+  incidence2 = incidence_matrix([[1, 2]])
   @test polyhedral_fan(f, incidence1, I3) isa PolyhedralFan{T}
   F1 = polyhedral_fan(f, incidence1, I3)
   F1NR = polyhedral_fan(f, incidence1, I3; non_redundant=true)
@@ -46,7 +46,7 @@
     @test maximal_cones(F1) isa SubObjectIterator{Cone{T}}
     @test dim.(maximal_cones(F1)) == [2, 2]
     @test _check_im_perm_rows(ray_indices(maximal_cones(F1)), incidence1)
-    @test _check_im_perm_rows(IncidenceMatrix(maximal_cones(F1)), incidence1)
+    @test _check_im_perm_rows(incidence_matrix(maximal_cones(F1)), incidence1)
     @test _check_im_perm_rows(maximal_cones(IncidenceMatrix, F1), incidence1)
     @test number_of_maximal_cones(F1) == 2
     @test lineality_space(F2) isa SubObjectIterator{RayVector{T}}
@@ -60,7 +60,7 @@
     @test rays.(cones(F2, 2)) == [[], []]
     @test isnothing(cones(F2, 1))
     @test _check_im_perm_rows(ray_indices(cones(F1, 2)), incidence1)
-    @test _check_im_perm_rows(IncidenceMatrix(cones(F1, 2)), incidence1)
+    @test _check_im_perm_rows(incidence_matrix(cones(F1, 2)), incidence1)
     @test _check_im_perm_rows(cones(IncidenceMatrix, F1, 2), incidence1)
 
     II = ray_indices(maximal_cones(NFsquare))
@@ -74,7 +74,7 @@
     @test f_vector(NFsquare) == [4, 4]
     @test rays(F1NR) == collect(eachrow(I3))
     @test _check_im_perm_rows(ray_indices(maximal_cones(F1NR)), incidence1)
-    @test _check_im_perm_rows(IncidenceMatrix(maximal_cones(F1NR)), incidence1)
+    @test _check_im_perm_rows(incidence_matrix(maximal_cones(F1NR)), incidence1)
     @test n_rays(F2NR) == 0
     @test lineality_dim(F2NR) == 1
     RMLF2NR = rays_modulo_lineality(F2NR)
@@ -82,7 +82,7 @@
     @test RMLF2NR[:rays_modulo_lineality] == collect(eachrow(R))
     @test lineality_space(F2NR) == collect(eachrow(L))
     @test _check_im_perm_rows(ray_indices(maximal_cones(F2NR)), incidence2)
-    @test _check_im_perm_rows(IncidenceMatrix(maximal_cones(F2NR)), incidence2)
+    @test _check_im_perm_rows(incidence_matrix(maximal_cones(F2NR)), incidence2)
 
     C = positive_hull(f, identity_matrix(ZZ, 0))
     pf = polyhedral_fan(C)
@@ -138,7 +138,7 @@ end
 end
 
 @testset "Star Subdivision" begin
-  f = polyhedral_fan(IncidenceMatrix([[1, 2, 3, 4]]), [1 0 0; 1 1 0; 1 1 1; 1 0 1])
+  f = polyhedral_fan(incidence_matrix([[1, 2, 3, 4]]), [1 0 0; 1 1 0; 1 1 1; 1 0 1])
   @test is_pure(f)
   @test is_fulldimensional(f)
   v0 = [1; 0; 0]
@@ -151,7 +151,7 @@ end
   @test number_of_maximal_cones(sf1) == 3
   @test number_of_maximal_cones(sf2) == 4
 
-  ff = polyhedral_fan(IncidenceMatrix([[1], [2, 3]]), [1 0 0; -1 0 0; 0 1 0])
+  ff = polyhedral_fan(incidence_matrix([[1], [2, 3]]), [1 0 0; -1 0 0; 0 1 0])
   @test !is_pure(ff)
   @test !is_fulldimensional(ff)
   w0 = [1; 0; 0]
diff --git a/test/PolyhedralGeometry/polyhedron.jl b/test/PolyhedralGeometry/polyhedron.jl
index 405e7a31cbfb..3556449f31dd 100644
--- a/test/PolyhedralGeometry/polyhedron.jl
+++ b/test/PolyhedralGeometry/polyhedron.jl
@@ -118,16 +118,16 @@
                      convex_hull.(Ref(f), [[-1 -1; -1 1], [1 -1; 1 1], [-1 -1; 1 -1], [-1 1; 1 1]]))
     @test _check_im_perm_rows(vertex_indices(faces(square, 1)),
                               [[1, 3], [2, 4], [1, 2], [3, 4]])
-    @test ray_indices(faces(square, 1)) == IncidenceMatrix(4, 0)
+    @test ray_indices(faces(square, 1)) == incidence_matrix(4, 0)
     @test _check_im_perm_rows(vertex_and_ray_indices(faces(square, 1)),
                               [[2, 4], [1, 3], [1, 2], [3, 4]])
-    @test _check_im_perm_rows(IncidenceMatrix(faces(square, 1)),
+    @test _check_im_perm_rows(incidence_matrix(faces(square, 1)),
                               [[1, 3], [2, 4], [1, 2], [3, 4]])
     @test _check_im_perm_rows(faces(IncidenceMatrix, square, 1),
                               [[1, 3], [2, 4], [1, 2], [3, 4]])
     @test _check_im_perm_rows(facet_indices(vertices(square)),
                               [[1, 3], [2, 3], [1, 4], [2, 4]])
-    @test _check_im_perm_rows(IncidenceMatrix(vertices(square)),
+    @test _check_im_perm_rows(incidence_matrix(vertices(square)),
                               [[1, 3], [2, 3], [1, 4], [2, 4]])
     @test _check_im_perm_rows(vertices(IncidenceMatrix, square),
                               [[1, 3], [2, 3], [1, 4], [2, 4]])
@@ -137,7 +137,7 @@
     @test _check_im_perm_rows(vertex_indices(faces(Pos, 1)), [[1], [1], [1]])
     @test _check_im_perm_rows(ray_indices(faces(Pos, 1)), [[1], [2], [3]])
     @test _check_im_perm_rows(vertex_and_ray_indices(faces(Pos, 1)), [[1, 4], [2, 4], [3, 4]])
-    @test _check_im_perm_rows(IncidenceMatrix(faces(Pos, 1)), [[1, 4], [2, 4], [3, 4]])
+    @test _check_im_perm_rows(incidence_matrix(faces(Pos, 1)), [[1, 4], [2, 4], [3, 4]])
     @test _check_im_perm_rows(faces(IncidenceMatrix, Pos, 1), [[1, 4], [2, 4], [3, 4]])
     @test isnothing(faces(Q2, 0))
     v = vertices(minkowski_sum(Q0, square))
@@ -166,17 +166,17 @@
       end
       @test _check_im_perm_rows(ray_indices(facets(S, Pos)), [[2, 3], [1, 3], [1, 2]])
       @test _check_im_perm_rows(vertex_and_ray_indices(facets(S, Pos)), [[2, 3, 4], [1, 3, 4], [1, 2, 4]])
-      @test _check_im_perm_rows(IncidenceMatrix(facets(S, Pos)), [[2, 3, 4], [1, 3, 4], [1, 2, 4]])
+      @test _check_im_perm_rows(incidence_matrix(facets(S, Pos)), [[2, 3, 4], [1, 3, 4], [1, 2, 4]])
       @test _check_im_perm_rows(vertex_indices(facets(S, Pos)), [[1], [1], [1]])
     end
     @test _check_im_perm_rows(facets(IncidenceMatrix, Pos),
         [[2, 3, 4], [1, 3, 4], [1, 2, 4]]
     )
     @test _check_im_perm_rows(facet_indices(vertices(Pos)), [[1, 2, 3]])
-    @test _check_im_perm_rows(IncidenceMatrix(vertices(Pos)), [[1, 2, 3]])
+    @test _check_im_perm_rows(incidence_matrix(vertices(Pos)), [[1, 2, 3]])
     @test _check_im_perm_rows(vertices(IncidenceMatrix, Pos), [[1, 2, 3]])
     @test _check_im_perm_rows(facet_indices(rays(Pos)), [[1, 3], [2, 3], [1, 2]])
-    @test _check_im_perm_rows(IncidenceMatrix(rays(Pos)), [[2, 3], [1, 3], [1, 2]])
+    @test _check_im_perm_rows(incidence_matrix(rays(Pos)), [[2, 3], [1, 3], [1, 2]])
     @test _check_im_perm_rows(rays(IncidenceMatrix, Pos), [[1, 3], [2, 3], [1, 2]])
     @test facets(Pair, Pos) isa SubObjectIterator{Pair{Matrix{T},T}}
     @test facets(Pos) isa SubObjectIterator{AffineHalfspace{T}}
@@ -606,7 +606,7 @@
           @test affine_inequality_matrix(facets(S, D)) == matrix(R, hcat(-b, vcat(A...)))
           @test halfspace_matrix_pair(facets(S, D)) == (A=matrix(R, vcat(A...)), b=b)
 
-          @test ray_indices(facets(S, D)) == IncidenceMatrix(12, 0)
+          @test ray_indices(facets(S, D)) == incidence_matrix(12, 0)
           @test _check_im_perm_rows(vertex_indices(facets(S, D)), [
             [1, 3, 5, 9, 10],
             [1, 2, 3, 4, 6],
@@ -635,7 +635,7 @@
             [9, 10, 14, 17, 19],
             [15, 17, 18, 19, 20],
           ])
-          @test _check_im_perm_rows(IncidenceMatrix(facets(S, D)), [
+          @test _check_im_perm_rows(incidence_matrix(facets(S, D)), [
             [1, 3, 5, 9, 10],
             [1, 2, 3, 4, 6],
             [1, 2, 5, 7, 8],
@@ -666,9 +666,9 @@
         [9, 10, 14, 17, 19],
         [15, 17, 18, 19, 20],
       ])
-      @test facet_indices(rays(D)) == IncidenceMatrix(0, 12)
-      @test IncidenceMatrix(rays(D)) == IncidenceMatrix(0, 12)
-      @test rays(IncidenceMatrix, D) == IncidenceMatrix(0, 12)
+      @test facet_indices(rays(D)) == incidence_matrix(0, 12)
+      @test incidence_matrix(rays(D)) == incidence_matrix(0, 12)
+      @test rays(IncidenceMatrix, D) == incidence_matrix(0, 12)
       @test _check_im_perm_rows(facet_indices(vertices(D)), [
         [1, 2, 3],
         [2, 3, 6],
@@ -691,7 +691,7 @@
         [9, 11, 12],
         [4, 9, 12],
       ])
-      @test _check_im_perm_rows(IncidenceMatrix(vertices(D)), [
+      @test _check_im_perm_rows(incidence_matrix(vertices(D)), [
         [1, 2, 3],
         [2, 3, 6],
         [1, 2, 7],
diff --git a/test/PolyhedralGeometry/scalar_types.jl b/test/PolyhedralGeometry/scalar_types.jl
index 802b1794a306..450e0ccc6280 100644
--- a/test/PolyhedralGeometry/scalar_types.jl
+++ b/test/PolyhedralGeometry/scalar_types.jl
@@ -47,7 +47,7 @@
   @test length(lattice_points(sd)) == 11
 
   let pc = polyhedral_complex(
-      E, IncidenceMatrix(facets(sd)), vertices(sd); non_redundant=true
+      E, incidence_matrix(facets(sd)), vertices(sd); non_redundant=true
     )
     @test issetequal(maximal_polyhedra(pc), faces(sd, 2))
   end
diff --git a/test/PolyhedralGeometry/setup_tests.jl b/test/PolyhedralGeometry/setup_tests.jl
index 6d2ba9757c99..c6cc8613922f 100644
--- a/test/PolyhedralGeometry/setup_tests.jl
+++ b/test/PolyhedralGeometry/setup_tests.jl
@@ -9,7 +9,7 @@ if !isdefined(Main, :_prepare_scalar_types)
 end
 
 function _check_im_perm_rows(inc::IncidenceMatrix, o)
-  oinc = IncidenceMatrix(o)
+  oinc = incidence_matrix(o)
   nr, nc = size(inc)
   (nr, nc) == size(oinc) &&
     issetequal(Polymake.row.(Ref(inc), 1:nr),
diff --git a/test/PolyhedralGeometry/subdivision_of_points.jl b/test/PolyhedralGeometry/subdivision_of_points.jl
index 2608f0e46f6f..eebf7eae44b8 100644
--- a/test/PolyhedralGeometry/subdivision_of_points.jl
+++ b/test/PolyhedralGeometry/subdivision_of_points.jl
@@ -2,7 +2,7 @@
   C = cube(2)
   square_weights = [0, 0, 1, 2]
   square_max_cells = [[1, 2, 3], [2, 3, 4]]
-  square_incidence = IncidenceMatrix(square_max_cells)
+  square_incidence = incidence_matrix(square_max_cells)
 
   square_by_weights = subdivision_of_points(C, square_weights)
   square_by_cells = subdivision_of_points(C, square_max_cells)
@@ -21,9 +21,9 @@
 
   moaepts = [4 0 0; 0 4 0; 0 0 4; 2 1 1; 1 2 1; 1 1 2]
   fulldim_moaepts = moaepts[:, 2:3]
-  moaeimreg0 = IncidenceMatrix(1, 6)
+  moaeimreg0 = incidence_matrix(1, 6)
   moaeimreg0[1, :] = [1, 1, 1, 0, 0, 0]
-  moaeimnonreg0 = IncidenceMatrix([
+  moaeimnonreg0 = incidence_matrix([
     [4, 5, 6], [1, 4, 2], [2, 4, 5], [2, 3, 5], [3, 5, 6], [1, 3, 6], [1, 4, 6]
   ])
 
diff --git a/test/PolyhedralGeometry/types.jl b/test/PolyhedralGeometry/types.jl
index f54e31a349b7..d079c17df422 100644
--- a/test/PolyhedralGeometry/types.jl
+++ b/test/PolyhedralGeometry/types.jl
@@ -1,12 +1,21 @@
 @testset "types" begin
   @testset "IncidenceMatrix" begin
-    im = IncidenceMatrix([[1, 2, 3], [4, 5, 6]])
+    im = incidence_matrix([[1, 2, 3], [4, 5, 6]])
     @test nrows(im) == 2
     @test ncols(im) == 6
     @test row(im, 1) isa Set{Int}
     @test row(im, 1) == Set{Int}([1, 2, 3])
     @test column(im, 2) isa Set{Int}
     @test column(im, 2) == Set{Int}([1])
+
+    @testset "IncidenceMatrix constructions" begin
+      @test incidence_matrix([true true true false false false; false false false true true true]) == im
+      @test incidence_matrix([1 1 1 0 0 0; 0 0 0 1 1 1]) == im
+      @test incidence_matrix(4, 2) == incidence_matrix([0 0; 0 0; 0 0; 0 0])
+      @test incidence_matrix(im) == im
+      @test incidence_matrix(3, 8, [[1, 2, 3], [4, 5, 6]]) == incidence_matrix([1 1 1 0 0 0 0 0; 0 0 0 1 1 1 0 0; 0 0 0 0 0 0 0 0])
+      @test_throws ArgumentError incidence_matrix([1 1 1 0 0 0; 0 0 0 1 2 1])
+    end
   end
 
   a = [1, 2, 3]
diff --git a/test/Serialization/PolyhedralGeometry.jl b/test/Serialization/PolyhedralGeometry.jl
index 0d8b8e1f4c03..9d04515719d3 100644
--- a/test/Serialization/PolyhedralGeometry.jl
+++ b/test/Serialization/PolyhedralGeometry.jl
@@ -87,7 +87,7 @@ using Oscar: _integer_variables
     end
 
     @testset "PolyhedralComplex" begin
-      IM = IncidenceMatrix([[1,2,3],[1,3,4]])
+      IM = incidence_matrix([[1,2,3],[1,3,4]])
       vr = [0 0; 1 0; 1 1; 0 1]
       PC = polyhedral_complex(IM, vr)
       test_save_load_roundtrip(path, PC) do loaded
@@ -156,7 +156,7 @@ using Oscar: _integer_variables
 
     @testset "SubdivisionOfPoints" begin
       moaepts = [4 0 0; 0 4 0; 0 0 4; 2 1 1; 1 2 1; 1 1 2]
-      moaeimnonreg0 = IncidenceMatrix([[4,5,6],[1,4,2],[2,4,5],[2,3,5],[3,5,6],[1,3,6],[1,4,6]])
+      moaeimnonreg0 = incidence_matrix([[4,5,6],[1,4,2],[2,4,5],[2,3,5],[3,5,6],[1,3,6],[1,4,6]])
       MOAE = subdivision_of_points(moaepts, moaeimnonreg0)
       test_save_load_roundtrip(path, MOAE) do loaded
         @test number_of_maximal_cells(MOAE) == number_of_maximal_cells(loaded)
diff --git a/test/Serialization/TropicalGeometry.jl b/test/Serialization/TropicalGeometry.jl
index 17f01e907466..7d80d1c430cf 100644
--- a/test/Serialization/TropicalGeometry.jl
+++ b/test/Serialization/TropicalGeometry.jl
@@ -10,7 +10,7 @@
             end
 
             @testset "Embedded" begin
-                IM = IncidenceMatrix([[1,2],[1,3],[1,4]])
+                IM = incidence_matrix([[1,2],[1,3],[1,4]])
                 VR = [0 0; 1 0; -1 0; 0 1]
                 PC = polyhedral_complex(QQFieldElem, IM, VR)
                 TC = tropical_curve(PC)