**Abstract** : Isomanifolds are the generalization of isosurfaces to arbitrary dimension and codimension, i.e. submanifolds of R d defined as the zero set of some multivariate multivalued smooth function f : R^d → R^{d−m} where m is the intrinsic dimension of the manifold. A natural way to approximate a smooth isomanifold M is to consider its Piecewise-Linear (PL) approximation or meshM based on a triangulation T of the ambient space R d. In this paper, we present a simple algorithm to construct such an approximation for arbitrary m and d up to a given precision D. The complexity of our algorithm is polynomial in d and δ, and exponential in m. Since, by previous results,M is O(D 2)-close and isotopic to M when δ = Ω(d 2.5), our algorithm constructs a geometrically close and topologically correct PL-approximation of isomanifolds of low dimensions in polynomial time. The algorithm is practical and can handle cases that are far ahead of the state-of-the-art. Combining this algorithm with dimensionality reduction techniques, the dependency on d in the sizes of the output sample and mesh can be completely removed with high probability. The crux of our algorithm is to use for the ambient triangulation T a regular triangulation from a particular family. This family consists of Freudenthal-Kuhn triangulations and their images through affine mappings. It also includes Coxeter triangulations of typeÃ d. We introduce an elegant and very compact data structure to implicitly store the full facial structure of such triangulations. This data structure allows to retrieve the faces or the cofaces of a simplex of any dimension in an output sensitive way, which is essential for our application and is of independent interest.