Changeset 26652

Timestamp:

11/22/21 07:03:07 (3 years ago)

Author:

Mathieu Morlighem

Message:

CHG: more work on julia solver

Location:

issm/trunk-jpl/src/jl/solve

Files:

• analyses.jl (modified) (2 diffs)
• elements.jl (modified) (1 diff)
• modules.jl (modified) (1 diff)
• nodes.jl (modified) (1 diff)
• solutionsequences.jl (added)
• solve.jl (modified) (1 diff)

issm/trunk-jpl/src/jl/solve/analyses.jl

@@ -80 +80 @@
 
 end#}}}
+function InputUpdateFromSolution(analysis::StressbalanceAnalysis,ug::Vector{Float64},element::Tria) #{{{
 
+        #Get dofs for this finite element
+        doflist = GetDofList(element,GsetEnum)
+
+        #Get solution vector for this element
+        numdof   = 3*2
+        values = Vector{Float64}(undef,numdof)
+        for i in 1:numdof values[i]=ug[doflist[i]] end
+
+        #Now split solution vector into x and y components
+        numnodes = 3
+        vx= Vector{Float64}(undef,numnodes)
+        vy= Vector{Float64}(undef,numnodes)
+        for i in 1:numnodes 
+                vx[i]=values[2*i-1] 
+                vy[i]=values[2*i] 
+                @assert isfinite(vx[i])
+                @assert isfinite(vy[i])
+        end
+
+        AddInput(element,VxEnum,vx,P1Enum)
+        AddInput(element,VyEnum,vy,P1Enum)
+end#}}}
 function Core(analysis::StressbalanceAnalysis,femmodel::FemModel)# {{{
 
@@ -95 +118 @@
 
 
-function solutionsequence_nonlinear(femmodel::FemModel,analysis::Analysis,maxiter::Int64,restol::Float64,reltol::Float64,abstol::Float64) # {{{
-
-        #Initialize number of iterations
-        count = 0
-        converged = false
-
-        #Get existing solution
-        ug = GetSolutionFromInputsx(analysis,femmodel)
-
-        print(ug)
-
-        #Loop until we reach convergence
-        while(~converged)
-
-                error("not implemented yet")
-
-                #Increase count
-                count += 1
-                if(count>maxiter)
-                        println("   maximum number of nonlinear iterations (",maxiter,") exceeded")
-                        converged = true
-                end
-        end
-
-        error("STOP")
-
-end# }}}

issm/trunk-jpl/src/jl/solve/elements.jl

@@ -21 +21 @@
         else
                 error("size ",size(data,1)," not supported yet");
+        end
+end #}}}
+function AddInput(element::Tria,inputenum::IssmEnum,data::Vector{Float64},interpolation::IssmEnum) #{{{
+        if interpolation==P1Enum
+                @assert length(data)==3
+                SetTriaInput(element.inputs,inputenum,P1Enum,element.vertexids,data)
+        else
+                error("interpolation ", interpolation, " not supported yet");
         end
 end #}}}

issm/trunk-jpl/src/jl/solve/modules.jl

@@ -134 +134 @@
 
 end#}}}
+function InputUpdateFromSolutionx(analysis::Analysis,ug::Vector{Float64},femmodel::FemModel) #{{{
+
+        #Go through elements and plug in solution
+        for i=1:length(femmodel.elements)
+                InputUpdateFromSolution(analysis,ug,femmodel.elements[i])
+        end
+
+        return ug
+
+end#}}}
+function Reducevectorgtofx(ug::Vector{Float64},nodes::Vector{Node}) #{{{
+
+        #Get size of output vector
+        fsize = NumberOfDofs(nodes,FsetEnum)
+
+        #Initialize output vector
+        uf = Vector{Float64}(undef,fsize)
+
+        #Go through elements and plug in solution
+        for i=1:length(nodes)
+                VecReduce(nodes[i],ug,uf)
+        end
+
+        return ug
+
+end#}}}

issm/trunk-jpl/src/jl/solve/nodes.jl

@@ -114 +114 @@
 
 end# }}}
+function VecReduce(node::Node,ug::Vector{Float64},uf::Vector{Float64}) #{{{
+
+        for i=1:node.gsize
+                if node.fdoflist[i]!=-1
+                        uf[node.fdoflist[i]] = ug[node.gdoflist[i]]
+                end
+        end
+
+end# }}}
 
 #Nodes functions

issm/trunk-jpl/src/jl/solve/solve.jl

@@ -10 +10 @@
 include("./femmodel.jl")
 include("./analyses.jl")
+include("./solutionsequences.jl")
 include("./modules.jl")