Directional is a C++ geometry processing library focused on directional-field processing. The functionality is based on the definitions and taxonomy surveyed theoretically in 1 and 2, and in many newer papers in the literature, cited within context in this tutorial and in the code. Directional contains tools to edit, analyze, and visualize directional fields of various degrees, orders, and symmetries.
Directional has gone through a major paradigm shift in version 3.0: the discretization is abstracted from face-based fields into general discrete tangent bundles that can represent a rich class of directional-field representations. As a result, many of the algorithm of the library can now work seamlessly on several representations without change (as an example, power fields can be computed on either vertex-based or face-based fields with the same function). The library comprises two basic elements:
-
Classes representing the basic tangent bundle and field structures. They are built with inheritance so that functions can be written polymorphically and consequently algorithms can be applied to several representations that have the minimal set of required properties.
-
Standalone functions, implementing directional-field algorithms, that take these classes as parameters.
Our paradigm avoids buffed classes with a complicated nested hierarchy; instead, the member functions in the classes are minimal,
and only used to implement the essential properties of a geometric object (for instance, the connection between tangent spaces).
The structure of these standalone functions follows the general philosophy of libigl: the library is header only, where each header contains a set of functions closely related (for instance, the precomputation and computation of some directional quantity over a mesh). For the most part, one header contains only one function. The atomic data structures are, for the most part, simple matrices in Eigen,
The visualization is done through a specialized class DirectionalViewer
, on the basis of libigl viewer, with many extended options that allow the rendering of directional fields.
The header files contain documentation of the parameters to each function and their required composition; in this tutorial we will mostly tie the functionality of Directional to the theoretical concepts of directional fields and the methods to process and visualize them.
This tutorial comprises an exhaustive set of examples that demonstrates the capabilities of Directional, where every subchapter entails a single concept. The tutorial code can be installed by going into the tutorial
folder from the main Directional folder, and typing the following instructions in a terminal:
mkdir build
cd build
cmake -DCMAKE_BUILD_TYPE=Release ../
make
This will build all tutorial chapters in the build
folder. The necessary dependencies will be appended and built automatically. To build in windows, use the cmake-gui ..
options instead of the last two commands, and create the project using Visual Studio, with the proper tutorial subchapter as the "startup project".
To access a single example, say 202_Sampling
, go to the build
subfolder, and the executable will be there. Command-line arguments are never required; the data is read from the shared
folder directly for each example.
Most examples require a modest amount of user interaction; the instructions of what to do are given in the command-line output upon execution.
In the continuum, directional fields are represented as elements of tangent spaces, where a tangent space is a linear space attached to each point of a manifold. The set of all such tangent spaces is called a tangent bundle.
In the discrete setting, a tangent bundle is a finite set of
The central data structure of Directional is the discrete tangent bundle, which is a finite and discrete graph TangentBundle
. Discrete tangent bundles can supply one of more of the following properties:
-
Intrinsic parameterization: vectors in a single tangent space are represented with coordinates
$\left(x_1,\cdots,x_d\right)$ in an arbitrary, but fixed basis of dimension$d$ . The basis itself does not need to be specified for many algorithms, it can stay abstract. -
Adjacency: a tangent bundle is represented as a graph where the tangent spaces are nodes, and the graph edges are adjacency relations, that encode the local environment of any tangent space. This in fact encodes the combinatorial and topological structure of the underlying discrete manifold, akin to the smooth structure in the continuum.
-
Connection: a connection defines the notion of parallelity in two adjacent tangent spaces, which encodes a metric structure on the underlying manifold, and allows to compute derivatives and curvature. Specifically, a connection is a directed-edge-based quantity on the tangent-bundle graph, given as a change of coordinates between the bases of the source tangent space to the target tangent space. This can be represented as an orthogonal matrix.
-
metric the metric on the bundle is supplied in two quantities:
connectionMass
is the weight of each edge$E_{TB}$ , andtangentSpaceMass
packs the weights of an inner product on vectors on the bundle. That is, it has either$V_{TB}$ or$V_{TB}+E_{TB}$ that are the non-zero components of the vector mass matrix. -
Cycles:
$G_{TB}$ can be equipped with a notion of cycles$F_{TB}$ that are simply-connected closed chains of spaces. Given a connection, the holonomy of a cycle is the failure of a vector to return to itself following a set of parallel transports around the cycle. Concretely, it is encoded in the logarithm the product of connection matrices. In$d=2$ , holonomy, which is then a single angle, is equivalent to the curvature of the cycle modulo$2\pi$ . There are two types of cycles, which define the topology of the underlying manifold: local cycles (aking to "faces" in$G_{TB}$ ), which are small closed loops, and global cycles, which can be homological generators and boundary loops. The singularities of fields are defined on the local cycles. -
Cochain complex: vector fields are often objects of differential geometry where the underlying manifold is equipped with notions of gradient, curl, divergence, and other vector calculus identities. A cochain complex is defined by a gradient operator
$G$ and a curl operator$C$ where$C \cdot G=0$ . The "cochain" aspect can be purely algebraic and not rely on any explicit nodes in$G_{TB}$ ; we call$G\cdot f$ for any "scalar" function$f$ exact (or conservative) fields and fields where$C\cdot v=0$ closed (or curl-free) fields. A cochain complex is enough to define deRham cohomology with the correct dimensions$ker(C)/im(G)$ . However to extract the explicit harmonic fields we need the next property. The combination of a metric and a cochain complex allows for a well-defined notion of Helmholtz-Hodge decomposition: any vector field$v$ can be decomposed into exact part$Gf$ , coexact part$M^{-1}C^Tg$ and harmonic part$h$ as follows [~boksebeld_2022]:
The inner product also introduces the discrete divergence operator
Oftentimes the above intrinsic quantities are enough for all algorithms; nevertheless, even for reasons of input, output, and visualization, a TangentBundle
will contain the following, embedding-based extrinsic quantities:
-
Sources and normals: point locations and their normals (the codimensional directions of the manifold), which define the tangent planes of the manifold in the Euclidean space. Note that this doesn't mean it's a watertight mesh; that could also encode a pointcloud, for instance.
-
extrinsic to intrinsic and back: functionality that projects extrinsic vectors into intrinsic space (might lose information), or produces the extrinsic representation of an intrinsic directional.
-
Cycles sources and normals. Like (7), but for the cycles themselves, This would be where singularities are
located
in space, and might just be a visualization quantity.
Two main types are currently implemented in directional: IntrinsicFaceTangentBundle
implements face-based tangent spaces, where the fields are tangent to the natural plane supporting triangles of surface meshes, and IntrinsicVertexTangentBundle
implements vertex-based intrinsic tangent spaces, which parameterize the cone environment of the
For example, this is how IntrinsicFaceTangentBundle
implements the above quantities:
- Intrinsic parameterization: a local basis in every face.
- Adjacency: dual (inner) edges between any two faces.
- Connection: the rotation matrix between the bases of any two adjacent faces.
- Cycles: the local cycles are around vertex
$1$ -rings, where singularities are then defined as (dual) vertex values, and global cycles are dual loops of generators and boundaries. - Cochain complex: the classical FEM face-based gradient and curl.
- Inner product: the natural face-based mass matrix (just a diagonal matrix of face areas).
- Sources are face barycenters, and normals are just face normals.
- The projection to the supporting plane of the face and encoding in local coordinates.
- Vertices and vertex normals (area-weighted from adjacent faces).
Some of the choices above can of course be variated to different flavors of face-based fields (for instance, the metric culminating in the mass weights). IntrinsicFaceTangentBundle
wraps around a an orientable input triangle mesh in a single connected-component. There are no general limitations on the genus or the boundaries. If your input comprises several connected components altogether, you should use several tangent bundles.
The representation of a directional field is its encoding in each discrete tangent plane. The most important element is the number of vectors in each tangent plane, which we denote as the degree of the field CartesianField
. Directional currently supports the following variants of Cartesian fields1:
-
Raw - a vector of
$d\times N$ entries represents an intrinsic$1^N$ -vector (a directional with$N$ independent vectors in each tangent plane) a dimension-dominant ordering:$(X_{1,1},\ldots, X_{1,d}),(X_{1,2},\ldots,X_{2,d}),\ldots (X_{N,1},\ldots, X_{N,d})$ per face. Vectors are assumed to be ordered in counterclockwise order in most Directional functions that process raw fields. the memory complexity is then$dN|V_{TB}|$ for the entire directional field. A Cartesian Field indicates being a raw-field by settingCartesianField::fieldType
todirectional::RAW_FIELD
. -
Power Field - This is a unique type to
$d=2$ . It encodes an$N$ -rotational-symmetric ($N$ -RoSy) object as a single complex number$y=u^N$ relative to the local basis in the tangent space, where the$N$ -RoSy is the set of roots$u \cdot e^{\frac{2\pi i k}{N}}, k \in [0,N-1]$ . The magnitude is also encoded this way, though it may be neglected in some applications. The memory complexity is then$2|V_{TB}|$ . -
PolyVector - Also unique to
$d=2$ , this is a generalization of power fields that represents an$N$ -directional object in a tangent space as the coefficients$a$ of a monic complex polynomial$f(z)=z^N+\sum_{i=0}^{N-1}{a_iz^i}$ , which roots$u$ are the encoded$1^N$ -vector field. In case where the field is an$N$ -RoSy, all coefficients but$a_0$ are zero. Note: A PolyVector that represents a perfect$N$ -RoSy would have all$a_i=0,\ \forall i>0$ , but$a_0$ would have opposite sign from the power-field representation of the same$N$ -RoSy. This is since the power field represents$u^N$ directly, whereas a PolyVector represents the coefficients of$z^N-U^N$ in this case. The memory complexity is$2N|V_{TB}|$ .
Directional provides a number of conversion functions to switch between different representations. Each of the functions is of the form rep1_to_rep2
, where rep1
and rep2
are the representation names in the above list. e.g., polyvector_to_raw()
. Some possible combinations are given by composing two functions in sequence. However, note that not every conversion is possible; for instance, it is not possible to convert from PolyVectors to power fields, as they do not possess the same power of expression. Converting into the more explicit raw representation is often needed for I/O and visualization.
Directional uses a specialized class called DirectionalViewer
which inherits and extends the libigl Viewer
class, augmenting it with functionality that pertains to directional fields. A mesh is then stored with its accompanying geometric quantities: the field, edge, vertex, and face-based scalar data, isolines, and more, that we detail in the tutorial per chapter in context. Like libigl viewer, Directional supports independent multiple meshes, each with its own set of quantities. Internally, the visualization schemes carete sub-meshes which serve as layers on the original meshes: arrows for glyphs, bars for edge highlights, etc. In practice this is encapsulated from the user and does not need to be controlled directly.
The most basic operation on directional fields is reading them from a file and drawing them in the most explicit way. In this example, a mesh and a field are read from a file and visualized as follows:
directional::readOFF(TUTORIAL_SHARED_PATH "/bumpy.off",mesh);
directional::read_raw_field(TUTORIAL_SHARED_PATH "/bumpy.rawfield", mesh, N, field);
directional::read_singularities(TUTORIAL_SHARED_PATH "/bumpy.sings", field);
directional::DirectionalViewer viewer;
viewer.set_mesh(mesh);
viewer.set_field(field);
The field is read in raw format (see File Formats), which is detailed in the Introduction. The field is face-based, and the singularities are consequently vertex-based,
The singularities and glyphs (and most other properties) can be toggled by functions of the type DirectionalViewer::toggle_field()
and DirectionalViewer::toggle_singularities()
.
Glyph Rendering on a mesh with singularities.
This example shows a Cartesian field computed (with the power field method described in Example 301 on either a vertex-based tangent bundle, or a face-based tangent bundle, to highlight the flexibility of choosing a discretization. The relevant code segments are:
directional::readOBJ(TUTORIAL_SHARED_PATH "/elephant.obj", mesh);
powerFaceField.init(mesh, POWER_FIELD, N);
powerVertexField.init(mesh, POWER_FIELD, N);
...
directional::power_field(mesh, constFaces, constVectors, Eigen::VectorXd::Constant(constFaces.size(),-1.0), N, powerFaceField);
directional::power_field(mesh, constVertices, constVectors, Eigen::VectorXd::Constant(constVertices.size(),-1.0), N, powerVertexField);
//computing power fields
directional::power_to_raw(powerFaceField, N, rawFaceField,true);
directional::power_to_raw(powerVertexField, N, rawVertexField,true);
One can see the stages of computing a field: first reading a mesh (readOBJ()
), then initializing the approxiate tangent bundle with the mesh (powerFace/VertexField.init()
), and then computing this power field on top of this representation (power_field()
). The field is converted to raw representation in power_to_raw()
) for later visualization.
Power fields on a face-based tangent bundle (left) and vertex-based (right)
This example demonstrates the editing paradigm in Directional, based on libigl picking. A face and a vector within the face are chosen, and clicking on a new direction for the vector changes it. Note the different colors for glyphs and selected faces. The specificiation of selection is done via the following code:
directionalViewer->set_selected_faces(selectedFaces);
directionalViewer->set_selected_vector(currF, currVec);
Editing several vectors on a single face.
Vector fields on surfaces are commonly visualized by tracing [streamlines] (https://en.wikipedia.org/wiki/Streamlines,_streaklines,_and_pathlines). Directional supports the seeding and tracing of streamlines, for all types of directionals. The seeds for the streamlines are initialized using DirectionalViewer::init_streamlines()
, and the lines are traced using DirectionalViewer::streamlines_next()
. Each call to DirectionalViewer::advance_streamlines()
extends each line by one triangle, allowing interactive rendering of the traced lines. The streamline have the same colors as the initial glyphs, where the colors fade into white as the streamline advance.
Interactive streamlines tracing.
It is possible to set and visualize scalar quantities on meshes at different discretization locations: either face based quantities that appear as flat colors per face, vertex-based (pointwise) quantities that interpolate linearly on faces, appearing smooth, and edge-based (integrated) quantities, that appear as flat quantities on a diamond mesh associates with each edge (taking a DirectionalViewer::set_X_data()
functions, that also allow the setting of the viewable range of the function (the rest is clipped).
Face-, Vertex- and edge-based data on a mesh, with a field as a layer of (white) glyphs.
On big meshes, it might appear cumbersome to view all glyphs on every face. It is possible to only view the glyphs on a subsample of faces, by using the sparsity
parameter in DirectionalViewer::set_field()
.
Dense and Sparse views of a field as glyphs.
In the following sections, we show some effects of working with different representations and converting between them.
One of the fundamental operations in directional-field processing is matching. That is, defining which vectors in tangent space
Given a raw field (in assumed CCW order in every tangent space), it is possible to devise the rotation angles
principal matching is done through the function principal_matching()
, that accepts a Cartesian field as a parameter, and computes the following:
- The matching on each (directed) TB-graph edge. This is stored in the
matching
member of the field class. - The indices of the cycles. The singular local cycles are are stored in the corresponding
singLocalCycles
andsingIndices
of the field class.
The singularities are computed as the index of each local cycle from the effort around it. The index of a cycle is the amount of rotations a directional object undergoes around a cycle. A directional must return to itself after a cycle, and therefore the index is an integer principal_matching()
does not update singularities around boundary or generator loops.
A Field is shown with singularities, and a single face is shown with the principal matching to its neighbors (in multiple colors).
This is an educational example that demonstrates the loss of information when generating a Cartesian field from rotation angles, and then trying to retrieve them back by principal matching. This causes low valence cycles and undersampling cause aliasing in the perceived field. There are three modes seen in the example:
-
In the polar mode, the user can prescribe the index of a singularity directly, and compute the field with index prescription (see example 401). With this, the rotation angles between adjacent faces can be arbitrarily large, and appear as noise in the low valence cycles.
-
In the principal-matching mode, the rotations are reconstructed from the field, without prior knowledge of the polar-prescribed rotations from the previous mode. The large rotation between adjacent faces is lost, which gives rise to a "singularity party": many perceived singularities or a lower index.
-
In the interpolation mode, the field is interpolated on the free faces (white) from the constrained faces (red), keeping the red band fixed from the polar mode. We see a field that is smooth in the Cartesian sense, with more uniformly-placed singularities.
Alternating: the polar mode, the principal-matching mode, and the Cartesian mode.
Given a matching (in this case, principal matching), it is possible to "comb" the field. That is, re-index each face (keeping the CCW order), so that the vector indexing aligns perfectly with the matching to the neighbors---then, the new matching on the dual edges becomes trivially zero. This operation is important in order to prepare a directional field for integration, for instance. In the presence of singularities, the field can only be combed up to a forest of paths that connect between singularities, also known as seams. Note that such paths do not necessarily cut the mesh into a simply-connected patch, but may only connects subgroups of singularities with indices adding up to an integer; as a trivial example, a 1-vector field is always trivially combed, even in the presence of integral singularities, and the set of seams is zero. The combing is done through the function directional::combing()
. The matching in the output combedField
is already set to the trivial matching in the combed regions, and the correct matching across the seam.
Colored indices of directionals, alternating between combed (with seams) and uncombed) indexing.
The Cartesian representation is a meta-category for representation of vectors in explicit coordinates, either
This representation is offered in 4, but they did not give it a specific name (the method in general is called "globally optimal"). We use the name "power fields" coined in 5.
A power field representation uses a complex basis in each tangent plane of a discrete tangent bundle, and represents an
By prescribing constraints
where directional::power_field()
. It is possible to softly prescribe the constraints
where
If the set
Hard (left) and soft (right) aligned constraints (yellow on red faces) interpolated to the rest of the mesh. Note the singularities that are discovered through principal matching.
A Polyvector field 3 is a generalization of power fields that allows to represent independent vectors in each tangent space, invariant to their order. The representation is as the coefficient set
where the roots
With the function polyvector_field()
one can solve the linear Polyvector problem in its full capacity; the input is a set of prescribed constraints
So that the set power_field()
is implemented by calling polyvector_field()
. This tutorial examples allows interacting with the alignment weights.
Top: Sharp-edge constraints (left; note sometimes more than one per face), Hard (middle) and soft (right) solution. Bottom: dominant-weighted smoothness (left), alignment (middle) and rotational symmetry (right).
This functionality only works with face-based fields via IntrinsicFaceTangentBundle
.
Vector-field guided surface parameterization is based on the idea of designing the candidate gradients
of the parameterization functions (which are tangent vector fields on the surface) instead of the functions themselves. Thus, vector-set fields (
A field that has zero PolyCurl everywhere is locally (away from singularities) integrable into curl_matching()
with that minimizes curl instead of rotation.
PolyCurl is iteratively reduced from an initial PolyVector field, Top: fields for iteration 0 (original), 10, and 50. Bottom: PolyCurl plots. The color is the norm of the vector of the roots of the PolyCurl. Its maximal value (infinity norm) appears below.
This functionality only works with face-based fields via IntrinsicFaceTangentBundle
.
Two tangent vectors
Conjugate vector fields are very important in architectural geometry: their integral lines form, informally speaking, an infinitesimal planar quad mesh. As such, the finite quad mesh that results from discretizing conjugate networks is a good candidate for consequent planarity parameterization 12.
Finding a conjugate vector field that satisfies given directional constraints is a standard problem in architectural geometry, which can be tackled by deforming a
A smooth $2^2$-PolyVector field (left) is deformed to become a conjugate field (right). Top: fields Bottom: conjugacy plots.
Polar fields are represented using angles. These angles may encode the rotation from some given basis on a tangent space (and so it is a "logarithmic" representation, when compared to Cartesian methods), or an angle difference between two neighboring tangent spaces (in the sense of deviation from parallel transport). The former usually requires integer variables for directional-field design. The latter does not, but state-of-the-art methods require the prescription of indices around independent dual cycles in the mesh. Currently, Directional supports the latter.
The notion of encoding rotation angles on dual edges, as means to encode deviation from parallel transport between adjacent tangent planes, appeared in several formats in the literature 13,8. The formulation and notation we use in Directional is that of Trivial Connections 13. Trivial connection solves for a single rotation angle
The basis cycles form the cycles around which curvatures (and singularities) have to be prescribed on the mesh. The sum on basis cycles is described in a sparse matrix directional::index_prescription()
, which can also accept a solver for precomputation, for the purpose of prefactoring
Indices are prescribed on several vertex singularities, and on a generator loop, to match the index theorem.
The full details of the method implemented in this chapter can be found in a technical report 14. Many of the therotical ideas for general IntrinsicFaceTangentBundle
.
Consider a seam edge
where
Seamless
In this example, we demonstrate the computation of such a integration, both permutationally, and fully seamless. The computed function is a
directional::IntegrationData intData(N);
directional::setup_integration(rawField, intData, meshCut, combedField);
intData.verbose=true;
intData.integralSeamless=false;
directional::integrate(combedField, intData, meshCut, cutUVRot ,cornerWholeUV);
cutUVRot=cutUVRot.block(0,0,cutUVRot.rows(),2);
intData.integralSeamless = true;
directional::integrate(combedField, intData, meshCut, cutUVFull,cornerWholeUV);
cutUVFull=cutUVFull.block(0,0,cutUVFull.rows(),2);
The data structure containing all relevant information about the integration is IntegrationData
. It contains some parameters that can be tuned to control the integration. Several relevant ones are:
double lengthRatio; //global scaling of functions
//Flags
bool integralSeamless; //Whether to do full translational seamless.
bool roundSeams; //Whether to round seams or round singularities
bool verbose; //output the integration log.
bool localInjectivity; //Enforce local injectivity; might result in failure!
lengthRatio
encodes a global scale for the Poisson problem (scaling the fields uniformly), where the ratio is measured against the bounding box diagonal. Some of the other parameters are demonstrated in the other examples in this chapter. The integrator takes the original (whole) mesh, and generates a cut-mesh (in VMeshCut,FMeshCut
) of disc-topology. The singularities are on the boundary of this mesh, and the function can consequently be defined without branching ambiguity on its vertices, with the appropriate permutation and translation across the cut seams.
Left: directional field. Center: permutationally-seamless integration. Right: fully-seamless integration.]
Directional can handle integration for every
In this example we demonstrate the results for lengthRatio
, and all fully seamless. Note that the density of the isolines increases with
Left to right: $N=2,4,7,11$. Top: field. Bottom: integer isolines.
It is possible to choose whether to round the seam jumps
Left to right (bottom is zoom-in of top): Field, rounding only seams (leading to the right singularity being offset), and rounding singularity function values directly.]
It is possible to constrain the functions to have linear relations between them, which reduce the degrees of freedom. This is done by inputting a matrix intData
through the field linRed
, and there are pre-made funtcions to set it, such as set_triangular_symmetry(int N)
.
Left to right: $6$-directional fields with a singularity, $N=6$ with only sign symmetry (the three lines don't always meet at all intersections), and the same with added triangular symmetry, where the intersections are enforced.
This subchapter demonstrates how we can create an actual polygonal mesh from the arrangement of isolines on the surface. This creates an arrangement of lines (in exact numbers) on every triangle, and stitches them together across triangles to get a complete conforming mesh. The new mesh is given in libhedra format18 of
The meshing unit is independent from the integration unit, and can be potentially used with external functions; one should fill the MeshFunctionIsolinesData
structure with the input, and call mesh_function_isolines()
.
The input is given as functions on the vertices of the whole (original) mesh in vertexNFunction
, where the user must supply two sparse matrices: orig2CutMat
to produce values per corner of the cut mesh, and the equivalent exactOrig2CutMatInteger
to do the same in exact numbers. This ensures that the values across seams are perfectly matching in rationals, and the robustness of the meshing.
There is a natural transition between the integrator and the mesher, which is done by calling setup_mesh_function_isolines()
, which translates the IntegratorData
to the meshing data.
The meshing function requires CGAL as a dependency, which is operated through the libigl CGAL dependency mechanism.
Left to right: polygonal meshes of the arrangements of isolines from the N=4,7,11 examples ($N=2$ is not yet supported) in Example 502. The screenshots are from Meshlab19]
This chapter introduces representations for fields that use high-order interpolants in order to work with fields with smoother bases in the finite-element paradigm.
This examples only works with IntrinsicFaceTangentBundle
by design, since the algorithm targets face-based fields.
Directional fields can be used with subdivision surfaces in a manner which is structure preserving. That is, the subdivision of a coarse directional field into a fine directional field subdivides a coarse gradient into a fine gradient, and its coarse curl into fine curl. The challenge of doing it for piecewise-constant fields is worked by 20. We optimize for a curl free field rawFieldCoarse
, subdivide it into rawFieldFine
using subdivide_field()
, and compute a seamless parameterization for both. The coarse field is optimized for being curl-free fast (using the PolyCurl optimization; see Example 303, and then the fine field is curl free by design, with the same singularities, which makes for an efficient process.
Top Left to right: coarse curl-reduced directional field, curl plot, and parameterization. Bottom: subdivided fine results.
Directional is a budding project, and there are many algorithms in the state-of-the-art that we look forward to implement, with the help of volunteer researchers and practitioners from the field. Prominent examples of desired implementations are:
-
Support for 3D fields, particularly Octahedral fields 21, both in tet meshes and with the boundary-element method.
-
A discrete exterior calculus framework.
-
Differential operators and Hodge decomposition.
-
Support for tensor fields.
-
Advanced and better visualization techniques.
Footnotes
-
Amir Vaxman, Marcel Campen, Olga Diamanti, Daniele Panozzo, David Bommes, Klaus Hildebrandt, Mirela Ben-Chen, Directional Field Synthesis, Design, and Processing, 2016. ↩ ↩2
-
Fernando de Goes, Mathieu Desbrun, Yiying Tong, Vector Field Processing on Triangle Meshes, 2016. ↩
-
Olga Diamanti, Amir Vaxman, Daniele Panozzo, Olga Sorkine-Hornung, Designing N-PolyVector Fields with Complex Polynomials, 2014. ↩ ↩2 ↩3 ↩4
-
Felix Knöppel, Keenan Crane, Ulrich Pinkall, and Peter Schröder, Globally Optimal Direction Fields, 2013. ↩ ↩2
-
Omri Azencot, Etienne Corman, Mirela Ben-Chen, Maks Ovsjanikov, Consistent Functional Cross Field Design for Mesh Quadrangulation, 2017. ↩
-
Christopher Brandt, Leonardo Scandolo, Elmar Eisemann, and Klaus Hildebrandt, Modeling n-Symmetry Vector Fields using Higher-Order Energies, 2018. ↩
-
Merel Meekes, Amir Vaxman, Unconventional Patterns on Surfaces, 2021. ↩
-
Nicolas Ray, Bruno Vallet, Wan Chiu Li, Bruno Lévy, N-Symmetry Direction Field Design, 2008. ↩ ↩2
-
David Bommes, Henrik Zimmer, Leif Kobbelt, Mixed-integer quadrangulation, 2009. ↩ ↩2
-
Daniele Panozzo, Enrico Puppo, Marco Tarini, Olga Sorkine-Hornung, Frame Fields: Anisotropic and Non-Orthogonal Cross Fields, 2014. ↩
-
Olga Diamanti, Amir Vaxman, Daniele Panozzo, Olga Sorkine-Hornung, Integrable PolyVector Fields, 2015. ↩
-
Yang Liu, Weiwei Xu, Jun Wang, Lifeng Zhu, Baining Guo, Falai Chen, Guoping Wang, General Planar Quadrilateral Mesh Design Using Conjugate Direction Field, 2008. ↩
-
Keenan Crane, Mathieu Desbrun, Peter Schröder, Trivial Connections on Discrete Surfaces, 2010. ↩ ↩2 ↩3
-
Amir Vaxman Directional Technical Reports: Seamless Integration, 2021. ↩
-
Merel Meekes, Amir Vaxman, [Unconventional Patterns on Surfaces] (https://webspace.science.uu.nl/~vaxma001/Unconventional_Patterns_on_Surfaces.pdf), 2021. ↩ ↩2 ↩3
-
Felix Kälberer, Matthias Nieser, Konrad Polthier, QuadCover - Surface Parameterization using Branched Coverings, 2007 ↩ ↩2
-
Ashish Myles, Nico Pietroni, Denis Zorin, Robust Field-aligned Global Parametrization, 2014. ↩
-
Amir Vaxman and others, libhedra: geometric processing and optimization of polygonal meshes. ↩
-
Paolo Cignoni, Marco Callieri, Massimiliano Corsini, Matteo Dellepiane, Fabio Ganovelli, Guido Ranzuglia, , MeshLab: an Open-Source Mesh Processing Tool. ↩
-
Bram Custers, Amir Vaxman, Subdivision Directional Fields, 2020. ↩
-
Justin Solomon, Amir Vaxman, David Bommes, Boundary Element Octahedral Fields in Volumes, 2017. ↩