#+TITLE: Mesh.org
#+AUTHOR: screwtape
* gmsh sphere geometry utility
#+BEGIN_EXAMPLE
  I just got into 3D printing. I'm not  this doesn't  produce;   your slicer
  interested   in downloading  others'  probably can add that automatically.
  stuff   in general,  I want to  make  The  program   is basically  a  poor
  stuff. Rather than thinking  through  man's   macro  system  inside   CLOS
  extrusions, I had the idea to simply  methods   so as to use class   scope
  csg-sum   unit spheres  at different  allocated    indices  to  create   a
  scales  and  positions  (using  gmsh  combined  geometry  file output.  It
  which uses oce). The lowest surfaces  adds spheres don't worry about it.
  of a sphere generally need  supports

#+END_EXAMPLE

** Geometry of two overlapping spheres
#+name: twospheres
#+begin_src lisp :noweb yes :results output
  <<mesh>>

  (defvar *aball* (make-instance 'sphere :x 0 :y 0 :z 0 :r 10 :lc 5))
  (defvar *bball* (make-instance 'sphere :x 10 :y 0 :z 0 :r 10 :lc 5))


  (geo-beginning *aball*)
  (fill-forms *aball*)
  (fill-forms *bball*)

  (mapc 'eval (forms *aball*))
  (mapc 'eval (forms *bball*))
#+end_src

#+RESULTS: twospheres
#+begin_example
  Mesh.Algorithm = 6;

  lc = 5;
  Point(1) = {0,0,0,lc};
  Point(2) = {10,0,0,lc};
  Point(3) = {0,10,0,lc};
  Circle(1) = {2,1,3};
  Point(4) = {-10,0,0,lc};
  Point(5) = {0,-10,0,lc};
  Circle(2) = {3,1,4};
  Circle(3) = {4,1,5};
  Circle(4) = {5,1,2};
  Point(6) = {0,0,-10,lc};
  Point(7) = {0,0,10,lc};
  Circle(5) = {3,1,6};
  Circle(6) = {6,1,5};
  Circle(7) = {5,1,7};
  Circle(8) = {7,1,3};
  Circle(9) = {2,1,7};
  Circle(10) = {7,1,4};
  Circle(11) = {4,1,6};
  Circle(12) = {6,1,2};
  Curve Loop(13) = {2,8,-10};
  Surface(14) = {13};
  Curve Loop(15) = {10,3,7};
  Surface(16) = {15};
  Curve Loop(17) = {-8,-9,1};
  Surface(18) = {17};
  Curve Loop(19) = {-11,-2,5};
  Surface(20) = {19};
  Curve Loop(21) = {-5,-12,-1};
  Surface(22) = {21};
  Curve Loop(23) = {-3,11,6};
  Surface(24) = {23};
  Curve Loop(25) = {-7,4,9};
  Surface(26) = {25};
  Curve Loop(27) = {-4,12,-6};
  Surface(28) = {27};
  Surface Loop(29) = {28,26,16,14,20,24,22,18};
  Volume(30) = {29};

  Physical Surface(1) = {28,26,16,14,20,24,22,18};
  Physical Volume(2) = {30};

  Point(8) = {10,0,0,lc};
  Point(9) = {20,0,0,lc};
  Point(10) = {10,10,0,lc};
  Circle(31) = {9,8,10};
  Point(11) = {0,0,0,lc};
  Point(12) = {10,-10,0,lc};
  Circle(32) = {10,8,11};
  Circle(33) = {11,8,12};
  Circle(34) = {12,8,9};
  Point(13) = {10,0,-10,lc};
  Point(14) = {10,0,10,lc};
  Circle(35) = {10,8,13};
  Circle(36) = {13,8,12};
  Circle(37) = {12,8,14};
  Circle(38) = {14,8,10};
  Circle(39) = {9,8,14};
  Circle(40) = {14,8,11};
  Circle(41) = {11,8,13};
  Circle(42) = {13,8,9};
  Curve Loop(43) = {32,38,-40};
  Surface(44) = {43};
  Curve Loop(45) = {40,33,37};
  Surface(46) = {45};
  Curve Loop(47) = {-38,-39,31};
  Surface(48) = {47};
  Curve Loop(49) = {-41,-32,35};
  Surface(50) = {49};
  Curve Loop(51) = {-35,-42,-31};
  Surface(52) = {51};
  Curve Loop(53) = {-33,41,36};
  Surface(54) = {53};
  Curve Loop(55) = {-37,34,39};
  Surface(56) = {55};
  Curve Loop(57) = {-34,42,-36};
  Surface(58) = {57};
  Surface Loop(59) = {58,56,46,44,50,54,52,48};
  Volume(60) = {59};

  Physical Surface(3) = {58,56,46,44,50,54,52,48};
  Physical Volume(4) = {60};

#+end_example

#+BEGIN_EXAMPLE
  Output to a geo file can be opened in gmsh, and then mesh>3D'd.   Everything
  could be better!
  
#+END_EXAMPLE

** Source
#+name: mesh
#+begin_src lisp
  (defpackage mesh (:use cl) (:export sphere geo-beginning fill-forms))
  (in-package mesh)

  (defclass 3d ()
    ((points :initform nil :accessor points)
     (lc :initarg :lc :reader lc)
     (complexes :initform nil :accessor complexes)
     (phys :initform nil :accessor phys)
     (surface :accessor surface)
     (volume :accessor volume)
     (next-point :allocation :class)
     (next-phy :allocation :class)
     (next-complex :allocation :class))
    (:default-initargs :lc 5))

  (defmethod new-pnt ((obj 3d)) "
  (new-pnt 3D)
  Adds a class-scope new point.
  "
    (let ((idx (incf (slot-value obj 'next-point))))
      (if (null (points obj))
	  (setf (points obj) (list idx))
	(nconc (points obj) (list idx)))
      (values idx)))

  (defmethod pof ((obj 3d) n) "
  (pof 3d x)
  Gets the locally scoped (1- x) indexed point
  "
    (nth (1- n) (points obj)))


  (defmethod cof ((obj 3d) n) "
  (cof 3d x)
  Gets the locally scoped (1- x) indexed complex
  "
    (nth (1- n) (complexes obj)))

  (defmethod new-phy ((obj 3d)) "
  (new-phy 3d)
  Adds a class-scope physical object
  "
    (let ((idx (incf (slot-value obj 'next-phy))))
      (if (null (phys obj))
	  (setf (phys obj) (list idx))
	(nconc (phys obj) (list idx)))
      (values idx)))

  (defmethod new-cpx ((obj 3d)) "
  (new-cpx 3d)
  Adds a class-scope complex object
  "
    (let ((idx (incf (slot-value obj 'next-complex))))
      (if (null (complexes obj))
	  (setf (complexes obj) (list idx))
	(nconc (complexes obj) (list idx)))
      (values idx)))

  (defmethod shared-initialize :after ((obj 3d) names &rest args)
    (declare (ignore names args))
    (unless (slot-boundp obj 'next-point)
     (setf (slot-value obj 'next-point) 0
	   (slot-value obj 'next-complex) 0 
	   (slot-value obj 'next-phy) 0)))

  (defclass sphere (3d)
    ((x :initarg :x :accessor x)
     (y :initarg :y :accessor y)
     (z :initarg :z :accessor z)
     (r :initarg :r :accessor r)
     (sphere-forms :type list :initform (list) :accessor forms))
    (:documentation "
  (make-instance mesh:sphere :x 3 :y 4 :z 5 :r 10 :lc 3)
  Instantiates an object for the printing of a gmsh geometry of 
  said sphere equivalent to gmsh's unit sphere
  (but with class scope indices)
  "))

  (defun blank () '(format t "~%"))

  (defmethod add-point ((obj 3d) x y z)
   (let ((p (new-pnt obj)))
   `(format t "Point(~d) = {~d,~d,~d,lc};~%" ,p ,x ,y ,z)))

  (defmethod add-circle ((obj 3d) a b c &aux (p (pof obj a)) (q (pof obj b)) (r (pof obj c)))
   (let ((s (new-cpx obj)))
    `(format t "Circle(~d) = {~d,~d,~d};~%" ,s ,p ,q ,r)))

  (defun sgn (n) (if (plusp n) 1 -1))

  (defmethod add-curvel ((obj 3d) a b c &aux 
			  (p (* (sgn a) (cof obj (abs a))))
			  (q (* (sgn b) (cof obj (abs b)))) 
			  (r (* (sgn c) (cof obj (abs c)))))
   (let ((s (new-cpx obj)))
    `(format t "Curve Loop(~d) = {~d,~d,~d};~%" ,s ,p ,q ,r)))

  (defmethod add-surface ((obj 3d) x)
   (let ((s (new-cpx obj)))
    `(format t "Surface(~d) = {~d};~%" ,s ,(cof obj x))))

  (defmethod add-volume ((obj 3d) x)
   (let ((s (new-cpx obj)))
    `(format t "Volume(~d) = {~d};~%" ,s ,(cof obj x))))

  (defmethod add-surfl ((obj 3d) &rest args)
   (let ((s (new-cpx obj)) (new-args (mapcar (lambda (x) (cof obj x)) args)))
    `(format t "Surface Loop(~d) = {~@{~d~^,~}};~%" ,s ,@new-args)))

  (defmethod phy-surf ((obj 3d) &rest args)
   (let ((s (new-phy obj)) (new-args (mapcar (lambda (x) (cof obj x)) args)))
    `(format t "Physical Surface(~d) = {~@{~d~^,~}};~%" ,s ,@new-args)))

  (defmethod phy-vol ((obj 3d) vol)
   (let ((s (new-phy obj)))
    `(format t "Physical Volume(~d) = {~d};~%" ,s ,(cof obj vol))))

  (defmethod geo-beginning ((obj sphere) &aux (mesh-algo 6) (lc (lc obj)))
    (setf (forms obj)
	  (list `(format t "~a = ~d;~%" "Mesh.Algorithm" ,mesh-algo)
		(blank)
		`(format t "lc = ~d;~%" ,lc))))

  (defmethod fill-forms ((obj sphere)
			 &aux (x (x obj)) (y (y obj)) (z (z obj))
			 (r (r obj)) (lc (lc obj))) "
  (fill-forms sphere)
  Prints the Gmsh geometry of a unit sphere at (x sphere) .. 
  of radius/lc (r sphere) (lc sphere)
  Exactly equivalent to the the unit sphere from gmsh.
  The point / extrusion / physical counters are class scope
  so subsequent mesh::3d instances can be printed in the same
  geometry. 'complexes' is an erroneous name.
  "
   (setf (forms obj)
    (nconc (forms obj) 
     (list (add-point obj x y z)
      (add-point obj (+ x r) y z)
      (add-point obj x (+ y r) z)
      (add-circle obj 2 1 3)
      (add-point obj (- x r) y z)
      (add-point obj x (- y r) z)
      (add-circle obj 3 1 4) ; 10
      (add-circle obj 4 1 5)
      (add-circle obj 5 1 2)
      (add-point obj x y (- z r))
      (add-point obj x y (+ z r)) ; 14
      (add-circle obj 3 1 6)
      (add-circle obj 6 1 5)
      (add-circle obj 5 1 7)
      (add-circle obj 7 1 3)
      (add-circle obj 2 1 7)
      (add-circle obj 7 1 4)
      (add-circle obj 4 1 6)
      (add-circle obj 6 1 2) ; 22
      (add-curvel obj 2 8 -10)
      (add-surface obj 13)
      (add-curvel obj 10 3 7)
      (add-surface obj 15)
      (add-curvel obj -8 -9 1)
      (add-surface obj 17) ; 28
      (add-curvel obj -11 -2 5)
      (add-surface obj 19)
      (add-curvel obj -5 -12 -1)
      (add-surface obj 21)
      (add-curvel obj -3 11 6)
      (add-surface obj 23)
      (add-curvel obj -7 4 9)
      (add-surface obj 25)
      (add-curvel obj -4 12 -6)
      (add-surface obj 27) ; 38
      (add-surfl obj 28 26 16 14 20 24 22 18)
      (add-volume obj 29)
      (blank)
      (phy-surf obj 28 26 16 14 20 24 22 18)
      (phy-vol obj 30)
      (blank)))))

#+end_src