source: issm/oecreview/Archive/22819-23185/ISSM-23172-23173.diff

Index: ../trunk-jpl/src/c/modules/SetActiveNodesLSMx/SetActiveNodesLSMx.cpp
===================================================================
--- ../trunk-jpl/src/c/modules/SetActiveNodesLSMx/SetActiveNodesLSMx.cpp        (revision 23172)
+++ ../trunk-jpl/src/c/modules/SetActiveNodesLSMx/SetActiveNodesLSMx.cpp        (revision 23173)
@@ -55,8 +55,9 @@
                xDelete<IssmDouble>(mask);
        }
 }/*}}}*/
-void GetMaskOfIceVerticesLSMx(FemModel* femmodel){/*{{{*/
 
+void GetMaskOfIceVerticesLSMx0(FemModel* femmodel){/*{{{*/
+
        /*Initialize vector with number of vertices*/
        int numvertices=femmodel->vertices->NumberOfVertices();
        if(numvertices==0)  return;
@@ -82,3 +83,42 @@
        delete vec_mask_ice;
        xDelete<IssmDouble>(mask_ice);
 }/*}}}*/
+void GetMaskOfIceVerticesLSMx(FemModel* femmodel){/*{{{*/
+
+       femmodel->SetCurrentConfiguration(LevelsetAnalysisEnum);
+
+       /*Create vector on gset*/
+       int configuration_type;
+       femmodel->parameters->FindParam(&configuration_type,ConfigurationTypeEnum);
+       int gsize=femmodel->nodes->NumberOfDofs(configuration_type,GsetEnum);
+       if(gsize==0)  return;
+       Vector<IssmDouble>* vec_mask_ice=new Vector<IssmDouble>(gsize);
+
+       /*Fill vector with values: */
+       for(int i=0;i<femmodel->elements->Size();i++){
+               Element* element=xDynamicCast<Element*>(femmodel->elements->GetObjectByOffset(i));
+
+               if(element->IsIceInElement()){
+                       int numnodes = element->GetNumberOfNodes();
+                       int  gsize_local=GetNumberOfDofs(element->nodes,numnodes,GsetEnum,NoneEnum);
+                       int* glist_local=GetGlobalDofList(element->nodes,numnodes,GsetEnum,NoneEnum);
+                       IssmDouble* ones = xNew<IssmDouble>(gsize_local);
+                       for(int n=0;n<gsize_local;n++) ones[n] = 1.;
+                       vec_mask_ice->SetValues(gsize_local,glist_local,ones,INS_VAL);
+                       xDelete<IssmDouble>(ones);
+                       xDelete<int>(glist_local);
+               }
+       }
+
+       /*Assemble vector and serialize */
+       vec_mask_ice->Assemble();
+       IssmDouble* mask_ice=vec_mask_ice->ToMPISerial();
+       delete vec_mask_ice;
+       for(int i=0;i<femmodel->elements->Size();i++){
+               Element* element=xDynamicCast<Element*>(femmodel->elements->GetObjectByOffset(i));
+               element->InputUpdateFromSolutionOneDof(mask_ice,IceMaskNodeActivationEnum);
+       }
+
+       /*Clean up and return*/
+       xDelete<IssmDouble>(mask_ice);
+}/*}}}*/
Index: ../trunk-jpl/src/c/modules/SetActiveNodesLSMx/SetActiveNodesLSMx.h
===================================================================
--- ../trunk-jpl/src/c/modules/SetActiveNodesLSMx/SetActiveNodesLSMx.h  (revision 23172)
+++ ../trunk-jpl/src/c/modules/SetActiveNodesLSMx/SetActiveNodesLSMx.h  (revision 23173)
@@ -8,5 +8,6 @@
 #include "../../classes/classes.h"
 
 void SetActiveNodesLSMx(FemModel* femmodel);
+void GetMaskOfIceVerticesLSMx0(FemModel* femmodel);
 void GetMaskOfIceVerticesLSMx(FemModel* femmodel);
-#endif  /* _UPDATESPCSX_H */
+#endif  /* _SETACTIVENODESLSMX_H*/
Index: ../trunk-jpl/src/c/classes/FemModel.cpp
===================================================================
--- ../trunk-jpl/src/c/classes/FemModel.cpp     (revision 23172)
+++ ../trunk-jpl/src/c/classes/FemModel.cpp     (revision 23173)
@@ -404,7 +404,7 @@
                        VerticesDofx(vertices,parameters);
 
                        if(VerboseMProcessor()) _printf0_("      detecting active vertices\n");
-                       GetMaskOfIceVerticesLSMx(this);
+                       GetMaskOfIceVerticesLSMx0(this);
                }
 
                if(VerboseMProcessor()) _printf0_("      resolving node constraints\n");
