Ray tracing is a set of techniques used in computer graphics to render photorealistic inside views of a scene in the Euclidean space . It considers that the light propagates along with straight lines. To ray trace a scene, rays are launched from each image pixel: where the light arrives. If these primary rays intersect some surface, the color is computed using direct and indirect light contributions. The direct illumination considers the rays coming directly from light sources and the indirect gathers the light rays bounced by other surfaces: secondary rays.
Ray tracing works intrinsically in space geometry since light travels along lines minimizing lengths (rays). Thus, unlike in rasterization, no tricks to compute information based on the global physical behavior of light are required. In particular, ray tracing could be extended to not necessarily be restricted to Euclidean geometry tracing straight lines. We believe Riemannian geometry is the right setting for this extension since it generalizes the required concepts: metric and rays.
Riemannian space construction has topology as a global geometric constraint. In , the authors define Riemannian shading to explore Nvidia RTX GPUs to visualize Non-Euclidean spaces endowed with non-trivial topology (dating back to Thurston’s geometrization conjecture). Instead, this work firstly focuses on time/user-dependent Riemannian metric constructions on the classic space to explore special effects like warping, mirages  and scene deformation . This implies spaces with trivial topology, but with generic Riemannian geometry. Then, we define Riemannian ray tracing to visualize scenes endowed with such effects. Our first results use graphs of functions and diffeomorphisms to construct the metrics, allowing the modeling of expressive effects. Finally, gaussian bump functions restrict these constructions providing local deformations .
That approach opens many geometric questions concerning new ways of representing rays and space. Thus, the goal of this paper is to establish the basis of a new line of research based on employing Riemannian geometry in ray tracing techniques. We believe curved rays can advance the state of the art in many areas, not restricted to rendering only.
2. Related works
In the literature, an expressive number of extensions and modifications of ray tracing were studied to increase efficiency, photorealism, functionality, among others. Almost all of these works focus on tracing straight line rays. Only a few approaches use nonlinear rays; we relate some of these.
Bar  used a nonlinear approach to ray trace smooth surfaces efficiently. The surfaces are considered to be bent flat planes, which are generated by a deformation of the space. Such deformation corresponds to a coordinate system, where the intersection problem considers the surface being flat, and the rays being curved.
In visualization of physics, Stam and Languénou  presented a ray tracing algorithm integrating perturbations with the basic equations from geometrical optics that govern the propagation of rays in a medium varying its index of refraction continuously. The curvature of the rays depends on the air’s index of refraction. Nonlinear ray tracing was also considered  in visualizing relativistic effects, the geometric behavior of nonlinear dynamical systems, and the movement of charged particles in a force field (e.g., electron movement).
In visualization of mathematics, ray tracing appeared recently in the visualization of Thurston geometries [2, 8, 5]. Berger et al.  used ray tracing to visualize Non-Euclidean spaces. However, this work does not provide real-time immersive and interactive visualization, and it was restricted to Euclidean and hyperbolic manifolds, where rays are modeled by straight lines.
Velho et. al [8, 5] introduced a framework, implemented on top of Nvidia RTX GPUs, for real-time immersive and interactive visualization of Euclidean, Hyperbolic, Spherical, Nil, Sol, and geometries: the six more interesting Thurston geometries. They introduced a novel ray tracing model for Riemannian manifolds focusing primarily on the topology of the underlying manifolds rather than the geometry. In this work, we consider only space which has a trivial topology, however, we endow it with complicated and interesting geometries.
3. Riemannian geometry
Riemannian geometry  was developed from the study of surfaces in . There is a natural metric in each tangent plane of a surface inherited from the scalar product of . This implies a distance measure in . Geodesics are curves in minimizing distance. We use the terms geodesic and ray interchangeably.
This section presents the definitions for Riemannian metrics in and their geodesics. Our ray tracing will operate in this setting. We follow a simplified version of the exposition in  which was inspired in the notation of Carmo .
3.1. Parametric manifold
A parametric -manifold is a topological space which can be parameterized by , that is, there is a homeomorphism from to — the chart. In the case of many parameterizations, the change of charts must smooth. Informally, the parametric manifold is a continuous deformation of the Euclidean space. Figure 1 shows a schematic view of a -manifold.
Let be a point in the -manifold , and be a parameterization of . The tangent space in
is the vector space generated by the vectors— the derivatives of the coordinates curves of x at . As admits a linear vector field structure, we can define a scalar product on each tangent space of . If such procedure is “smooth” the result is a metric on .
3.2. Riemannian metric
A Riemannian metric in is a map that attributes to each point a scalar product in the tangent space , such that in coordinates, , the function is smooth. This function does not depend on the coordinate system. We may use instead of .
Expressing two tangent vectors , , at a point , in terms of the basis , and , we obtain the scalar product:
The metric tensorin is fully determined by the matrix and generalizes the classical Euclidean inner product of . The pair is a Riemannian manifold. The chart x bridges and , so that the pull-back of the metric to results in a Riemannian manifold isometric to . This work investigates the design and visualization of Riemannian metrics on .
We now investigate the geodesics (rays) in the Riemannian manifold . In the Euclidean metric, these curves have null acceleration (second derivative). Here, the analogous concept of “acceleration” is defined using the covariant derivative. This is a way to calculate the derivative of tangent fields along curves in manifolds.
A geodesic in a Riemannian manifold is a curve with null covariant derivative:
This differs from the classical by the addition of , which includes the Christoffel symbols of .
To linearize System 3.2, we add new variables being the first derivatives , obtaining thus the geodesic flow of :
To compute the geodesic flow in a Riemannian manifold we need its Christoffel symbols at each point. These can be computed in terms of the metric coefficients and its derivatives  by the formula:
Where is the inverse matrix of .
3.4. Exponential map
Let be a point in — the observer. A key idea in the ray tracing algorithm is to trace rays from . This is also behind the exponential map: a natural application encoding all the rays leaving :
is the geodesic satisfying and . Intuitively, the exponential map takes a vector tangent at a point , and runs a geodesic from towards the direction until the parameter (see Figure 1).
A natural concern arises: Is well defined? The well-known Hopf–Rinow theorem  gives us an answer. It states that is defined at the whole tangent space iff the Riemannian manifold is a complete metric space. In our case, diffeomorphic to , then is defined in the whole space. Additionally, the theorem say that every point in can be connected to by a geodesic. Thus, is surjective implying that every scene in can be ray traced using this map.
However, may not be injective: when there exists points in connected to by many geodesics. Visualizing such manifolds using rasterization is infeasible.
The Hadamard theorem  states that the exponential map will be injective if, additionally, has non-positive sectional curvature. If we consider being the paraboloid, given by the graph of , the exponential map is not injective: there exists pairs of points in connected by many geodesics.
4. Riemannian shading and illumination
In computer graphics, Shading is the process of assigning a color to a pixel. Classic approaches to perform such tasks are not suited for Riemannian manifolds due to the nonlinear nature of their rays.
Consider a viewing -sphere centered in an observer point in a -manifold . We give a color for each ray direction in the the observer field by tracing a ray; carries the image. We call this procedure Riemannian shading. Specifically, the sphere is centered at the origin of . For each direction , we attribute a color by launching a ray from towards using the exponential map . If intersects a scene object, we define an RGB color. This is the Riemannian shading , where is a color space.
The Riemannian illumination of a point in an embedded surface comes from direct geodesics connecting to the light sources and indirect geodesics connecting other Riemannian-illuminated points to . The radiant intensity at is modeled using the Lambertian reflectance, which depends on the Riemannian metric, or more generally from a BRDF. For a point in , we compute the Riemannian shading using ray tracing and Riemannian illumination.
Classical ray tracing  approximates physical illumination. Our Riemannian ray tracing model can be also used to compute a Riemannian shading function for local or global illumination. This work uses pseudo-color based on properties of the space, such as a coordinated point, to define the Riemannian shading function.
4.1. Ray Marching
To launch a ray we need to solve the geodesic flow (Eq 3.3) and subsequently compute the intersection of the given ray with the scene objects. In general, our rays will be curved and are solutions of the geodesic flow which we have resorted to numerical integration methods.
We use Euler’s method for integration in the parameterization image, since there we use the Euclidean metric. The Euler’s numerical integration method approximates a ray starting at in the direction by a polygonal :
where is the integration step and is the ray satisfying and . We use Equation 3.3 to compute .
The polygonal is used to ray trace a scene in . When the scene is accurately rendered. For the ray-scene intersection, we check the for each segment given by the polygonal approximation given by Equation 4.1 (i.e., ray marching).
For now, we are testing the intersection with each segment sequentially. In future works, we pretend to use many segments of the polygonal to increase the GPU occupancy and optimize the intersection computations. This idea was suggested by Ingo Wald during GTC 2020.
We could use the Runge–Kutta method to approximate solutions of the geodesic flow. Instead, we consider the Euler method since it uses fewer computations giving GPU performance.
5. Examples of Riemannian metrics
This section presents two examples of Riemannian metrics in . The first is given by pushing the metric from the graph of three-dimensional functions, the other comes from deformations of given by diffeomorphisms.
5.1. Graph of a function
Let be a smooth function. We construct a Riemannian metric in which is the pullback of the metric of the graph of :
is parameterized in the obvious way by . The tangent space at a point is generated by the vectors:
Where is the partial derivative of in the standard direction . We use the Euclidean metric of to induce a Riemannian metric in .
Let be a point in , we compute the metric at . Replacing the expression of Equation 5.1 in , we obtain the coefficients of the Riemannian metric:
which are smooth functions of , then is Riemannian. In terms of measure, the metric makes the volume vary beautifully because its volume formula, used to compute integrals in , has the form . Then, it only matches with the standard if .
We present the geodesics of . A direct calculation using the software Maple implies in , with , , being all distinct, and , if . Substituting these in Equation 3.4, and using software Maple, we obtain a simple formula for the Christoffel symbols of :
Where is the second order partial derivative in the direction , and . Replacing the Christoffel symbols of in this Equation 3.4 we obtain.
For a concrete example, consider being the graph of the polynomial function . The Christoffel symbols of this parameterization are computed using Equation 5.3. Replacing the result in Equation 5.4 we obtain a geodesic flow. Figure 2 gives an inside view of the space endowed with the geometry pulled back from by the parameterization . In the left image, we positioned the camera looking towards direction . On the right image, the camera points at the direction . The scene is composed of a regular grid in distorted by .
Let , given by , be a diffeomorphism. We put a metric in based in the deformations provided by . Our base manifold in this case is parameterized by , and its associated base is given by:
We pullback the Euclidean metric of through the differential of , this will provide a new metric on .
The metric at a point is completely determined by the metric tensor . Thus subtituing Equation 5.5 in , we obtain the following formula:
As are smooth functions, is Riemannian.
We now compute the geodesics of . By Equation 3.2 we must first calculate the Christoffel symbols. Replacing the metric in Equation 5.6 in Christoffel symbols Equation 3.4 and using software Maple, we obtain a simple formula:
Where is the Jacobian matrix of , and is without the -line and the -column: the -minor of . The term in the equation is the -coefficient of the inverse of . Thus we obtain a beautiful formula.
Let be a diffeomorphism of . Then, the Christoffel symbols of are given by:
Theorem 5.1 says that the Christoffel symbols of can be computed using only the Jacobian and Hessian operators of the parameterization . For the best of our knowledge, this is the first time such a formula appears. Equation 5.7 will be very useful when the parameterization be a composition of parameterizations.
Replacing the Christoffel symbols in the Eq. 3.2, we obtain its geodesic flow:
For an explicit example of a deformation of , we consider the following twist map
The Christoffel symbols of this parameterization are computed using Equation 5.7. Replacing the result in Equation 5.8 we obtain a geodesic flow. We give an inside view of this “twisted” by tracing rays using Riemannian shading. Figure 3 shows an immersive view of a regular grid in distorted by .
6. Deforming the space
This section describes two convenient ways of accumulating local deformations of . The first is an extension of Subsection 5.1 consisting of summing functions. The second deals with the composition of deformations. In both cases, the functions and deformations are modeled using Gaussians. Many techniques in Computer Graphics depends on these functions due to their flexibility in modeling and smoothness.
6.1. Local metric deformations
We model graphs and deformations using the three-dimensional Gaussian function
Where is the amplitude, is the center, and , , are the spreads.
Using Subsection 5.1, we visualize the graph . Figure 4 shows inside views of by varying the parameters of . The first image presents the case of zero amplitude. The second image illustrates a small deformation given by a small amplitude and spreads values. The third and fourth images have small amplitude with lager spreads, and vice-versa. In all the cases, the scene is composed of a regular grid and the ray tracing has a small (black) fog component. This is why in the fourth image there is a black region in the center of the deformation.
We present a formula to deform a neighborhood of a point in . Consider a Gaussian function centered at , choosing small spreads values for , the formula
deforms the neighborhood of at the direction .
For a more refined formula, replace by a vector field . The resulted formula deforms a neighborhood of following the instructions of .
6.2. Accumulating deformations
We present two ways of accumulating deformations: summing the functions or the deformations presented above, and composing diffeomorphism.
Let be functions taking values in . The graph of the function is given by
In other words, to compute the metric we only need to store the partial derivatives of each function . The pair is the desired Riemannian manifold.
We now compute the geodesic flow of . The geodesic Equation 3.2 is given by the Christoffel symbols of , which are presented as by Eq. 5.3 in terms of . Therefore, the symbols can be computed using the first and second partial derivatives of the functions ’s.
To show the effects in visualizing graphs of summed functions, consider a regular grid embedded in (Left image in Figure 5). Adding a unique Gaussian function in the summation (middle image in Figure 5). Summing two Gaussian functions (right image in Figure 5).
Let be a sequence of diffeomorphisms of , we pretend to compute the Riemannian metric on that carries the deformations of the composition . Then we derive the geodesic flow of the resulting Riemannian manifold for a possible ray tracing.
Equation 5.6 states that the metric is given by , where . In other words, to compute the metric coefficients , we only need the Jacobian matrix of which is easily obtained using the chain rule formula .
To compute the geodesic flow of the Riemannian manifold , we need its Christoffel symbols, which are provided by the Jacobian and the Hessian operator of (from Theorem 5.1) using the formula . Specifically, to compute , we invert the Jacobian of using the formula
It is not direct to compute the Hessian component of , where is the -coordinate of . Firstly, we consider the composition of only two diffeomorphisms , then we extend it to the composition of diffeomorphisms.
The Hessian matrix of the coordinate of , that is , is given by
The general formula is obtained through an inductive argument. We simplify the notation by considering , for each . Therefore
From a computational point of view, this is an important formula because it computes in terms of the Jacobians and Hessians of each diffeomorphism separately. Therefore, during the computations only the partial derivatives of each diffeomorphism must be stored, the partial derivatives of the compositions are derived from them.
-  (1986) Ray tracing deformed surfaces. ACM SIGGRAPH Computer Graphics 20 (4), pp. 287–296. Cited by: §1, §2.
-  (2014) An image-space algorithm for immersive views in 3-manifolds and orbifolds. Visual Computer. Cited by: §2.
-  (1992) Riemannian geometry. Birkhäuser. Cited by: §3.3, §3.4, §3.4, §3, §3.
-  (1995) Nonlinear ray tracing: visualizing strange worlds. The Visual Computer 11 (5), pp. 263–274. Cited by: §2.
-  (2020) Visualization of nil, sol, and sl2 geometries. In In preparation, Cited by: §2, §2, §3.
-  (2019)(Website) External Links: Cited by: §1.
-  (1996) Ray tracing in non-constant media. In Rendering Techniques’ 96, pp. 225–234. Cited by: §1, §2.
-  (2020) Immersive visualization of the classical non-euclidean spaces using real-time ray tracing in vr. In In Proceedings of the 46th Graphics Interface, Cited by: §1, §2, §2.
-  (1979) An improved illumination model for shaded display. In ACM SIGGRAPH Computer Graphics, Vol. 13, pp. 14. Cited by: §1, §4.