function [xy,trilist,outeredge] = circle_t6_meshgen(h,curveit)
% [xy,trilist,outeredge] = circle_t6_meshgen(h,curveit)
% This function should generate a mesh of t6 elements lying inside
% the unit disk.
% The input arguments are h: a positive scalar between 0 and 1
% indicating the size
% (max edge length) of the desired elements, and curveit,
% which is either 0 or 1.  0 means that the elements lying on the
% the boundary should not be curved; 1 means that they should be curved
% so that their midpoint nodes lie on the circle itself.  Thus, the
% former is sub-parametric whereas the latter is isoparametric.
%
% The output arguments are as follows.  xy is an n-by-2 array listing
% the (x,y) coordinates of the n nodes of the mesh (one node per row).
% trilist is an m-by-6 array listing the six node of each triangle
% in the mesh (one triangle per row).  Each entry in trilist is an
% integer between 1 and n (i.e., a row index into array xy).  
%
% The six nodes of each triangle should be arranged according to the
% following pattern:
%  1
%  2 3
%  4 5 6
% and with this orientation (i.e., the ray 4->1 should be counterclockwise
% from the ray 4->6).
%
% Finally, outeredge is a p-by-2 array indicating which triangle edges are on
% the boundary (unit circle).  No triangle should have more than one
% edge on the boundary.  Each row in outeredge indicates one boundary
% edge.  The first entry in the row is a triangle number (between 1 and m,
% where m is the number of rows in trilist) and the second is an integer
% 1 to 3, where 1 indicates the edge 1-2-4, 2 indicates the edge 4-5-6,
% and 3 indicates the edge 1-3-6.
%
% Some hints how to write this routine: First select a set of well-spaced
% points of distance approximately h from each other.  This can be done
% by placing down a grid of points, deleting the points outside the unit circle
% or that fall too close to the circle, and then putting down more points
% evenly spaced on the unit circle.
% Next, you have to generate triangles for those points.  The delaunay
% function might be helpful.  Next, you have to generate all the midpoint
% nodes.  For this purpose, you need to match up every pair of triangles
% that has a common edge.  This can be done either by setting up a sparse
% matrix (using the sparse function) or by making arrays of edge endpoints
% and then sorting them (using either sort or sortrows).
% Another way to write this routine is to generate a completely regular
$ mesh of a hexagon, and then distort the mesh.