Index: /issm/trunk-jpl/src/jl/core/analyses/masstransportanalysis.jl
===================================================================
--- /issm/trunk-jpl/src/jl/core/analyses/masstransportanalysis.jl	(revision 26737)
+++ /issm/trunk-jpl/src/jl/core/analyses/masstransportanalysis.jl	(revision 26738)
@@ -22,5 +22,5 @@
 	numdof = 1
 	for i in 1:md.mesh.numberofvertices
-		push!(nodes,Node(i,i,true,numdof,-ones(Int64,numdof), ones(Int64,numdof), -ones(Int64,numdof), zeros(numdof)))
+		push!(nodes,Node(i,i,true,true,numdof,-ones(Int64,numdof), ones(Int64,numdof), -ones(Int64,numdof), zeros(numdof)))
 	end
 end#}}}
@@ -70,4 +70,7 @@
 end #}}}
 function CreateKMatrix(analysis::MasstransportAnalysis,element::Tria)# {{{
+
+	#Return if there is no ice in this element
+	if(!IsIceInElement(element)) return end
 
 	#Internmediaries
@@ -136,4 +139,7 @@
 end #}}}
 function CreatePVector(analysis::MasstransportAnalysis,element::Tria)# {{{
+
+	#Return if there is no ice in this element
+	if(!IsIceInElement(element)) return end
 
 	#Internmediaries
@@ -265,2 +271,5 @@
 	AddInput(element, BaseEnum,    newbase,    P1Enum)
 end#}}}
+function UpdateConstraints(analysis::MasstransportAnalysis, femmodel::FemModel) #{{{
+	SetActiveNodesLSMx(femmodel)
+end#}}}
Index: /issm/trunk-jpl/src/jl/core/analyses/stressbalanceanalysis.jl
===================================================================
--- /issm/trunk-jpl/src/jl/core/analyses/stressbalanceanalysis.jl	(revision 26737)
+++ /issm/trunk-jpl/src/jl/core/analyses/stressbalanceanalysis.jl	(revision 26738)
@@ -27,5 +27,5 @@
 	numdof = 2
 	for i in 1:md.mesh.numberofvertices
-		push!(nodes,Node(i,i,true,numdof,-ones(Int64,numdof), ones(Int64,numdof), -ones(Int64,numdof), zeros(numdof)))
+		push!(nodes,Node(i,i,true,true,numdof,-ones(Int64,numdof), ones(Int64,numdof), -ones(Int64,numdof), zeros(numdof)))
 	end
 end#}}}
@@ -91,4 +91,6 @@
 function CreateKMatrix(analysis::StressbalanceAnalysis,element::Tria)# {{{
 
+	if(!IsIceInElement(element)) return end
+
 	#Internmediaries
 	numnodes = 3
@@ -153,4 +155,6 @@
 end #}}}
 function CreatePVector(analysis::StressbalanceAnalysis,element::Tria)# {{{
+
+	if(!IsIceInElement(element)) return end
 
 	#Internmediaries
@@ -271,2 +275,5 @@
 	AddInput(element, VelEnum, vel, P1Enum)
 end#}}}
+function UpdateConstraints(analysis::StressbalanceAnalysis, femmodel::FemModel) #{{{
+	SetActiveNodesLSMx(femmodel)
+end#}}}
Index: /issm/trunk-jpl/src/jl/core/elements.jl
===================================================================
--- /issm/trunk-jpl/src/jl/core/elements.jl	(revision 26737)
+++ /issm/trunk-jpl/src/jl/core/elements.jl	(revision 26738)
@@ -46,4 +46,19 @@
 	end
 end #}}}
+function InputUpdateFromVector(element::Tria, vector::Vector{Float64}, enum::IssmEnum, layout::IssmEnum) #{{{
+
+	lidlist = element.vertexids
+	data = Vector{Float64}(undef, 3)
+
+	if(layout==VertexSIdEnum)
+		for i in 1:3
+			data[i] = vector[element.vertices[i].sid]
+			@assert isfinite(data[i])
+		end
+		SetTriaInput(element.inputs, enum, P1Enum, lidlist, data)
+	else
+		error("layout ", layout, " not supported yet");
+	end
+end #}}}
 function Update(element::Tria, inputs::Inputs, index::Int64, md::model, finiteelement::IssmEnum) #{{{
 
@@ -108,4 +123,16 @@
 	InputServe!(element,input)
 	return input
+
+end # }}}
+function GetInputListOnNodes!(element::Tria, vector::Vector{Float64}, enum::IssmEnum) # {{{
+
+	#Get Input first 
+	input = GetInput(element, enum)
+
+	#Get value at each vertex (i.e. P1 Nodes)
+	gauss=GaussTria(P1Enum)
+	for i in 1:gauss.numgauss
+		vector[i] = GetInputValue(input, gauss, i)
+	end
 
 end # }}}
