Volumes and surface areas of solids
Measure any closed mesh: surface area from its triangles, volume and centre of mass from signed tetrahedra, and why the result goes wrong for open meshes or coordinates far from the origin.
Surface area
Volume from signed tetrahedra
var shapes = new (string Name, Mesh M)[]{ ("box 4 x 3 x 2.5", Mesh.Box(4, 3, 2.5)), ("L-shaped room", Mesh.Prism(new[] { (0.0, 0.0), (4.0, 0.0), (4.0, 2.0), (2.0, 2.0), (2.0, 4.0), (0.0, 4.0) }, 3)), ("pyramid", Mesh.Pyramid(4, 3))}; foreach (var (name, m) in shapes) Console.WriteLine($"{name}: area {Area(m):0.###} m2, volume {Volume(m, new V3(0, 0, 0)):0.###} m3, from O = (50, -20, 7): {Volume(m, new V3(50, -20, 7)):0.###} m3"); static double Area(Mesh m) => m.F.Sum(t => m.Normal(t).Length / 2); static double Volume(Mesh m, V3 o) => m.F.Sum(t => V3.Dot(m.V[t.A] - o, V3.Cross(m.V[t.B] - o, m.V[t.C] - o)) / 6); sealed class Mesh{ public readonly List<V3> V = new(); public readonly List<(int A, int B, int C)> F = new(); // counter-clockwise seen from outside public static Mesh Box(double lx, double ly, double lz) => Prism(new[] { (0.0, 0.0), (lx, 0.0), (lx, ly), (0.0, ly) }, lz); // A counter-clockwise footprint extruded up by h (floor and roof fanned from corner 0) public static Mesh Prism((double X, double Y)[] poly, double h) { var m = new Mesh(); int n = poly.Length; foreach (var (x, y) in poly) m.V.Add(new V3(x, y, 0)); foreach (var (x, y) in poly) m.V.Add(new V3(x, y, h)); for (int k = 1; k < n - 1; k++) { m.F.Add((0, k + 1, k)); // floor, facing down m.F.Add((n, n + k, n + k + 1)); // roof, facing up } for (int i = 0; i < n; i++) { int j = (i + 1) % n; m.F.Add((i, j, j + n)); // each wall as two triangles m.F.Add((i, j + n, i + n)); } return m; } public static Mesh Pyramid(double b, double h) { var m = new Mesh(); m.V.AddRange(new V3[] { new(0, 0, 0), new(b, 0, 0), new(b, b, 0), new(0, b, 0), new(b / 2, b / 2, h) }); m.F.AddRange(new[] { (0, 1, 4), (1, 2, 4), (2, 3, 4), (3, 0, 4), (0, 2, 1), (0, 3, 2) }); return m; } // (b - a) x (c - a): points out of the solid, and its length is twice the triangle's area public V3 Normal((int A, int B, int C) t) => V3.Cross(V[t.B] - V[t.A], V[t.C] - V[t.A]);} readonly record struct V3(double X, double Y, double Z){ public static V3 operator +(V3 a, V3 b) => new(a.X + b.X, a.Y + b.Y, a.Z + b.Z); public static V3 operator -(V3 a, V3 b) => new(a.X - b.X, a.Y - b.Y, a.Z - b.Z); public static V3 operator *(V3 v, double k) => new(v.X * k, v.Y * k, v.Z * k); public static double Dot(V3 a, V3 b) => a.X * b.X + a.Y * b.Y + a.Z * b.Z; public static V3 Cross(V3 a, V3 b) => new(a.Y * b.Z - a.Z * b.Y, a.Z * b.X - a.X * b.Z, a.X * b.Y - a.Y * b.X); public double Length => Math.Sqrt(Dot(this, this)); public V3 Normalized() => this * (1 / Length); static string F(double v) => Math.Abs(v) < 5e-4 ? "0" : v.ToString("0.###"); public override string ToString() => $"({F(X)}, {F(Y)}, {F(Z)})";}box 4 x 3 x 2.5: area 59 m2, volume 30 m3, from O = (50, -20, 7): 30 m3 L-shaped room: area 72 m2, volume 36 m3, from O = (50, -20, 7): 36 m3 pyramid: area 44.844 m2, volume 16 m3, from O = (50, -20, 7): 16 m3
Centre of mass
var footprint = new[] { (0.0, 0.0), (4.0, 0.0), (4.0, 2.0), (2.0, 2.0), (2.0, 4.0), (0.0, 4.0) };var slab = Mesh.Prism(footprint, 0.25); // an L-shaped precast slab double vol = 0;var moment = new V3(0, 0, 0);foreach (var t in slab.F){ V3 a = slab.V[t.A], b = slab.V[t.B], c = slab.V[t.C]; // tetrahedron (origin, a, b, c) double v = V3.Dot(a, V3.Cross(b, c)) / 6; vol += v; moment += (a + b + c) * (v / 4); // its centroid is (0 + a + b + c) / 4} Console.WriteLine($"volume {vol:0.###} m3, mass {vol * 2400:0} kg at 2400 kg/m3");Console.WriteLine($"centre of mass {moment * (1 / vol)}"); sealed class Mesh{ public readonly List<V3> V = new(); public readonly List<(int A, int B, int C)> F = new(); // counter-clockwise seen from outside public static Mesh Box(double lx, double ly, double lz) => Prism(new[] { (0.0, 0.0), (lx, 0.0), (lx, ly), (0.0, ly) }, lz); // A counter-clockwise footprint extruded up by h (floor and roof fanned from corner 0) public static Mesh Prism((double X, double Y)[] poly, double h) { var m = new Mesh(); int n = poly.Length; foreach (var (x, y) in poly) m.V.Add(new V3(x, y, 0)); foreach (var (x, y) in poly) m.V.Add(new V3(x, y, h)); for (int k = 1; k < n - 1; k++) { m.F.Add((0, k + 1, k)); // floor, facing down m.F.Add((n, n + k, n + k + 1)); // roof, facing up } for (int i = 0; i < n; i++) { int j = (i + 1) % n; m.F.Add((i, j, j + n)); // each wall as two triangles m.F.Add((i, j + n, i + n)); } return m; } public static Mesh Pyramid(double b, double h) { var m = new Mesh(); m.V.AddRange(new V3[] { new(0, 0, 0), new(b, 0, 0), new(b, b, 0), new(0, b, 0), new(b / 2, b / 2, h) }); m.F.AddRange(new[] { (0, 1, 4), (1, 2, 4), (2, 3, 4), (3, 0, 4), (0, 2, 1), (0, 3, 2) }); return m; } // (b - a) x (c - a): points out of the solid, and its length is twice the triangle's area public V3 Normal((int A, int B, int C) t) => V3.Cross(V[t.B] - V[t.A], V[t.C] - V[t.A]);} readonly record struct V3(double X, double Y, double Z){ public static V3 operator +(V3 a, V3 b) => new(a.X + b.X, a.Y + b.Y, a.Z + b.Z); public static V3 operator -(V3 a, V3 b) => new(a.X - b.X, a.Y - b.Y, a.Z - b.Z); public static V3 operator *(V3 v, double k) => new(v.X * k, v.Y * k, v.Z * k); public static double Dot(V3 a, V3 b) => a.X * b.X + a.Y * b.Y + a.Z * b.Z; public static V3 Cross(V3 a, V3 b) => new(a.Y * b.Z - a.Z * b.Y, a.Z * b.X - a.X * b.Z, a.X * b.Y - a.Y * b.X); public double Length => Math.Sqrt(Dot(this, this)); public V3 Normalized() => this * (1 / Length); static string F(double v) => Math.Abs(v) < 5e-4 ? "0" : v.ToString("0.###"); public override string ToString() => $"({F(X)}, {F(Y)}, {F(Z)})";}volume 3 m3, mass 7200 kg at 2400 kg/m3 centre of mass (1.667, 1.667, 0.125)
Keep O close to the solid
var box = Mesh.Box(4, 3, 2.5);var site = new V3(512_340.37, 6_401_220.81, 35.2); // survey coordinatesvar placed = new Mesh();placed.V.AddRange(box.V.Select(p => p + site));placed.F.AddRange(box.F); Console.WriteLine($"near the origin: {Volume(box, new V3(0, 0, 0))} m3");Console.WriteLine($"at site, O = origin: {Volume(placed, new V3(0, 0, 0))} m3");Console.WriteLine($"at site, O = a corner: {Volume(placed, placed.V[0])} m3"); static double Volume(Mesh m, V3 o) => m.F.Sum(t => V3.Dot(m.V[t.A] - o, V3.Cross(m.V[t.B] - o, m.V[t.C] - o)) / 6); sealed class Mesh{ public readonly List<V3> V = new(); public readonly List<(int A, int B, int C)> F = new(); // counter-clockwise seen from outside public static Mesh Box(double lx, double ly, double lz) => Prism(new[] { (0.0, 0.0), (lx, 0.0), (lx, ly), (0.0, ly) }, lz); // A counter-clockwise footprint extruded up by h (floor and roof fanned from corner 0) public static Mesh Prism((double X, double Y)[] poly, double h) { var m = new Mesh(); int n = poly.Length; foreach (var (x, y) in poly) m.V.Add(new V3(x, y, 0)); foreach (var (x, y) in poly) m.V.Add(new V3(x, y, h)); for (int k = 1; k < n - 1; k++) { m.F.Add((0, k + 1, k)); // floor, facing down m.F.Add((n, n + k, n + k + 1)); // roof, facing up } for (int i = 0; i < n; i++) { int j = (i + 1) % n; m.F.Add((i, j, j + n)); // each wall as two triangles m.F.Add((i, j + n, i + n)); } return m; } public static Mesh Pyramid(double b, double h) { var m = new Mesh(); m.V.AddRange(new V3[] { new(0, 0, 0), new(b, 0, 0), new(b, b, 0), new(0, b, 0), new(b / 2, b / 2, h) }); m.F.AddRange(new[] { (0, 1, 4), (1, 2, 4), (2, 3, 4), (3, 0, 4), (0, 2, 1), (0, 3, 2) }); return m; } // (b - a) x (c - a): points out of the solid, and its length is twice the triangle's area public V3 Normal((int A, int B, int C) t) => V3.Cross(V[t.B] - V[t.A], V[t.C] - V[t.A]);} readonly record struct V3(double X, double Y, double Z){ public static V3 operator +(V3 a, V3 b) => new(a.X + b.X, a.Y + b.Y, a.Z + b.Z); public static V3 operator -(V3 a, V3 b) => new(a.X - b.X, a.Y - b.Y, a.Z - b.Z); public static V3 operator *(V3 v, double k) => new(v.X * k, v.Y * k, v.Z * k); public static double Dot(V3 a, V3 b) => a.X * b.X + a.Y * b.Y + a.Z * b.Z; public static V3 Cross(V3 a, V3 b) => new(a.Y * b.Z - a.Z * b.Y, a.Z * b.X - a.X * b.Z, a.X * b.Y - a.Y * b.X); public double Length => Math.Sqrt(Dot(this, this)); public V3 Normalized() => this * (1 / Length); static string F(double v) => Math.Abs(v) < 5e-4 ? "0" : v.ToString("0.###"); public override string ToString() => $"({F(X)}, {F(Y)}, {F(Z)})";}near the origin: 30 m3 at site, O = origin: 29.99652267806232 m3 at site, O = a corner: 30 m3
| Quantity | From a mesh |
|---|---|
Show a hint
Show the solution
var footprint = new[] { (0.0, 0.0), (4.0, 0.0), (4.0, 2.0), (2.0, 2.0), (2.0, 4.0), (0.0, 4.0) };var room = Mesh.Prism(footprint, 3); double fromMesh = room.F.Sum(t => V3.Dot(room.V[t.A], V3.Cross(room.V[t.B], room.V[t.C])) / 6); double shoelace = 0;for (int i = 0; i < footprint.Length; i++){ var (x1, y1) = footprint[i]; var (x2, y2) = footprint[(i + 1) % footprint.Length]; shoelace += (x1 * y2 - x2 * y1) / 2;} double walls = room.F.Select(room.Normal).Where(n => Math.Abs(n.Z) < 1e-9).Sum(n => n.Length / 2); Console.WriteLine($"volume from the mesh: {fromMesh:0.###} m3");Console.WriteLine($"floor area x height: {shoelace} x 3 = {shoelace * 3} m3");Console.WriteLine($"wall area: {walls:0.###} m2"); sealed class Mesh{ public readonly List<V3> V = new(); public readonly List<(int A, int B, int C)> F = new(); // counter-clockwise seen from outside public static Mesh Box(double lx, double ly, double lz) => Prism(new[] { (0.0, 0.0), (lx, 0.0), (lx, ly), (0.0, ly) }, lz); // A counter-clockwise footprint extruded up by h (floor and roof fanned from corner 0) public static Mesh Prism((double X, double Y)[] poly, double h) { var m = new Mesh(); int n = poly.Length; foreach (var (x, y) in poly) m.V.Add(new V3(x, y, 0)); foreach (var (x, y) in poly) m.V.Add(new V3(x, y, h)); for (int k = 1; k < n - 1; k++) { m.F.Add((0, k + 1, k)); // floor, facing down m.F.Add((n, n + k, n + k + 1)); // roof, facing up } for (int i = 0; i < n; i++) { int j = (i + 1) % n; m.F.Add((i, j, j + n)); // each wall as two triangles m.F.Add((i, j + n, i + n)); } return m; } public static Mesh Pyramid(double b, double h) { var m = new Mesh(); m.V.AddRange(new V3[] { new(0, 0, 0), new(b, 0, 0), new(b, b, 0), new(0, b, 0), new(b / 2, b / 2, h) }); m.F.AddRange(new[] { (0, 1, 4), (1, 2, 4), (2, 3, 4), (3, 0, 4), (0, 2, 1), (0, 3, 2) }); return m; } // (b - a) x (c - a): points out of the solid, and its length is twice the triangle's area public V3 Normal((int A, int B, int C) t) => V3.Cross(V[t.B] - V[t.A], V[t.C] - V[t.A]);} readonly record struct V3(double X, double Y, double Z){ public static V3 operator +(V3 a, V3 b) => new(a.X + b.X, a.Y + b.Y, a.Z + b.Z); public static V3 operator -(V3 a, V3 b) => new(a.X - b.X, a.Y - b.Y, a.Z - b.Z); public static V3 operator *(V3 v, double k) => new(v.X * k, v.Y * k, v.Z * k); public static double Dot(V3 a, V3 b) => a.X * b.X + a.Y * b.Y + a.Z * b.Z; public static V3 Cross(V3 a, V3 b) => new(a.Y * b.Z - a.Z * b.Y, a.Z * b.X - a.X * b.Z, a.X * b.Y - a.Y * b.X); public double Length => Math.Sqrt(Dot(this, this)); public V3 Normalized() => this * (1 / Length); static string F(double v) => Math.Abs(v) < 5e-4 ? "0" : v.ToString("0.###"); public override string ToString() => $"({F(X)}, {F(Y)}, {F(Z)})";}volume from the mesh: 36 m3 floor area x height: 12 x 3 = 36 m3 wall area: 48 m2
Remember
Surface area is the sum of ½ |(b − a) × (c − a)| over the triangles. Volume is the sum of signed tetrahedra (a − O) · ((b − O) × (c − O)) / 6; the parts outside the solid cancel. It needs a closed mesh with consistent winding; a negative total means the mesh is inside out. The centre of mass is the volume-weighted average of the tetrahedra’s centroids. Put O near the solid: far-away coordinates lose precision even in double.