source: issm/oecreview/Archive/21337-21723/ISSM-21573-21574.diff@ 21726

../trunk-jpl/src/c/classes/FemModel.cpp

@@ -2953 +2953 @@
 }
 /*}}}*/
 void FemModel::ReMesh(void){/*{{{*/
-        
+
         /*Variables*/
         IssmDouble *newx                        = NULL;
         IssmDouble *newy                        = NULL;
@@ -2972 +2972 @@
         /*Partitioning the new mesh. Maybe ElementsAndVerticesPartitioning.cpp could be modified to set this without iomodel.*/
         this->ElementsAndVerticesPartitioning(newnumberofvertices,newnumberofelements,elementswidth,newelementslist,&my_elements,&my_vertices);
 
-        if(this->loads->Size()!=0){
-                _error_("not supported yet");
-        }
+        if(this->loads->Size()!=0) _error_("not supported yet");
 
         /*Create vertices*/
         Vertices* new_vertices=new Vertices();
@@ -3008 +3006 @@
                 if(analysis_enum==StressbalanceAnalysisEnum) this->CreateConstraints(newnumberofvertices,newnumberofelements,nodecounter,constraintcounter,newx,newy,my_vertices,new_constraints);
                 this->UpdateElements(newnumberofelements,newelementslist,my_elements,nodecounter,i,new_elements);
 
-
-                if(IssmComm::GetRank()==2){
-                        _printf_("New constraints: CPU 1\n");
-                        //this->constraints->DeepEcho();
-                        new_constraints->DeepEcho();
-                        printarray(my_vertices,newnumberofvertices,1);
-                        _error_("stop!");
-                }
-
-
                 if(new_nodes->Size()) nodecounter=new_nodes->MaximumId();
                 constraintcounter = new_constraints->NumberOfConstraints();
                 /*Make sure nodecounter is at least 0 (if no node exists, maxid will be -1*/
@@ -3069 +3057 @@
         /*Finally: interpolate all inputs and insert them into the new elements.*/
         this->InterpolateInputs(new_vertices,new_elements);
 
-        if(IssmComm::GetRank()==2){
-                _printf_("New constraints: CPU 1\n");
-                //this->constraints->DeepEcho();
-                new_constraints->DeepEcho();
-                _error_("stop!");
-        }
-
-
         /*Delete old structure and set new pointers*/
         delete this->vertices;          this->vertices          = new_vertices;
         delete this->elements;          this->elements          = new_elements;
@@ -3630 +3610 @@
         Vector<IssmDouble>* vspcvxflag  = new Vector<IssmDouble>(numberofnodes_analysistype);
         Vector<IssmDouble>* vspcvyflag  = new Vector<IssmDouble>(numberofnodes_analysistype);
         
-        for(int i=0;i<numberofnodes_analysistype;i++){ 
-                vspcvx->SetValue(i,0.,INS_VAL);
-                vspcvy->SetValue(i,0.,INS_VAL);
-                vspcvxflag->SetValue(i,0.,INS_VAL);
-                vspcvyflag->SetValue(i,0.,INS_VAL);
-        }
-
         for(int i=0;i<this->constraints->Size();i++){
                 SpcStatic* spc                  = xDynamicCast<SpcStatic*>(this->constraints->GetObjectByOffset(i));
                 int dof                                 = spc->GetDof();
@@ -3644 +3617 @@
                 IssmDouble spcvalue     = spc->GetValue(); 
                 int nodeindex                   = node-1;
                 
-                if(IssmComm::GetRank()==0)      
-                  {//itapopo
-                 if(dof==0) {//vx
-                        if(IssmComm::GetRank()==0) _printf_("Node: " << nodeindex << "\tdof: " << dof << "\tspcvalue: " << spcvalue << "\n") ;
+                if(dof==0) {//vx
                         vspcvx->SetValue(nodeindex,spcvalue,INS_VAL);
                         vspcvxflag->SetValue(nodeindex,1,INS_VAL);
                 } 
@@ -3655 +3625 @@
                         vspcvy->SetValue(nodeindex,spcvalue,INS_VAL);
                         vspcvyflag->SetValue(nodeindex,1,INS_VAL);
                 }
-
-                  }//itapopo
         }
 
-        #ifdef _HAVE_PETSC_
-                _printf_("Tem PETSC\n");
-        #endif
-
         /*Assemble*/
         vspcvx->Assemble();
         vspcvy->Assemble();
         vspcvxflag->Assemble();
         vspcvyflag->Assemble();
 
-        IssmDouble pvalue1,pvalue2;
-        vspcvx->GetValue(&pvalue1,0);
-        vspcvxflag->GetValue(&pvalue2,0);
-        if(IssmComm::GetRank()==0) _printf_("Node: " << 0 << "\tspcvalue: " << pvalue1 << "\tflag: " << pvalue2 << "\n") ;
-        vspcvx->GetValue(&pvalue1,5);
-        vspcvxflag->GetValue(&pvalue2,5);
-        if(IssmComm::GetRank()==0) _printf_("Node: " << 5 << "\tspcvalue: " << pvalue1 << "\tflag: " << pvalue2 << "\n") ;
-        
-        
-        
         /*Serialize*/
         spcvx            = vspcvx->ToMPISerial();
         spcvy            = vspcvy->ToMPISerial();
@@ -3701 +3655 @@
         InterpFromMeshToMesh2dx(&newspcvxflag,elementslist,x,y,nods_data,nels_data,spcvxflag,M_data,N_data,newx,newy,N_interp,NULL);
         InterpFromMeshToMesh2dx(&newspcvyflag,elementslist,x,y,nods_data,nels_data,spcvyflag,M_data,N_data,newx,newy,N_interp,NULL);
 
-        if(IssmComm::GetRank()==1){
-                _printf_("CPU 1:\n");
-                _printf_("Old spcvx:\n");
-                printarray(spcvx,nods_data,1);
-                _printf_("Old spcvxflag:\n");
-                printarray(spcvxflag,nods_data,1);
-                _error_("Stop!!");
-        }
-
-
-
-
-
         int count                                       = 0;
         IssmDouble eps                          = 1.e-8;
 
         /*Now, insert the interpolated constraints in the data set (constraints)*/
         for(int i=0;i<newnumberofvertices;i++){
-                if(my_vertices[i])
-                /*spcvx*/
-                if(!xIsNan<IssmDouble>(newspcvx[i]) && newspcvxflag[i]>(1-eps)){
-                        constraints->AddObject(new SpcStatic(constraintcounter+count+1,nodecounter+i+1,0,newspcvx[i],StressbalanceAnalysisEnum));
-                        //add count'th spc, on node i+1, setting dof 1 to vx.
-                        count++;
+                if(my_vertices[i]){ 
+                        /*spcvx*/
+                        if(!xIsNan<IssmDouble>(newspcvx[i]) && newspcvxflag[i]>(1-eps)){
+                                constraints->AddObject(new SpcStatic(constraintcounter+count+1,nodecounter+i+1,0,newspcvx[i],StressbalanceAnalysisEnum));
+                                //add count'th spc, on node i+1, setting dof 1 to vx.
+                                count++;
+                        }
                 }
         }
         count=0;
         for(int i=0;i<newnumberofvertices;i++){
-                if(my_vertices[i])
-                /*spcvy*/
-                if(!xIsNan<IssmDouble>(newspcvy[i]) && newspcvyflag[i]>(1-eps) ){
-                        constraints->AddObject(new SpcStatic(constraintcounter+count+1,nodecounter+i+1,1,newspcvy[i],StressbalanceAnalysisEnum)); 
-                        //add count'th spc, on node i+1, setting dof 1 to vx.
-                        count++;
+                if(my_vertices[i]){
+                        /*spcvy*/
+                        if(!xIsNan<IssmDouble>(newspcvy[i]) && newspcvyflag[i]>(1-eps) ){
+                                constraints->AddObject(new SpcStatic(constraintcounter+count+1,nodecounter+i+1,1,newspcvy[i],StressbalanceAnalysisEnum)); 
+                                //add count'th spc, on node i+1, setting dof 1 to vx.
+                                count++;
+                        }
                 }
         }