@@ -197,4 +224,16 @@
 	end
 end#}}}
+function IsIceInElement(element::Tria) #{{{
+	#We consider that an element has ice if at least one of its nodes has a negative level set
+
+	input=GetInput(element, MaskIceLevelsetEnum)
+
+	if GetInputMin(input)<0
+		return true
+	else
+		return false
+	end
+
+end#}}}
 function GetIcefrontCoordinates!(element::Tria, xyz_front::Matrix{Float64}, xyz_list::Matrix{Float64}, levelsetenum::IssmEnum) #{{{
 
Index: /issm/trunk-jpl/src/jl/core/inputs.jl
===================================================================
--- /issm/trunk-jpl/src/jl/core/inputs.jl	(revision 26737)
+++ /issm/trunk-jpl/src/jl/core/inputs.jl	(revision 26738)
@@ -109,4 +109,9 @@
 
 end#}}}
+function GetInputMin(input::ElementInput) #{{{
+
+	return minimum(input.element_values)
+
+end#}}}
 function GetInputValue(input::ElementInput,gauss::GaussTria,i::Int64) #{{{
 
Index: /issm/trunk-jpl/src/jl/core/modules.jl
===================================================================
--- /issm/trunk-jpl/src/jl/core/modules.jl	(revision 26737)
+++ /issm/trunk-jpl/src/jl/core/modules.jl	(revision 26738)
@@ -49,16 +49,18 @@
 		#Configure objects
 		ConfigureObjectx(elements, nodes[i], vertices, parameters, inputs, i)
-
-		#Constrain and Number nodes
-		SpcNodesx(nodes[i], constraints[i], parameters)
-		NodesDofx(nodes[i], parameters)
 	end
 
 	#Build FemModel
-	return FemModel(analyses, elements, vertices,
-						 Vector{Node}(undef,0), nodes,
-						 parameters, inputs,
-						 Vector{Constraint}(undef,0), constraints,
-						 results)
+	femmodel = FemModel(analyses, elements, vertices,
+							  Vector{Node}(undef,0), nodes,
+							  parameters, inputs,
+							  Vector{Constraint}(undef,0), constraints,
+							  results)
+
+	println("      detecting active vertices")
+	GetMaskOfIceVerticesLSMx0(femmodel)
+
+	return femmodel
+
 end# }}}
 function CreateElements(elements::Vector{Tria},md::model) #{{{
@@ -205,4 +207,12 @@
 
 end#}}}
+function InputUpdateFromVectorx(femmodel::FemModel, vector::Vector{Float64}, enum::IssmEnum, layout::IssmEnum)
+
+	#Go through elements and plug in solution
+	for i=1:length(femmodel.elements)
+		InputUpdateFromVector(femmodel.elements[i], vector, enum, layout)
+	end
+
+end#}}}
 function InputDuplicatex(femmodel::FemModel, oldenum::IssmEnum, newenum::IssmEnum) #{{{
 	DuplicateInput(femmodel.inputs, oldenum, newenum)
@@ -266,6 +276,6 @@
 		pe = CreatePVector(analysis,femmodel.elements[i])
 
-		AddToGlobal!(Ke,Kff,Kfs)
-		AddToGlobal!(pe,pf)
+		if(!isnothing(Ke)) AddToGlobal!(Ke,Kff,Kfs) end
+		if(!isnothing(pe)) AddToGlobal!(pe,pf) end
 	end
 
@@ -320,2 +330,59 @@
 
 end# }}}
