Changeset 5998

Timestamp:

09/23/10 15:02:51 (15 years ago)

Author:

Mathieu Morlighem

Message:

merged the two solvers

Location:

issm/trunk/src/m/solvers

Files:

• solver_diagnostic_nonlinear.m (modified) (1 diff)
• solver_diagnostic_nonlinear_kff.m (deleted)
• solver_diagnostic_nonlinear_kgg.m (deleted)

issm/trunk/src/m/solvers/solver_diagnostic_nonlinear.m

@@ -7 +7 @@
         %Branch on partitioning schema requested
         kffpartitioning=femmodel.parameters.Kff;
+        min_mechanical_constraints=femmodel.parameters.MinMechanicalConstraints;
 
-        if ~kffpartitioning,
-                femmodel=solver_diagnostic_nonlinear_kgg(femmodel,conserve_loads);
-        else
-                femmodel=solver_diagnostic_nonlinear_kff(femmodel,conserve_loads);
+        %keep a copy of loads for now
+        loads=femmodel.loads;
+
+        %initialize solution vector
+        converged=0; count=1;
+
+        %Start non-linear iteration using input velocity: 
+        ug=GetSolutionFromInputs(femmodel.elements, femmodel.nodes, femmodel.vertices, loads, femmodel.materials, femmodel.parameters);
+        uf=Reducevectorgtof( ug, femmodel.nodesets,femmodel.parameters);
+
+        %Update the solution to make sure that vx and vxold are similar
+        [femmodel.elements,femmodel.materials]=InputUpdateFromSolution(femmodel.elements,femmodel.nodes,femmodel.vertices,loads,femmodel.materials,femmodel.parameters,ug);
+
+        while(~converged),
+
+                %save pointer to old velocity
+                old_ug=ug;
+                old_uf=uf;
+
+                if kffpartitioning,
+                        disp('-------------- file: solver_diagnostic_nonlinear.m line: 31'); 
+                        [K_gg,K_ff,K_fs,p_g,p_f,kmax]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,loads,femmodel.materials,femmodel.parameters);
+                        p_f = Reduceload( p_f, K_fs, femmodel.ys);
+                else
+                        disp('-------------- file: solver_diagnostic_nonlinear.m line: 35'); 
+                        [K_gg,K_ff,K_fs,p_g,p_f,kmax]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,loads,femmodel.materials,femmodel.parameters);
+                        [K_ff, K_fs] = Reducematrixfromgtof( K_gg, femmodel.nodesets,femmodel.parameters); 
+                        p_f = Reduceloadfromgtof( p_g, K_fs, femmodel.ys, femmodel.nodesets,femmodel.parameters);
+                end
+
+                uf=Solver(K_ff,p_f,old_uf,femmodel.parameters);
+                ug= Mergesolutionfromftog( uf, femmodel.ys, femmodel.nodesets,femmodel.parameters); 
+
+                [femmodel.elements,femmodel.materials]=InputUpdateFromSolution(femmodel.elements,femmodel.nodes,femmodel.vertices,loads,femmodel.materials,femmodel.parameters,ug);
+                [loads,constraints_converged,num_unstable_constraints] =PenaltyConstraints( femmodel.elements,femmodel.nodes,femmodel.vertices,loads, femmodel.materials,femmodel.parameters);
+
+                displaystring(femmodel.parameters.Verbose,'%s%i','      number of unstable constraints: ',num_unstable_constraints);
+
+                %Figure out if convergence have been reached
+                converged=convergence(K_ff,p_f,uf,old_uf,femmodel.parameters);
+
+                %add convergence status into  status
+                [femmodel.elements loads]=InputUpdateFromConstant(femmodel.elements,femmodel.nodes,femmodel.vertices,loads,femmodel.materials,femmodel.parameters,double(converged),ConvergedEnum);
+
+                %rift convergence
+                if ~constraints_converged,
+                        if converged,
+                                if num_unstable_constraints <= min_mechanical_constraints,
+                                        converged=1;
+                                else 
+                                        converged=0;
+                                end
+                        end
+                end
+
+                %increase count
+                count=count+1;
+                if(converged==1)break;
+                        if(count>max_nonlinear_iterations),
+                                displaystring(femmodel.parameters.Verbose,'%s%i%s','      maximum number of iterations ',max_nonlinear_iterations,' exceeded');
+                        end
+                end
+
+        end
+
+        %deal with loads:
+        if conserve_loads==false,
+                femmodel.loads=loads;
         end
 end