Changeset 15448

Timestamp:

07/05/13 15:45:34 (12 years ago)

Author:

Mathieu Morlighem

Message:

BUG: fixed spc pb for quadratic elements. Only works for StaticSpc for now

Location:

issm/trunk-jpl/src/c

Files:

• classes/Constraints/SpcStatic.cpp (modified) (1 diff)
• modules/ModelProcessorx/DiagnosticHoriz/CreateConstraintsDiagnosticHoriz.cpp (modified) (4 diffs)

issm/trunk-jpl/src/c/classes/Constraints/SpcStatic.cpp

@@ -22 +22 @@
 SpcStatic::SpcStatic(int spc_sid,int spc_nodeid, int spc_dof,IssmDouble spc_value,int spc_analysis_type){
 
-        sid=spc_sid;
-        nodeid=spc_nodeid;
-        dof=spc_dof;
-        value=spc_value;
-        analysis_type=spc_analysis_type;
+        sid           = spc_sid;
+        nodeid        = spc_nodeid;
+        dof           = spc_dof;
+        value         = spc_value;
+        analysis_type = spc_analysis_type;
 
         return;

issm/trunk-jpl/src/c/modules/ModelProcessorx/DiagnosticHoriz/CreateConstraintsDiagnosticHoriz.cpp

@@ -18 +18 @@
         IssmDouble  stokesreconditioning;
         bool    isstokes,isl1l2,ismacayealpattyn;
+        int     fe_ssa;
    bool    spcpresent=false;
         int Mx,Nx;
@@ -50 +51 @@
         iomodel->Constant(&isl1l2,FlowequationIsl1l2Enum);
         iomodel->Constant(&ismacayealpattyn,FlowequationIsmacayealpattynEnum);
+        iomodel->Constant(&fe_ssa,FlowequationFeSsaEnum);
 
         /*Recover pointer: */
@@ -219 +221 @@
                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,spcvx[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
                                         count++;
+
                                 }
                                 else if (Mx==iomodel->numberofvertices+1) {
@@ -306 +309 @@
                                                 break;
                                         default: _error_("Vertex approximation " << EnumToStringx(reCast<int,IssmDouble>(vertices_type[i])) << " not supported");
+                                }
+                        }
+                }
+        }
+
+        /*SPC Quadratic elements*/
+        if(ismacayealpattyn&& fe_ssa==1){
+
+                int v1,v2,el1,el2;
+                bool my_edge;
+
+                if(Mx!=iomodel->numberofvertices) _error_("transient spc not supported yet");
+                CreateEdges(iomodel);
+
+                for(i=0;i<iomodel->numberofedges;i++){
+
+                        v1 = iomodel->edges[4*i+0]-1;
+                        v2 = iomodel->edges[4*i+1]-1;
+                        el1= iomodel->edges[4*i+2]-1;
+                        el2= iomodel->edges[4*i+3]-1;
+
+                        my_edge = iomodel->my_elements[el1];
+                        if(!my_edge && el2>=0){
+                                my_edge = iomodel->my_elements[el2];
+                        }
+
+                        if(my_edge){
+                                if(!xIsNan<IssmDouble>(spcvx[v1]) && !xIsNan<IssmDouble>(spcvx[v2])){
+                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+iomodel->numberofvertices+i+1,
+                                                                        1,(spcvx[v1]+spcvx[v2])/yts,DiagnosticHorizAnalysisEnum));
+                                        count++;
+                                }
+                                if(!xIsNan<IssmDouble>(spcvx[v1]) && !xIsNan<IssmDouble>(spcvx[v2])){
+                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+iomodel->numberofvertices+i+1,
+                                                                        2,(spcvy[v1]+spcvy[v2])/yts,DiagnosticHorizAnalysisEnum));
+                                        count++;
+                                }
+                                if(!xIsNan<IssmDouble>(spcvx[v1]) && !xIsNan<IssmDouble>(spcvx[v2])){
+                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+iomodel->numberofvertices+i+1,
+                                                                        3,(spcvz[v1]+spcvz[v2])/yts,DiagnosticHorizAnalysisEnum));
+                                        count++;
                                 }
                         }