Index: /issm/trunk-jpl/src/jl/solve/analyses/stressbalanceanalysis.jl
===================================================================
--- /issm/trunk-jpl/src/jl/solve/analyses/stressbalanceanalysis.jl	(revision 26657)
+++ /issm/trunk-jpl/src/jl/solve/analyses/stressbalanceanalysis.jl	(revision 26658)
@@ -66,4 +66,11 @@
 
 end #}}}
+function CreateKMatrix(analysis::StressbalanceAnalysis,element::Tria)# {{{
+
+	#Initialize Element matrix
+	Ke = ElementMatrix(element.nodes)
+
+	return Ke
+end #}}}
 function GetSolutionFromInputs(analysis::StressbalanceAnalysis,ug::Vector{Float64},element::Tria) #{{{
 
Index: /issm/trunk-jpl/src/jl/solve/elementmatrix.jl
===================================================================
--- /issm/trunk-jpl/src/jl/solve/elementmatrix.jl	(revision 26657)
+++ /issm/trunk-jpl/src/jl/solve/elementmatrix.jl	(revision 26658)
@@ -7,4 +7,55 @@
 end #}}}
 function ElementMatrix(nodes::Vector{Node})#{{{
-	error("not implemented yet")
+
+	#Get matrix size
+	nrows = NumberOfDofs(nodes,GsetEnum)
+
+	#Initialize element matrix with zeros
+	values = zeros(nrows,nrows)
+
+	#Get dof lists
+	gglobaldoflist=GetGlobalDofList(nodes,nrows,GsetEnum)
+	fglobaldoflist=GetGlobalDofList(nodes,nrows,FsetEnum)
+	sglobaldoflist=GetGlobalDofList(nodes,nrows,SsetEnum)
+
+	return ElementMatrix(nrows,gglobaldoflist,fglobaldoflist,sglobaldoflist,values)
 end#}}}
+function Base.show(io::IO, this::ElementMatrix)# {{{
+
+	println(io,"ElementMatrix:")
+	println(io,"   nrows: ",this.nrows)
+	println(io,"   gglobaldoflist: ",this.gglobaldoflist)
+	println(io,"   fglobaldoflist: ",this.fglobaldoflist)
+	println(io,"   sglobaldoflist: ",this.sglobaldoflist)
+	print(io,"   values: ")
+	display(this.values)
+end# }}}
+
+mutable struct ElementVector#{{{
+	nrows::Int64
+	fglobaldoflist::Vector{Int64}
+	values::Vector{Float64}
+end #}}}
+function ElementVector(nodes::Vector{Node})#{{{
+
+	#Get matrix size
+	nrows = NumberOfDofs(nodes,GsetEnum)
+
+	#Initialize element matrix with zeros
+	values = zeros(nrows,nrows)
+
+	#Get dof list
+	fglobaldoflist=GetGlobalDofList(nodes,nrows,FsetEnum)
+
+	return ElementVector(nrows,fglobaldoflist,values)
+end#}}}
+function Base.show(io::IO, this::ElementVector)# {{{
+
+	println(io,"ElementVector:")
+	println(io,"   nrows: ",this.nrows)
+	println(io,"   gglobaldoflist: ",this.gglobaldoflist)
+	println(io,"   fglobaldoflist: ",this.fglobaldoflist)
+	println(io,"   sglobaldoflist: ",this.sglobaldoflist)
+	print(io,"   values: ")
+	display(this.values)
+end# }}}
Index: /issm/trunk-jpl/src/jl/solve/nodes.jl
===================================================================
--- /issm/trunk-jpl/src/jl/solve/nodes.jl	(revision 26657)
+++ /issm/trunk-jpl/src/jl/solve/nodes.jl	(revision 26658)
@@ -95,15 +95,15 @@
 	elseif setenum==FsetEnum
 		for i=1:node.gsize
-			if  node.fdoflist[i]!=-1
+			#if  node.fdoflist[i]!=-1
 				count += 1
 				doflist[count] = node.fdoflist[i]
-			end
+			#end
 		end
 	elseif setenum==SsetEnum
 		for i=1:node.gsize
-			if  node.sdoflist[i]!=-1
+			#if  node.sdoflist[i]!=-1
 				count += 1
 				doflist[count] = node.sdoflist[i]
-			end
+			#end
 		end
 	else
@@ -112,4 +112,20 @@
 
 	return count
+
+end# }}}
+function GetGlobalDofList(nodes::Vector{Node},ndofs::Int64,setenum::IssmEnum) #{{{
+
+	#Allocate list
+	doflist = Vector{Int64}(undef,ndofs)
+
+	#Assign values 
+	count = 0
+	for i in 1:length(nodes)
+		count = GetDofList(nodes[i],doflist,count,setenum)
+	end
+	println(count," ",ndofs)
+	@assert count==ndofs
+
+	return doflist
 
 end# }}}