@@ -2632,7 +2632,7 @@
        delete this->constraints;       this->constraints       = new_constraints;
        delete this->materials;         this->materials = new_materials;
 
-       GetMaskOfIceVerticesLSMx(this);
+       GetMaskOfIceVerticesLSMx0(this);
 
        /*Insert MISMIP+ bed topography FIXME it could be stay in another place*/
        this->parameters->FindParam(&basalforcing_model,BasalforcingsEnum);
Index: ../trunk-jpl/src/c/classes/Elements/Tria.cpp
===================================================================
--- ../trunk-jpl/src/c/classes/Elements/Tria.cpp        (revision 23172)
+++ ../trunk-jpl/src/c/classes/Elements/Tria.cpp        (revision 23173)
@@ -3958,12 +3958,18 @@
 
        _assert_(fieldvalue==0.); //field value != 0 not implemented yet
 
+       /*Get field on vertices (we do not allow for higher order elements!!)*/
+       IssmDouble lsf[NUMVERTICES];
+       this->GetInputListOnVertices(&lsf[0],fieldenum);
+
        /*1. check that we do cross fieldvalue in this element*/
-       Input* input = inputs->GetInput(fieldenum);
-       if(!input) _error_("Cannot calculate distance to "<<EnumToStringx(fieldenum)<<", input not found");
-       IssmDouble minvalue = input->Min();
+       IssmDouble minvalue = lsf[0];
+       IssmDouble maxvalue = lsf[0];
+       for(int i=1;i<NUMVERTICES;i++){
+               if(lsf[i]>maxvalue) maxvalue = lsf[i];
+               if(lsf[i]<minvalue) minvalue = lsf[i];
+       }
        if(minvalue>fieldvalue) return;
-       IssmDouble maxvalue = input->Max();
        if(maxvalue<fieldvalue) return;
 
        /*2. Find coordinates of where levelset crosses 0*/
@@ -3987,8 +3993,7 @@
                }
        }
        else if(numiceverts==NUMVERTICES){ //NUMVERTICES ice vertices: calving front lies on element edge
-               IssmDouble lsf[NUMVERTICES];
-               this->GetInputListOnVertices(&lsf[0],fieldenum);
+
                for(int i=0;i<NUMVERTICES;i++){
                        if(lsf[indices[i]]==0.){
                                x[counter]=xyz_list[indices[i]][0];
Index: ../trunk-jpl/src/c/classes/Nodes.cpp
===================================================================
--- ../trunk-jpl/src/c/classes/Nodes.cpp        (revision 23172)
+++ ../trunk-jpl/src/c/classes/Nodes.cpp        (revision 23173)
@@ -348,11 +348,7 @@
        /*Grab max of all cpus: */
        ISSM_MPI_Allreduce((void*)&flag,(void*)&allflag,1,ISSM_MPI_INT,ISSM_MPI_MAX,IssmComm::GetComm());
 
-       if(allflag){
-               return true;
-       }
-       else{
-               return false;
-       }
+       if(allflag) return true;
+       else        return false;
 }
 /*}}}*/