+function UpdateConstraintsx(femmodel::FemModel, analysis::Analysis)# {{{
+
+	#First, see if the analysis needs to change constraints
+	UpdateConstraints(analysis, femmodel)
+
+	#Second, constraints might be time dependent
+	SpcNodesx(femmodel.nodes, femmodel.constraints, femmodel.parameters)
+
+	#Now, update degrees of freedoms
+	NodesDofx(femmodel.nodes, femmodel.parameters)
+
+end# }}}
+function SetActiveNodesLSMx(femmodel::FemModel) #{{{
+
+	#Check mask of each element to see if element is active
+	numnodes = 3
+	mask = Vector{Float64}(undef, numnodes)
+	for i in 1:length(femmodel.elements)
+		GetInputListOnNodes!(femmodel.elements[i], mask, IceMaskNodeActivationEnum)
+		for j in 1:numnodes
+			node = femmodel.elements[i].nodes[j]
+			if(mask[j]==1.) Activate!(node)
+			else             Deactivate!(node)
+			end
+		end
+
+	end
+end#}}}
+function GetMaskOfIceVerticesLSMx0(femmodel::FemModel) #{{{
+
+	#Initialize vector with number of vertices
+	numvertices=length(femmodel.vertices)
+   if(numvertices==0) return end
+
+	#Initialize vector
+	nbv = 3
+	onesvec      = ones(nbv)
+	vec_mask_ice = IssmVector(numvertices)
+	vertexids    = Vector{Int64}(undef, nbv)
+
+	#Assign values to vector
+	for i in 1:length(femmodel.elements)
+		if (IsIceInElement(femmodel.elements[i]))
+			for j in 1:nbv
+				vertexids[j] = femmodel.elements[i].vertices[j].sid
+			end
+			SetValues!(vec_mask_ice, nbv, vertexids, onesvec)
+		end
+	end
+	Assemble!(vec_mask_ice)
+
+	#Serialize vector
+	vec_mask_ice_serial = ToSerial(vec_mask_ice)
+
+	#Update IceMaskNodeActivationEnum in elements
+	InputUpdateFromVectorx(femmodel, vec_mask_ice_serial, IceMaskNodeActivationEnum, VertexSIdEnum)
+end#}}}
Index: /issm/trunk-jpl/src/jl/core/nodes.jl
===================================================================
--- /issm/trunk-jpl/src/jl/core/nodes.jl	(revision 26737)
+++ /issm/trunk-jpl/src/jl/core/nodes.jl	(revision 26738)
@@ -4,4 +4,5 @@
 	sid::Int64
 	indexingupdate::Bool
+	active::Bool
 	gsize::Int64
 	gdoflist::Vector{Int64}
@@ -24,4 +25,30 @@
 	println(io,"   svalues: ",this.svalues)
 end# }}}
+function Activate!(node::Node) #{{{
+
+	if(!node.active)
+		node.indexingupdate = true
+		node.active = true
+		for i in 1:node.gsize
+			node.fdoflist[i]  = +1
+			node.sdoflist[i]  = -1
+			node.svalues[i]   = 0.0
+		end
+	end
+
+end# }}}
+function Deactivate!(node::Node) #{{{
+
+	if(node.active)
+		node.indexingupdate = true
+		node.active = false
+		for i in 1:node.gsize
+			node.fdoflist[i]  = -1
+			node.sdoflist[i]  = +1
+			node.svalues[i]   = 0.0
+		end
+	end
+
+end# }}}
 function ApplyConstraint(node::Node,dof::Int8,value::Float64) #{{{
 
@@ -42,10 +69,10 @@
 
 	if setenum==GsetEnum
-		for i=1:node.gsize
+		for i in 1:node.gsize
 			node.gdoflist[i] = dofcount
 			dofcount += 1
 		end
 	elseif setenum==FsetEnum
-		for i=1:node.gsize
+		for i in 1:node.gsize
 			if  node.fdoflist[i]!=-1
 				@assert node.sdoflist[i]==-1
@@ -55,5 +82,5 @@
 		end
 	elseif setenum==SsetEnum
-		for i=1:node.gsize
+		for i in 1:node.gsize
 			if  node.sdoflist[i]!=-1
 				@assert node.fdoflist[i]==-1
Index: /issm/trunk-jpl/src/jl/core/solutionsequences.jl
===================================================================
--- /issm/trunk-jpl/src/jl/core/solutionsequences.jl	(revision 26737)
+++ /issm/trunk-jpl/src/jl/core/solutionsequences.jl	(revision 26738)
@@ -1,3 +1,6 @@
 function solutionsequence_linear(femmodel::FemModel,analysis::Analysis) # {{{
+
+	#First, update constraints in case the levelset has changed
+	UpdateConstraintsx(femmodel, analysis)
 
 	#Get new matrices
@@ -20,4 +23,7 @@
 end# }}}
 function solutionsequence_nonlinear(femmodel::FemModel,analysis::Analysis,maxiter::Int64,restol::Float64,reltol::Float64,abstol::Float64) # {{{
+
+	#First, update constraints in case the levelset has changed
+	UpdateConstraintsx(femmodel, analysis)
 
 	#Initialize number of iterations
Index: /issm/trunk-jpl/src/jl/core/toolkits.jl
===================================================================
--- /issm/trunk-jpl/src/jl/core/toolkits.jl	(revision 26737)
+++ /issm/trunk-jpl/src/jl/core/toolkits.jl	(revision 26738)
@@ -92,4 +92,9 @@
 
 end#}}}
+function ToSerial(vector::IssmVector)#{{{
+
+	return vector.vector
+
+end#}}}
 function Norm(x::IssmVector,type::Int64)#{{{