From efbc2b1043536e11bd0a39a6873cb2836744e6da Mon Sep 17 00:00:00 2001
From: steve <kd2cca@gmail.com>
Date: Wed, 24 Jul 2019 23:14:26 -0400
Subject: [PATCH] wip try rfunctions with dual contours

---
 examples/fea.jl |  3 ++-
 src/frep.jl     | 22 ++++++++++++++--------
 2 files changed, 16 insertions(+), 9 deletions(-)

diff --git a/examples/fea.jl b/examples/fea.jl
index f467597..34e1879 100644
--- a/examples/fea.jl
+++ b/examples/fea.jl
@@ -1,4 +1,5 @@
 using Descartes
+using Meshing
 
 # params
 beam_size = [50,10,10]
@@ -17,4 +18,4 @@ for i = 1:hole_ct
     global c = diff(c, h)
 end
 
-save("fea.ply",HomogenousMesh(c))
+save("fea.ply",HomogenousMesh(c, samples=(50,50,50), algorithm=DualContours()))
diff --git a/src/frep.jl b/src/frep.jl
index 9a4a614..88acf70 100644
--- a/src/frep.jl
+++ b/src/frep.jl
@@ -1,8 +1,14 @@
 # http://en.wikipedia.org/wiki/Function_representation
 
+# rvechev functions for smooth max and min
+rmax(x,y) = (x+y+sqrt((x-y)^2))/2
+rmax(x, y, z...) = rmax(x, rmax(y,z...))
+rmin(x,y) = (x+y-sqrt((x-y)^2))/2
+rmin(x, y, z...) = rmin(x, rmin(y,z...))
+
 function _radius(a,b,r)
     if abs(a-b) >= r
-        return min(a,b)
+        return rmin(a,b)
     else
         return b+r*sin(pi/4+asin((a-b)/(r*sqrt(2))))-r
     end
@@ -21,7 +27,7 @@ function FRep(p::Cylinder, v)
     x = v[1]*it[1,1]+v[2]*it[1,2]+v[3]*it[1,3]+it[1,4]
     y = v[1]*it[2,1]+v[2]*it[2,2]+v[3]*it[2,3]+it[2,4]
     z = v[1]*it[3,1]+v[2]*it[3,2]+v[3]*it[3,3]+it[3,4]
-    max(-z+p.bottom, z-p.height-p.bottom, sqrt(x*x + y*y) - p.radius)
+    rmax(-z+p.bottom, z-p.height-p.bottom, sqrt(x*x + y*y) - p.radius)
 end
 
 function FRep(p::Cuboid, v)
@@ -31,26 +37,26 @@ function FRep(p::Cuboid, v)
     z = v[1]*it[3,1]+v[2]*it[3,2]+v[3]*it[3,3]+it[3,4]
     dx, dy, dz = p.dimensions
     lbx, lby,lbz = p.lowercorner
-    max(-x+lbx, x-dx-lbx,
+    rmax(-x+lbx, x-dx-lbx,
         -y+lby, y-dy-lby,
         -z+lbz, z-dz-lbz)
 end
 
 function FRep(u::CSGUnion, v)
-    min(FRep(u.left, v),FRep(u.right, v))
+    rmin(FRep(u.left, v),FRep(u.right, v))
 end
 
 function FRep(u::CSGDiff, v)
-    max(FRep(u.left, v), -FRep(u.right, v))
+    rmax(FRep(u.left, v), -FRep(u.right, v))
 end
 
 function FRep(u::CSGIntersect, v)
-    max(FRep(u.left, v), FRep(u.right, v))
+    rmax(FRep(u.left, v), FRep(u.right, v))
 end
 
 function FRep(s::Shell, v)
     r = FRep(s.primitive, v)
-    max(r, -r-s.distance)
+    rmax(r, -r-s.distance)
 end
 
 function FRep(u::RadiusedCSGUnion, v)
@@ -85,5 +91,5 @@ end
 
 function FRep(p::LinearExtrude, v)
     r = FRep(s.primitive, v)
-    max(max(-z,z-p.height), r)
+    rmax(rmax(-z,z-p.height), r)
 end