NEW: adding transient spc

src/core/analyses/levelsetanalysis.jl

@@ -3,19 +3,41 @@ struct LevelsetAnalysis <: Analysis#{{{
 end #}}}
 
 #Model Processing
-function CreateConstraints(analysis::LevelsetAnalysis,constraints::Vector{Constraint},md::model) #{{{
+function CreateConstraints(analysis::LevelsetAnalysis,constraints::Vector{AbstractConstraint},md::model) #{{{
 
 	#load constraints from model
 	spclevelset = md.levelset.spclevelset
-	@assert size(spclevelset,1)==md.mesh.numberofvertices
-	@assert size(spclevelset,2)==1
 
 	count = 1
-	for i in 1:md.mesh.numberofvertices
-		if ~isnan(spclevelset[i])
-			push!(constraints,Constraint(count,i,1,spclevelset[i]))
-			count+=1
+	if size(spclevelset,1)==md.mesh.numberofvertices
+		@assert size(spclevelset,2)==1
+		for i in 1:md.mesh.numberofvertices
+			if ~isnan(spclevelset[i])
+				push!(constraints,Constraint(count,i,1,spclevelset[i]))
+				count+=1
+			end
 		end
+
+	elseif size(spclevelset,1)==md.mesh.numberofvertices+1
+
+		times = spclevelset[end,:]
+		N     = size(spclevelset, 2)
+
+		for i in 1:md.mesh.numberofvertices
+			spcpresent = false
+			for j in 1:N
+				if(~isnan(spclevelset[i,j]))
+					spcpresent = true
+					break
+				end
+			end
+			if spcpresent
+				push!(constraints, ConstraintTransient(count, i, 1, times, spclevelset[i,:]))
+				count+=1
+			end
+		end
+	else
+		error("Size of spclevelset not supported yet")
 	end
 
 	return nothing

src/core/analyses/masstransportanalysis.jl

@@ -3,7 +3,7 @@ struct MasstransportAnalysis <: Analysis#{{{
 end #}}}
 
 #Model Processing
-function CreateConstraints(analysis::MasstransportAnalysis,constraints::Vector{Constraint},md::model) #{{{
+function CreateConstraints(analysis::MasstransportAnalysis,constraints::Vector{AbstractConstraint},md::model) #{{{
 
 	#load constraints from model
 	spcthickness = md.masstransport.spcthickness

src/core/analyses/stressbalanceanalysis.jl

@@ -3,7 +3,7 @@ struct StressbalanceAnalysis <: Analysis#{{{
 end #}}}
 
 #Model Processing
-function CreateConstraints(analysis::StressbalanceAnalysis,constraints::Vector{Constraint},md::model) #{{{
+function CreateConstraints(analysis::StressbalanceAnalysis,constraints::Vector{AbstractConstraint},md::model) #{{{
 
 	#load constraints from model
 	spcvx = md.stressbalance.spcvx

src/core/constraints.jl

@@ -1,19 +1,65 @@
 #Constraint class definition
-struct Constraint #{{{
+abstract type AbstractConstraint end
+struct Constraint <: AbstractConstraint #{{{
 	id::Int64
 	nodeid::Int64
 	dof::Int8
 	value::Float64
 end# }}}
+struct ConstraintTransient <: AbstractConstraint #{{{
+	id::Int64
+	nodeid::Int64
+	dof::Int8
+	times::Vector{Float64}
+	values::Vector{Float64}
+end# }}}
 
 #Constraint functions
 function ConstrainNode(constraint::Constraint,nodes::Vector{Node},parameters::Parameters) #{{{
 
-	#Chase through nodes and find the node to which this SpcStatic apply
+	#Chase through nodes and find the node to which this SpcStatic applies
 	node = nodes[constraint.nodeid]
 
 	#Apply Constraint
 	ApplyConstraint(node, constraint.dof, constraint.value)
 
 	return nothing
 end# }}}
+function ConstrainNode(constraint::ConstraintTransient, nodes::Vector{Node}, parameters::Parameters) #{{{
+
+	#Chase through nodes and find the node to which this SpcTransient applies
+	node = nodes[constraint.nodeid]
+
+	#Find time in parameters
+	time = FindParam(Float64, parameters, TimeEnum)
+
+	# Now, go fetch value for this time:
+	if (time<=times[1])
+		value = values[1]
+		found = true
+	elseif (time>=times[end])
+		value = values[end]
+		found = true
+	else
+		for i in 1:length(times)-1
+			if (times[i]<=time && time<times[i+1])
+				alpha = (time-times[i])/(times[i+1]-times[i]);
+				@assert alpha>=0.0 && alpha<=1.0
+				value = (1-alpha)*values[i] + alpha*values[i+1];
+				found = true
+				break
+			end
+		end
+	end
+
+	if(!found) error("could not find time segment for constraint") end
+
+	#Apply Constraint
+	if isnan(value)
+		RelaxConstraint(node, constraint.dof)
+	else
+		ApplyConstraint(node, constraint.dof, constraint.value)
+	end
+
+	return nothing
+end# }}}

src/core/femmodel.jl

@@ -10,8 +10,8 @@ mutable struct FemModel #{{{
 	parameters::Parameters
 	inputs::Inputs
 
-	constraints::Vector{Constraint}
-	constraints_list::Vector{Vector{Constraint}}
+	constraints::Vector{AbstractConstraint}
+	constraints_list::Vector{Vector{AbstractConstraint}}
 
 	#loads::Vector{Loads}
 

src/core/modules.jl

@@ -45,7 +45,7 @@ function ModelProcessor(md::model, solutionstring::Symbol) #{{{
 	#Initialize analysis specific datasets
 	numanalyses = length(analyses)
 	nodes       = Vector{Vector{Node}}(undef,numanalyses)
-	constraints = Vector{Vector{Constraint}}(undef,numanalyses)
+	constraints = Vector{Vector{AbstractConstraint}}(undef,numanalyses)
 	for i in 1:numanalyses
 		analysis = analyses[i]
 		println("   creating datasets for analysis ", typeof(analysis))
@@ -431,7 +431,7 @@ function SetActiveNodesLSMx(femmodel::FemModel) #{{{
 	end
 	return nothing
 end#}}}
-function SpcNodesx(nodes::Vector{Node},constraints::Vector{Constraint},parameters::Parameters) #{{{
+function SpcNodesx(nodes::Vector{Node},constraints::Vector{AbstractConstraint},parameters::Parameters) #{{{
 
 	for i in 1:length(constraints)
 		ConstrainNode(constraints[i],nodes,parameters)

src/core/nodes.jl

@@ -167,6 +167,14 @@ function GetGlobalDofList(nodes::Vector{Node},ndofs::Int64,setenum::IssmEnum) #{
 	return doflist
 
 end# }}}
+function RelaxConstraint(node::Node, dof::Int8) #{{{
+
+	node.indexingupdate = true
+	node.fdoflist[dof]  = +1
+	node.sdoflist[dof]  = -1
+
+	return nothing
+end# }}}
 function VecReduce(node::Node,ug::Vector{Float64},uf::IssmVector) #{{{
 
 	for i=1:node.gsize

src/usr/classes.jl

@@ -280,7 +280,7 @@ end# }}}
 # }}}
 #Levelset{{{
 mutable struct Levelset
-	spclevelset::Vector{Float64}
+	spclevelset::Union{Vector{Float64},Matrix{Float64}}
 	stabilization::Int64
 	reinit_frequency::Int64
 	kill_icebergs::Int64