Index: /issm/trunk-jpl/src/m/model/mesh/meshnodensity.m
===================================================================
--- /issm/trunk-jpl/src/m/model/mesh/meshnodensity.m	(revision 11010)
+++ /issm/trunk-jpl/src/m/model/mesh/meshnodensity.m	(revision 11011)
@@ -9,10 +9,10 @@
 %
 %   Usage:
-%      md=setmesh(md,domainname)
-%   or md=setmesh(md,domainname,riftname);
+%      md=meshnodensity(md,domainname)
+%   or md=meshnodensity(md,domainname,riftname);
 %
 %   Examples:
-%      md=setmesh(md,'DomainOutline.exp');
-%      md=setmesh(md,'DomainOutline.exp','Rifts.exp');
+%      md=meshnodensity(md,'DomainOutline.exp');
+%      md=meshnodensity(md,'DomainOutline.exp','Rifts.exp');
 
 if (nargin==2),
Index: sm/trunk-jpl/src/m/model/mesh/setmesh.m
===================================================================
--- /issm/trunk-jpl/src/m/model/mesh/setmesh.m	(revision 11010)
+++ 	(revision )
@@ -1,85 +1,0 @@
-function md=setmesh(md,domainname,varargin)
-%MESH - create model mesh
-%
-%   This routine creates a model mesh using TriMesh and a domain outline, to within a certain resolution
-%   where md is a @model object, domainname is the name of an Argus domain outline file, 
-%   and resolution is a characteristic length for the mesh (same unit as the domain outline
-%   unit). Riftname is an optional argument (Argus domain outline) describing rifts.
-%
-%   Usage:
-%      md=setmesh(md,domainname,resolution)
-%   or md=setmesh(md,domainname,riftname, resolution)
-%
-%   Examples:
-%      md=setmesh(md,'DomainOutline.exp',1000);
-%      md=setmesh(md,'DomainOutline.exp','Rifts.exp',1500);
-
-%Figure out a characteristic area. Resolution is a node oriented concept (ex a 1000m  resolution node would 
-%be made of 1000*1000 area squares). 
-if (nargin==3),
-	resolution=varargin{1};
-	riftname='';
-end
-if (nargin==4),
-	riftname=varargin{1};
-	resolution=varargin{2};
-end
-
-%Check that mesh was not already run, and warn user: 
-if md.mesh.numberofelements~=0,
-	choice=input('This model already has a mesh. Are you sure you want to go ahead? (y/n)','s');
-	if ~strcmp(choice,'y')
-		disp('no meshing done ... exiting');
-		return
-	end
-end
-
-area=resolution^2;
-
-%Mesh using TriMesh
-if strcmp(riftname,''),
-	[md.mesh.elements,md.mesh.x,md.mesh.y,md.mesh.segments,md.mesh.segmentmarkers]=TriMesh(domainname,area,'yes');
-else
-	[elements,x,y,segments,segmentmarkers]=TriMesh(domainname,riftname,area,'yes');
-
-	%check that all the created nodes belong to at least one element
-	orphan=find(~ismember([1:length(x)],sort(unique(elements(:)))));
-	for i=1:length(orphan),
-		%get rid of the orphan node i
-		%update x and y
-		x=[x(1:orphan(i)-(i-1)-1); x(orphan(i)-(i-1)+1:end)];
-		y=[y(1:orphan(i)-(i-1)-1); y(orphan(i)-(i-1)+1:end)];
-		%update elements
-		pos=find(elements>orphan(i)-(i-1));
-		elements(pos)=elements(pos)-1;
-		%update segments
-		pos1=find(segments(:,1)>orphan(i)-(i-1));
-		pos2=find(segments(:,2)>orphan(i)-(i-1));
-		segments(pos1,1)=segments(pos1,1)-1;
-		segments(pos2,2)=segments(pos2,2)-1;
-	end
-
-	%plug into md
-	md.mesh.x=x;
-	md.mesh.y=y;
-	md.mesh.elements=elements;
-	md.mesh.segments=segments;
-	md.mesh.segmentmarkers=segmentmarkers;
-end
-
-%Fill in rest of fields:
-md.mesh.numberofelements=length(md.mesh.elements);
-md.mesh.numberofvertices=length(md.mesh.x);
-md.mesh.z=zeros(md.mesh.numberofvertices,1);
-md.mesh.vertexonboundary=zeros(md.mesh.numberofvertices,1); md.mesh.vertexonboundary(md.mesh.segments(:,1:2))=1;
-md.mesh.vertexonbed=ones(md.mesh.numberofvertices,1);
-md.mesh.vertexonsurface=ones(md.mesh.numberofvertices,1);
-md.mesh.elementonbed=ones(md.mesh.numberofelements,1);
-md.mesh.elementonsurface=ones(md.mesh.numberofelements,1);
-
-%Now, build the connectivity tables for this mesh.
-md.mesh.vertexconnectivity=NodeConnectivity(md.mesh.elements,md.mesh.numberofvertices);
-md.mesh.elementconnectivity=ElementConnectivity(md.mesh.elements,md.mesh.vertexconnectivity);
-
-%type of model
-md.mesh.dimension=2;
Index: /issm/trunk-jpl/src/m/model/mesh/triangle.m
===================================================================
--- /issm/trunk-jpl/src/m/model/mesh/triangle.m	(revision 11011)
+++ /issm/trunk-jpl/src/m/model/mesh/triangle.m	(revision 11011)
@@ -0,0 +1,85 @@
+function md=triangle(md,domainname,varargin)
+%TRIANGLE - create model mesh using the triangle package
+%
+%   This routine creates a model mesh using TriMesh and a domain outline, to within a certain resolution
+%   where md is a @model object, domainname is the name of an Argus domain outline file, 
+%   and resolution is a characteristic length for the mesh (same unit as the domain outline
+%   unit). Riftname is an optional argument (Argus domain outline) describing rifts.
+%
+%   Usage:
+%      md=triangle(md,domainname,resolution)
+%   or md=triangle(md,domainname,riftname, resolution)
+%
+%   Examples:
+%      md=triangle(md,'DomainOutline.exp',1000);
+%      md=triangle(md,'DomainOutline.exp','Rifts.exp',1500);
+
+%Figure out a characteristic area. Resolution is a node oriented concept (ex a 1000m  resolution node would 
+%be made of 1000*1000 area squares). 
+if (nargin==3),
+	resolution=varargin{1};
+	riftname='';
+end
+if (nargin==4),
+	riftname=varargin{1};
+	resolution=varargin{2};
+end
+
+%Check that mesh was not already run, and warn user: 
+if md.mesh.numberofelements~=0,
+	choice=input('This model already has a mesh. Are you sure you want to go ahead? (y/n)','s');
+	if ~strcmp(choice,'y')
+		disp('no meshing done ... exiting');
+		return
+	end
+end
+
+area=resolution^2;
+
+%Mesh using TriMesh
+if strcmp(riftname,''),
+	[md.mesh.elements,md.mesh.x,md.mesh.y,md.mesh.segments,md.mesh.segmentmarkers]=TriMesh(domainname,area,'yes');
+else
+	[elements,x,y,segments,segmentmarkers]=TriMesh(domainname,riftname,area,'yes');
+
+	%check that all the created nodes belong to at least one element
+	orphan=find(~ismember([1:length(x)],sort(unique(elements(:)))));
+	for i=1:length(orphan),
+		%get rid of the orphan node i
+		%update x and y
+		x=[x(1:orphan(i)-(i-1)-1); x(orphan(i)-(i-1)+1:end)];
+		y=[y(1:orphan(i)-(i-1)-1); y(orphan(i)-(i-1)+1:end)];
+		%update elements
+		pos=find(elements>orphan(i)-(i-1));
+		elements(pos)=elements(pos)-1;
+		%update segments
+		pos1=find(segments(:,1)>orphan(i)-(i-1));
+		pos2=find(segments(:,2)>orphan(i)-(i-1));
+		segments(pos1,1)=segments(pos1,1)-1;
+		segments(pos2,2)=segments(pos2,2)-1;
+	end
+
+	%plug into md
+	md.mesh.x=x;
+	md.mesh.y=y;
+	md.mesh.elements=elements;
+	md.mesh.segments=segments;
+	md.mesh.segmentmarkers=segmentmarkers;
+end
+
+%Fill in rest of fields:
+md.mesh.numberofelements=length(md.mesh.elements);
+md.mesh.numberofvertices=length(md.mesh.x);
+md.mesh.z=zeros(md.mesh.numberofvertices,1);
+md.mesh.vertexonboundary=zeros(md.mesh.numberofvertices,1); md.mesh.vertexonboundary(md.mesh.segments(:,1:2))=1;
+md.mesh.vertexonbed=ones(md.mesh.numberofvertices,1);
+md.mesh.vertexonsurface=ones(md.mesh.numberofvertices,1);
+md.mesh.elementonbed=ones(md.mesh.numberofelements,1);
+md.mesh.elementonsurface=ones(md.mesh.numberofelements,1);
+
+%Now, build the connectivity tables for this mesh.
+md.mesh.vertexconnectivity=NodeConnectivity(md.mesh.elements,md.mesh.numberofvertices);
+md.mesh.elementconnectivity=ElementConnectivity(md.mesh.elements,md.mesh.vertexconnectivity);
+
+%type of model
+md.mesh.dimension=2;
Index: /issm/trunk-jpl/src/m/utils/Mesh/roundmesh.m
===================================================================
--- /issm/trunk-jpl/src/m/utils/Mesh/roundmesh.m	(revision 11010)
+++ /issm/trunk-jpl/src/m/utils/Mesh/roundmesh.m	(revision 11011)
@@ -23,5 +23,5 @@
 
 %Call Bamg
-md=setmesh(md,'RoundDomainOutline.exp',resolution);
+md=triangle(md,'RoundDomainOutline.exp',resolution);
 %md=bamg(md,'domain','RoundDomainOutline.exp','hmin',resolution);