As mentioned in the last post, today’s post looks at how to draw hyperbolic geometry using the Poincaré disk model inside AutoCAD.
This is an interesting exercise, but probably won’t ultimately help us with the hyperbolic tessellations we’re aiming to create: it’s interesting as it will end up with us having a mechanism for mapping hyperbolic geometry onto an arbitrary circle inside AutoCAD, but as we’re probably going to end up creating straight-lined segments for the polygons we use to create our pattern, it’s mostly an exercise (and will probably be largely irrelevant to our end result). The tough part actually comes in the next post, when we move on to the tessellation algorithm.
To get started, I managed to find a straightforward C++ implementation of a set of functions to create hyperbolic geometry on a Poincaré disk. This library was developed by by David Coeurjolly, Research Director for the team focused on Multiresolution, Discrete and Combinatorial Models at LIRIS (Laboratoire d'InfoRmatique en Image et Systèmes d'information), part of the CNRS (Centre National de la Recherche Scientifique) in France. David’s implementation makes use of the Cairo graphics library, so there was a modest amount of work needed to generate geometry inside AutoCAD, instead. Indeed this was the main focus of my efforts: I fully admit I haven’t spent much time looking at the underlying mathematics (in particular, please don’t ask me to explain how the radius is calculated in the ComputeCircleParameters() method… ;-).
Here’s the resultant C# file containing the library class I’ve named HyperbolicEngine (for want of a better term):
using Autodesk.AutoCAD.ApplicationServices;
using Autodesk.AutoCAD.DatabaseServices;
using Autodesk.AutoCAD.Geometry;
using System;
namespace HyperbolicGeometry
{
public class HyperbolicEngine
{
// We'll hold some AutoCAD objects to provide quick access
// to the drawing database
Transaction _tr;
Database _db;
BlockTableRecord _btr;
// We'll also hold information about the target circle
double _radius;
Point3d _center;
Matrix3d _disp;
// Whether we'll draw curved or straight lines
bool _straight;
// Target layers for points, supporting lines and our
// actual geometry
ObjectId _ptLay, _lnLay, _hgLay;
// Constructor
public HyperbolicEngine(
Transaction tr,
Database db,
Circle cir,
bool straight,
ObjectId? ptLay = null,
ObjectId? lnLay = null,
ObjectId? hgLay = null
)
{
_tr = tr;
_db = db;
if (cir == null)
{
_radius = 400;
_center = Point3d.Origin;
}
else
{
_center = cir.Center;
_radius = cir.Radius;
}
_disp = Matrix3d.Displacement(_center - Point3d.Origin);
_btr =
(BlockTableRecord)tr.GetObject(
db.CurrentSpaceId, OpenMode.ForWrite
);
_straight = straight;
_ptLay = ptLay.HasValue ? ptLay.Value : db.Clayer;
_lnLay = lnLay.HasValue ? lnLay.Value : db.Clayer;
_hgLay = hgLay.HasValue ? hgLay.Value : db.Clayer;
}
// Functions to map a point from the unit disk to WCS
private Point3d WcsPoint(double x, double y)
{
Point3d pt = new Point3d(cX(x), cY(y), 0);
return pt.TransformBy(_disp);
}
private double cX(double x)
{
return x * _radius;
}
private double cY(double y)
{
return y * _radius;
}
/*
* Compute the circle containing points a and b, and
* orthogonal to the unit circle (a.k.a shortest path between
* a and b in the hyperbolic sense).
*
* Return true if the output is a circle. False if the
* shortest path is a straight line (diameter of the unit
* circle).
*/
private static bool ComputeCircleParameters(
Point2d a, Point2d b,
out double cx, out double cy, out double radius
)
{
double ax = a.X, ay = a.Y, bx = b.X, by = b.Y;
cx = 0; cy = 0; radius = 0;
if (Math.Abs(ax * by - ay * bx) < Tolerance.Global.EqualPoint)
return false;
cx =
(-1 / 2.0 *
(ay * bx * bx + ay * by * by -
(ax * ax + ay * ay + 1) * by + ay
)
/ (ax * by - ay * bx)
);
cy =
(1 / 2.0 *
(ax * bx * bx + ax * by * by -
(ax * ax + ay * ay + 1) * bx + ax
)
/ (ax * by - ay * bx)
);
radius =
-1 / 2.0 * Math.Sqrt(
ax * ax * ax * ax * bx * bx +
ax * ax * ax * ax * by * by -
2 * ax * ax * ax * bx * bx * bx -
2 * ax * ax * ax * bx * by * by +
2 * ax * ax * ay * ay * bx * bx +
2 * ax * ax * ay * ay * by * by -
2 * ax * ax * ay * bx * bx * by -
2 * ax * ax * ay * by * by * by +
ax * ax * bx * bx * bx * bx +
2 * ax * ax * bx * bx * by * by +
ax * ax * by * by * by * by -
2 * ax * ay * ay * bx * bx * bx -
2 * ax * ay * ay * bx * by * by +
ay * ay * ay * ay * bx * bx +
ay * ay * ay * ay * by * by -
2 * ay * ay * ay * bx * bx * by -
2 * ay * ay * ay * by * by * by +
ay * ay * bx * bx * bx * bx +
2 * ay * ay * bx * bx * by * by +
ay * ay * by * by * by * by -
2 * ax * ax * ax * bx -
2 * ax * ax * ay * by +
4 * ax * ax * bx * bx -
2 * ax * ay * ay * bx +
8 * ax * ay * bx * by -
2 * ax * bx * bx * bx -
2 * ax * bx * by * by -
2 * ay * ay * ay * by +
4 * ay * ay * by * by -
2 * ay * bx * bx * by -
2 * ay * by * by * by +
ax * ax - 2 * ax * bx +
ay * ay - 2 * ay * by +
bx * bx + by * by
)
/ (ax * by - ay * bx);
if (radius < 0.0)
radius =
1 / 2.0 * Math.Sqrt(
ax * ax * ax * ax * bx * bx +
ax * ax * ax * ax * by * by -
2 * ax * ax * ax * bx * bx * bx -
2 * ax * ax * ax * bx * by * by +
2 * ax * ax * ay * ay * bx * bx +
2 * ax * ax * ay * ay * by * by -
2 * ax * ax * ay * bx * bx * by -
2 * ax * ax * ay * by * by * by +
ax * ax * bx * bx * bx * bx +
2 * ax * ax * bx * bx * by * by +
ax * ax * by * by * by * by -
2 * ax * ay * ay * bx * bx * bx -
2 * ax * ay * ay * bx * by * by +
ay * ay * ay * ay * bx * bx +
ay * ay * ay * ay * by * by -
2 * ay * ay * ay * bx * bx * by -
2 * ay * ay * ay * by * by * by +
ay * ay * bx * bx * bx * bx +
2 * ay * ay * bx * bx * by * by +
ay * ay * by * by * by * by -
2 * ax * ax * ax * bx -
2 * ax * ax * ay * by +
4 * ax * ax * bx * bx -
2 * ax * ay * ay * bx +
8 * ax * ay * bx * by -
2 * ax * bx * bx * bx -
2 * ax * bx * by * by -
2 * ay * ay * ay * by +
4 * ay * ay * by * by -
2 * ay * bx * bx * by -
2 * ay * by * by * by +
ax * ax - 2 * ax * bx +
ay * ay - 2 * ay * by +
bx * bx + by * by
)
/ (ax * by - ay * bx);
return true;
}
// Draw the unit circle
public void DrawUnitCircle()
{
if (_tr != null && _btr != null)
{
Circle c = new Circle(_center, Vector3d.ZAxis, _radius);
c.LayerId = _hgLay;
_btr.AppendEntity(c);
_tr.AddNewlyCreatedDBObject(c, true);
}
}
// Draw a point in the Poincare disc
public void DrawPoint(Point2d a)
{
if (_tr != null && _btr != null)
{
DBPoint p = new DBPoint(WcsPoint(a.X, a.Y));
p.LayerId = _ptLay;
_btr.AppendEntity(p);
_tr.AddNewlyCreatedDBObject(p, true);
}
}
// Draw a hyperbolic line through a and b
public void DrawSupportingLine(
Point2d a, Point2d b, bool withPoints
)
{
if (withPoints)
{
DrawPoint(a);
DrawPoint(b);
}
double ax = a.X, ay = a.Y, bx = b.X, by = b.Y;
double cx, cy, r;
bool result =
ComputeCircleParameters(a, b, out cx, out cy, out r);
Entity ent;
if (!result)
{
// We project a and b points onto the unit circle
double theta = Math.Atan2(ay, ax);
double theta2 = Math.Atan2(by, bx);
ent =
new Line(
WcsPoint(Math.Cos(theta), Math.Sin(theta)),
WcsPoint(Math.Cos(theta2), Math.Sin(theta2))
);
}
else
{
ent =
new Circle(WcsPoint(cx, cy), Vector3d.ZAxis, r * _radius);
}
if (ent != null)
{
ent.LayerId = _lnLay;
if (_tr != null && _btr != null)
{
_btr.AppendEntity(ent);
_tr.AddNewlyCreatedDBObject(ent, true);
}
else
{
ent.Dispose();
}
}
}
// Draw a hyperbolic edge between a and b
public void DrawEdge(
Point2d a, Point2d b, bool withPoints,
bool straightSegments = false
)
{
if (withPoints)
{
DrawPoint(a);
DrawPoint(b);
}
double ax = a.X, ay = a.Y, bx = b.X, by = b.Y;
double cx = 0.0, cy = 0.0, r = 0.0;
bool result =
_straight ? false :
ComputeCircleParameters(a, b, out cx, out cy, out r);
// Near-aligned points -> straight segment
Entity ent = null;
if (!result)
{
ent = new Line(WcsPoint(ax, ay), WcsPoint(bx, by));
}
else
{
double theta = -Math.Atan2(-ay + cy, ax - cx);
double theta2 = -Math.Atan2(-by + cy, bx - cx);
// We recenter the angles to [0, 2Pi]
while (theta < 0)
{
theta += 2 * Math.PI;
}
while (theta2 < 0)
{
theta2 += 2 * Math.PI;
}
Point3d cen = WcsPoint(cx, cy);
double delta = Math.Abs(theta - theta2);
if (
(theta > theta2 && delta > Math.PI) ||
(theta < theta2 && delta < Math.PI)
)
{
ent = new Arc(cen, r * _radius, theta, theta2);
}
else
{
ent = new Arc(cen, r * _radius, theta2, theta);
}
}
if (ent != null)
{
ent.LayerId = _hgLay;
if (_tr != null && _btr != null)
{
_btr.AppendEntity(ent);
_tr.AddNewlyCreatedDBObject(ent, true);
}
else
{
ent.Dispose();
}
}
}
// Draw a hyperbolic triangle with vertices (a,b,c)
public void DrawTriangle(
Point2d a, Point2d b, Point2d c,
bool withPoints = false, bool withSupportLines = false
)
{
DrawPolygon(
new Point2d[] { a, b, c }, withPoints, withSupportLines
);
}
// Draw a hyperbolic polygon with the provided vertices
public void DrawPolygon(
Point2d[] pts, bool withPoints = false,
bool withSupportLines = false
)
{
if (withPoints)
{
foreach (Point2d pt in pts)
{
DrawPoint(pt);
}
}
int numPts = pts.Length;
if (withSupportLines && !_straight)
{
for (int i = 0; i < numPts; i++)
{
DrawSupportingLine(
pts[i], pts[(i + 1) % numPts], withPoints
);
}
}
for (int i = 0; i < numPts; i++)
{
DrawEdge(
pts[i], pts[(i + 1) % numPts], withPoints
);
}
}
}
}
Here’s the command class implementation:
using Autodesk.AutoCAD.ApplicationServices;
using Autodesk.AutoCAD.Colors;
using Autodesk.AutoCAD.DatabaseServices;
using Autodesk.AutoCAD.EditorInput;
using Autodesk.AutoCAD.Geometry;
using Autodesk.AutoCAD.Runtime;
using System;
namespace HyperbolicGeometry
{
public class Commands
{
[CommandMethod("HGT")]
public static void CreateHyperbolicTriangles()
{
Document doc =
Application.DocumentManager.MdiActiveDocument;
Database db = doc.Database;
Transaction tr =
doc.TransactionManager.StartTransaction();
using (tr)
{
LayerTable lt =
(LayerTable)tr.GetObject(
db.LayerTableId, OpenMode.ForWrite
);
// Points and supporting lines will be on layers that
// are grey and initially frozen
ObjectId ptLay = CreateLayer(tr, lt, "HG_PTS", 9, false);
ObjectId lnLay = CreateLayer(tr, lt, "HG_LNS", 9, false);
ObjectId htLay = CreateLayer(tr, lt, "HG_CNT", 6, true);
// Create our hyperbolic engine
HyperbolicEngine he =
new HyperbolicEngine(
tr, db, null, false, ptLay, lnLay, htLay
);
he.DrawUnitCircle();
he.DrawTriangle(
new Point2d(0.8, -0.3),
new Point2d(0.2, -0.6),
new Point2d(-0.61, -0.6),
true, true
);
tr.Commit();
}
}
[CommandMethod("HGS")]
public static void CreateHyperbolicSpiro()
{
Document doc =
Application.DocumentManager.MdiActiveDocument;
Database db = doc.Database;
Transaction tr =
doc.TransactionManager.StartTransaction();
using (tr)
{
LayerTable lt =
(LayerTable)tr.GetObject(
db.LayerTableId, OpenMode.ForWrite
);
// Points and supporting lines will be on layers that
// are grey and initially frozen
ObjectId ptLay = CreateLayer(tr, lt, "HG_PTS", 9, false);
ObjectId lnLay = CreateLayer(tr, lt, "HG_LNS", 9, false);
ObjectId htLay = CreateLayer(tr, lt, "HG_CNT", 6, true);
// Create our hyperbolic engine
HyperbolicEngine he =
new HyperbolicEngine(
tr, db, null, false, ptLay, lnLay, htLay
);
he.DrawUnitCircle();
int max = 13;
for (int k = 0; k < (max - 1); k++)
{
he.DrawTriangle(
new Point2d(1, 0),
new Point2d(
Math.Cos(k * 134 / (double)max),
Math.Sin(k * 134 / (double)max)
),
new Point2d(
Math.Cos((k + 1) * 134 / (double)max),
Math.Sin((k + 1) * 134 / (double)max)
),
true, true
);
}
tr.Commit();
}
}
private static ObjectId CreateLayer(
Transaction tr, LayerTable lt,
string lname, short col, bool on
)
{
if (lt.Has(lname))
return lt[lname];
LayerTableRecord ltr = new LayerTableRecord();
ltr.Color = Color.FromColorIndex(ColorMethod.ByAci, col);
ltr.Name = lname;
ltr.IsFrozen = !on;
ObjectId layId = lt.Add(ltr);
tr.AddNewlyCreatedDBObject(ltr, true);
return layId;
}
}
}
The code implements two commands, both of which have been lifted fairly directly from David’s test cases.
- HGT – creates a simple hyperbolic triangle
- HGS – creates a more complex set of hyperbolic triangles
Here are the results of running the HGT command:
And if you turn on the layer showing the supporting lines (HG_LNS), you’ll see the circles used to generated these arcs:
The HGS command creates more interesting geometry:
Which, of course, also has more interesting supporting lines along with it:
You can certainly see similarities between geometry on a Poincaré disk and the hypocycloids generated by Spiro.
Next time we’ll look at hyperbolic tessellations and how to generate them inside AutoCAD.