(define coarse-square-triangulation
;; This defines a triangulation of the unit square with two triangles,
;; five edges, and four vertices.
(let ((pll (make-Point 0.0 0.0)) ; lower-left
(plr (make-Point 1.0 0.0)) ; lower-right
(pul (make-Point 0.0 1.0)) ; upper-left
(pur (make-Point 1.0 1.0))); upper-right
(with-co-instantiation
((vll Boundary-vertex ; lower-left
point: pll
edges: (Vertex-edges<-list (list bottom diagonal left))
boundary-edges: (Boundary-vertex-boundary-edges<-list (list bottom left))
)
(vlr Boundary-vertex ; lower-right
point: plr
edges: (Vertex-edges<-list (list right bottom))
boundary-edges: (Boundary-vertex-boundary-edges<-list (list right bottom))
)
(vul Boundary-vertex ; upper-left
point: pul
edges: (Vertex-edges<-list (list left top))
boundary-edges: (Boundary-vertex-boundary-edges<-list (list left top))
)
(vur Boundary-vertex ; upper-right
point: pur
edges: (Vertex-edges<-list (list top diagonal right))
boundary-edges: (Boundary-vertex-boundary-edges<-list (list top right))
)
(bottom Boundary-edge
vertex-1: vll vertex-2: vlr
)
(top Boundary-edge
vertex-1: vul vertex-2: vur
)
(left Boundary-edge
vertex-1: vll vertex-2: vul
)
(right Boundary-edge
vertex-1: vlr vertex-2: vur
)
(diagonal Edge
vertex-1: vll vertex-2: vur
)
(tul Triangle ; upper-left
vertex-1: vll vertex-2: vur vertex-3: vul
edge-1: diagonal edge-2: top edge-3: left)
(tlr Triangle ; lower-right
vertex-1: vll vertex-2: vlr vertex-3: vur
edge-1: bottom edge-2: right edge-3: diagonal))
(instantiate Triangulation
vertices: (Triangulation-vertices<-list
(list vll vlr vul vur))
edges: (Triangulation-edges<-list
(list bottom top left right diagonal))
triangles: (Triangulation-triangles<-list
(list tul tlr))
has-polygonal-boundary?: #t))))
(define-generic (refine (o)))
;; We have temporarily replaced with-co-instantiation with instantiate
;; and explicitly initialized all fields so that fill-other-fields!
;; is not called. This doubles (!) the speed of (refine triangulation)
;; even with the space-filling-curve ordering of all the vertices
;; edges and triangles. I hope Christian can change the macro expansion
;; to be smarter so that we can use with-co-instantiation again.
;; Each old vertex introduces a new vertex at the same location.
;; We'll store this new vertex in Vertex-fine of each old vertex.
(define-handy-method (refine (v Vertex))
(set! fine (instantiate Vertex
index: #f
fine: '()
edges: '()
point: point)))
(define-handy-method (refine (v Boundary-vertex))
(set! fine (instantiate Boundary-vertex
index: #f
fine: '()
edges: '()
boundary-index: #f
boundary-edges: '()
point: point))) ; in the refine methods only
;; Each old edge has a new vertex in the middle, and
;; introduces two new edges. We'll store a list of the new vertex
;; and the two new edges in Edge-fine
(define-handy-method (refine (e Edge))
(let ((vertex-1 (Vertex-fine vertex-1))
(vertex-2 (Vertex-fine vertex-2)))
(let ((point-1 (Vertex-point vertex-1))
(point-2 (Vertex-point vertex-2)))
(let ((middle (instantiate Vertex
index: #f
fine: '()
edges: '()
point: (Point-scale 0.5
(Point-add point-1
point-2)))))
(let ((left (instantiate Edge
index: #f
fine: '()
vertex-1: vertex-1 ; existing vertex is vertex-1
vertex-2: middle))
(right (instantiate Edge
index: #f
fine: '()
vertex-1: vertex-2 ; existing vertex is vertex-1
vertex-2: middle)))
(Edge-fine-set! e (list middle left right))
;; add edge information to new vertices
(for-each
(lambda (e)
(Edge-for-each-vertex
(lambda (v)
(with-access
v (Vertex edges)
(set! edges (cons e edges))))
e))
(cdr (Edge-fine e))))))))
(define-handy-method (refine (e Boundary-edge))
(let ((vertex-1 (Vertex-fine vertex-1))
(vertex-2 (Vertex-fine vertex-2)))
(let ((point-1 (Vertex-point vertex-1))
(point-2 (Vertex-point vertex-2)))
(let ((middle (instantiate Boundary-vertex
index: #f
fine: '()
edges: '()
boundary-index: #f
boundary-edges: '()
point: (Point-scale 0.5
(Point-add point-1
point-2)))))
(let ((left (instantiate Boundary-edge
index: #f
fine: '()
boundary-index: #f
vertex-1: vertex-1 ; existing vertex is vertex-1
vertex-2: middle))
(right (instantiate Boundary-edge
index: #f
fine: '()
boundary-index: #f
vertex-1: vertex-2 ; existing vertex is vertex-1
vertex-2: middle)))
(Edge-fine-set! e (list middle left right))
;; add edge information to new vertices
(for-each
(lambda (e)
(Edge-for-each-vertex
(lambda (v)
(with-access
v (Boundary-vertex edges boundary-edges)
(set! edges (cons e edges))
(set! boundary-edges (cons e boundary-edges))))
e))
(cdr (Edge-fine e))))))))
;; each old triangle is decomposed into four new triangles and
;; generates three new edges
(define-handy-method (refine (t Triangle))
(let ((old-vertex-1 vertex-1)
(old-vertex-2 vertex-2)
(old-vertex-3 vertex-3)
(old-edge-1 edge-1)
(old-edge-2 edge-2)
(old-edge-3 edge-3))
(let ((vertex-1 (Vertex-fine old-vertex-1))
(vertex-2 (Vertex-fine old-vertex-2))
(vertex-3 (Vertex-fine old-vertex-3))
(mid-vertex-1 (car (Edge-fine old-edge-1)))
(mid-vertex-2 (car (Edge-fine old-edge-2)))
(mid-vertex-3 (car (Edge-fine old-edge-3))))
(let ((mid1 (instantiate Edge
index: #f
fine: '()
vertex-1: mid-vertex-1
vertex-2: mid-vertex-2))
(mid2 (instantiate Edge
index: #f
fine: '()
vertex-1: mid-vertex-2
vertex-2: mid-vertex-3))
(mid3 (instantiate Edge
index: #f
fine: '()
vertex-1: mid-vertex-3
vertex-2: mid-vertex-1)))
(let ((t1 (instantiate Triangle
index: #f
fine: '()
vertex-1: vertex-1
vertex-2: mid-vertex-1
vertex-3: mid-vertex-3
edge-1: (if (eq? (Edge-vertex-1 old-edge-1) old-vertex-1)
(cadr (Edge-fine old-edge-1))
(caddr (Edge-fine old-edge-1)))
edge-2: mid3
edge-3: (if (eq? (Edge-vertex-1 old-edge-3) old-vertex-1)
(cadr (Edge-fine old-edge-3))
(caddr (Edge-fine old-edge-3)))))
(t2 (instantiate Triangle
index: #f
fine: '()
vertex-1: vertex-2
vertex-2: mid-vertex-2
vertex-3: mid-vertex-1
edge-1: (if (eq? (Edge-vertex-1 old-edge-2) old-vertex-2)
(cadr (Edge-fine old-edge-2))
(caddr (Edge-fine old-edge-2)))
edge-2: mid1
edge-3: (if (eq? (Edge-vertex-1 old-edge-1) old-vertex-2)
(cadr (Edge-fine old-edge-1))
(caddr (Edge-fine old-edge-1)))))
(t3 (instantiate Triangle
index: #f
fine: '()
vertex-1: vertex-3
vertex-2: mid-vertex-3
vertex-3: mid-vertex-2
edge-1: (if (eq? (Edge-vertex-1 old-edge-3) old-vertex-3)
(cadr (Edge-fine old-edge-3))
(caddr (Edge-fine old-edge-3)))
edge-2: mid2
edge-3: (if (eq? (Edge-vertex-1 old-edge-2) old-vertex-3)
(cadr (Edge-fine old-edge-2))
(caddr (Edge-fine old-edge-2)))
))
(t4 (instantiate Triangle
index: #f
fine: '()
vertex-1: mid-vertex-1
vertex-2: mid-vertex-2
vertex-3: mid-vertex-3
edge-1: mid1
edge-2: mid2
edge-3: mid3)))
(let ((new-triangles (list t1 t2 t3 t4))
(new-edges (list mid1 mid2 mid3)))
(Triangle-fine-set! t (list new-triangles new-edges))
;; add new edge information to points; none of the new
;;edges are on the boundary
(for-each (lambda (e)
(Edge-for-each-vertex
(lambda (v)
(with-access
v (Vertex edges)
;; accumulate in a list
(set! edges (cons e edges))))
e))
new-edges)))))))
(define-method (refine (t Triangulation))
(let ((old-vertices (Triangulation-vertices->list
(Triangulation-vertices t)))
(old-edges (Triangulation-edges->list
(Triangulation-edges t)))
(old-triangles (Triangulation-triangles->list
(Triangulation-triangles t))))
(for-each refine old-vertices)
(for-each refine old-edges)
(for-each refine old-triangles)
(let ((new-Triangulation
(instantiate Triangulation
triangles: (Triangulation-triangles<-list
(mappend1 (lambda (t)
(car (Triangle-fine t)))
old-triangles))
edges: (Triangulation-edges<-list
(append (mappend1 (lambda (e)
(cdr (Edge-fine e)))
old-edges)
(mappend1 (lambda (t)
(cadr (Triangle-fine t)))
old-triangles)))
vertices: (Triangulation-vertices<-list
(append (map1 Vertex-fine
old-vertices)
(map1 (lambda (e)
(car (Edge-fine e)))
old-edges)))
has-polygonal-boundary?: (Triangulation-has-polygonal-boundary? t)
;; if the triangulation has a polygonal boundary,
;; then the refined triangulation has the same limits
;; as the original triangulation
limits: (and (Triangulation-has-polygonal-boundary? t)
(Triangulation-limits t)))))
;; convert temporary lists to internal data structures
(Triangulation-for-each-vertex
(lambda (v)
(with-access
v (Vertex edges)
(set! edges (Vertex-edges<-list edges))))
new-Triangulation)
(Triangulation-for-each-boundary-vertex
(lambda (v)
(with-access
v (Boundary-vertex boundary-edges)
(set! boundary-edges (Boundary-vertex-boundary-edges<-list boundary-edges))))
new-Triangulation)
;; convert fine fields on coarse triangulation
;; to what we need outside refine.
;; Vertex-fine is OK.
(Triangulation-for-each-edge
(lambda (e)
(with-access
e (Edge fine)
(set! fine (car fine)))) ; keep only the new vertex
t)
(Triangulation-for-each-triangle
(lambda (t)
(Triangle-fine-set! t '())) ; I don't believe we need anything here
t)
new-Triangulation)))
(define (Triangulation-consistent? t)
(Triangulation-for-each-triangle
(lambda (t)
(with-access
t (Triangle vertex-1 vertex-2 vertex-3 edge-1 edge-2 edge-3)
(for-each
(lambda (edge vertex-1 vertex-2)
(if (not (and (Edge-contains-vertex? edge vertex-1)
(Edge-contains-vertex? edge vertex-2)))
(for-each
display
(list "edge-k is not between vertex-k and vertex-k+1: ("
(Point-x (Vertex-point (Edge-vertex-1 edge))) ","
(Point-y (Vertex-point (Edge-vertex-1 edge))) ") ("
(Point-x (Vertex-point (Edge-vertex-2 edge))) ","
(Point-y (Vertex-point (Edge-vertex-2 edge))) ")\n"))))
(list edge-1 edge-2 edge-3)
(list vertex-1 vertex-2 vertex-3)
(list vertex-2 vertex-3 vertex-1))))
t)
(Triangulation-for-each-edge
(lambda (e)
(Edge-for-each-vertex
(lambda (v)
(if (not (Vertex-contains-edge? v e))
(for-each display (list "edge is not in vertex-edges: ("
(Point-x (Vertex-point v)) ","
(Point-y (Vertex-point v)) ")\n"))))
e))
t)
(Triangulation-for-each-boundary-edge
(lambda (e)
(Edge-for-each-vertex
(lambda (v)
(if (not (Boundary-vertex? v))
(for-each display (list "Boundary edge has a non-boundary vertex: "
(Point-x (Vertex-point v)) ","
(Point-y (Vertex-point v)) ")\n"))
(if (not (Boundary-vertex-contains-boundary-edge? v e))
(for-each display (list "Boundary edge is not in Vertex-boundary-edges: "
(Point-x (Vertex-point v)) ","
(Point-y (Vertex-point v)) ")\n")))))
e))
t)
(Triangulation-for-each-vertex
(lambda (v)
(Vertex-for-each-edge
(lambda (e)
(if (not (Edge-contains-vertex? e v))
(for-each
display
(list
"edge is in vertex-edges, but does not contain this vertex: ("
(Point-x (Vertex-point v)) ","
(Point-y (Vertex-point v)) ")\n"))))
v))
t)
(Triangulation-for-each-boundary-vertex
(lambda (v)
(Boundary-vertex-for-each-boundary-edge
(lambda (e)
(if (not (Boundary-edge? e))
(for-each display (list "Boundary vertex has a non-boundary edge in boundary-edges: "
(Point-x (Vertex-point v)) ","
(Point-y (Vertex-point v)) ")\n"))))
v))
t))