Changeset 17250

Timestamp:

02/10/14 15:37:11 (11 years ago)

Author:

jbondzio

Message:

BUG: minor bugs; transient coupling between stressbalance and lsm

Location:

issm/trunk-jpl/src/c

Files:

• analyses/LevelsetAnalysis.cpp (modified) (2 diffs)
• analyses/StressbalanceAnalysis.cpp (modified) (1 diff)
• cores/transient_core.cpp (modified) (1 diff)
• modules/SetActiveNodesLSMx/SetActiveNodesLSMx.cpp (modified) (2 diffs)

issm/trunk-jpl/src/c/analyses/LevelsetAnalysis.cpp

@@ -66 +66 @@
         /*parameters: */
         bool save_results;
-
-        /* extrapolate */
-        Analysis* analysis = new ExtrapolationAnalysis();
-        femmodel->parameters->SetParam(VxEnum,ExtrapolationVariableEnum);
-        analysis->Core(femmodel);
-        delete analysis;
 
         /*activate formulation: */
@@ -165 +159 @@
                                 /* Artificial Diffusion */
                                 element->ElementSizes(&hx,&hy,&hz);
-                                vel=sqrt(vx*vx + vy*vy );
+                                vel=sqrt(vx*vx + vy*vy) + 1e-14;
                                 h=sqrt( pow(hx*vx/vel,2) + pow(hy*vy/vel,2) ); //FIXME: is this correct?
 

issm/trunk-jpl/src/c/analyses/StressbalanceAnalysis.cpp

@@ -1044 +1044 @@
 void StressbalanceAnalysis::UpdateConstraints(FemModel* femmodel){/*{{{*/
         /*Default, do nothing*/
-        _printf0_("   Updating active and non-active nodes in Femmodel \n");
-        //SetActiveNodesLSMx(femmodel->elements);
+        _printf0_("   Updating active and non-active nodes for StressbalanceAnalysis \n");
+        SetActiveNodesLSMx(femmodel->elements);
         return;
 }/*}}}*/

issm/trunk-jpl/src/c/cores/transient_core.cpp

@@ -150 +150 @@
                 if(islevelset){
                         if(VerboseSolution()) _printf0_("   computing movement of ice boundaries\n");
+
+                        /* extrapolate required variables */
+                        Analysis* extanalysis = new ExtrapolationAnalysis();
+                        int vars[2] = {VxEnum, VyEnum};
+                        for(int iv=0;iv<2;iv++){
+                                femmodel->parameters->SetParam(vars[iv],ExtrapolationVariableEnum); 
+                                extanalysis->Core(femmodel);
+                        }
+                        delete extanalysis;     
+
                         analysis = new LevelsetAnalysis();
                         analysis->Core(femmodel);
                         delete analysis;
+
+                        /* update vertices included for next calculation */
+                        GetMaskOfIceVerticesLSMx(femmodel);
                 }
 

issm/trunk-jpl/src/c/modules/SetActiveNodesLSMx/SetActiveNodesLSMx.cpp

@@ -16 +16 @@
                 int         numnodes = element->GetNumberOfNodes();
                 IssmDouble *mask     = xNew<IssmDouble>(numnodes);
-                element->GetInputListOnNodes(&mask[0],IceMaskNodeActivationEnum);
-                
+                // include switch for elements with multiple different sets of nodes
+                switch(element->GetElementType()){
+                        case MINIEnum:case TaylorHoodEnum:{
+                                Input* input=element->GetInput(IceMaskNodeActivationEnum);
+                                if(!input) _error_("Input " << EnumToStringx(IceMaskNodeActivationEnum) << " not found in element");
+
+                                /* Start looping on the number of vertices: */
+                                Gauss* gauss=element->NewGauss();
+                                for(int iv=0;iv<element->NumberofNodesVelocity();iv++){
+                                        gauss->GaussNode(element->VelocityInterpolation(),iv);
+                                        input->GetInputValue(&mask[iv],gauss);
+                                }
+                                for(int iv=0;iv<element->NumberofNodesPressure();iv++){
+                                        gauss->GaussNode(element->PressureInterpolation(),iv);
+                                        input->GetInputValue(&mask[element->NumberofNodesVelocity()+iv],gauss);
+                                }
+                                delete gauss;
+                                break;
+                        }
+                        default:
+                                element->GetInputListOnNodes(&mask[0],IceMaskNodeActivationEnum);
+                                break;
+                }
+
                 for(int in=0;in<numnodes;in++){
                         Node* node=element->GetNode(in);
@@ -33 +55 @@
 void GetMaskOfIceVerticesLSMx(FemModel* femmodel){/*{{{*/
 
+        /* Intermediaries */
+        int i;
         /*Initialize vector with number of vertices*/
         int numvertices=femmodel->vertices->NumberOfVertices();
         Vector<IssmDouble>* vec_mask_ice=new Vector<IssmDouble>(numvertices); //vertices that have ice at next time step
-
+        for(i=0;i<numvertices;i++){
+                vec_mask_ice->SetValue(i,0.,INS_VAL);
+        }
         /*Fill vector with values: */
-        for(int i=0;i<femmodel->elements->Size();i++){
+        for(i=0;i<femmodel->elements->Size();i++){
                 Element* element=dynamic_cast<Element*>(femmodel->elements->GetObjectByOffset(i));
                 SetMaskOfIceElement(vec_mask_ice, element);