source: issm/oecreview/Archive/26740-27031/ISSM-27030-27031.diff@ 27032

../trunk-jpl/src/m/boundaryconditions/SetMLHOBC.py

@@ -1 @@
-import numpy as np
-
-
-def SetMLHOBC(md):
-    """
-    SETMLHOBC - Create the boundary conditions for stressbalance for MLHO: VxBase, VyBase, VxShear, VyShear
-
-       Usage:
-          md = SetIceShelfBC(md, varargin)
-
-       Example:
-          md = SetIceShelfBC(md)
-
-    """
-
-    #node on Dirichlet (boundary and ~icefront)
-    md.stressbalance.spcvx_base = md.stressbalance.spcvx
-    md.stressbalance.spcvy_base = md.stressbalance.spcvy
-    md.stressbalance.spcvx_shear = np.nan * md.stressbalance.spcvx
-    md.stressbalance.spcvy_shear = np.nan * md.stressbalance.spcvy
-
-    return md

../trunk-jpl/src/m/boundaryconditions/SetMLHOBC.m

@@ -1 @@
-function md=SetMLHOBC(md)
-%SETMLHOBC - Create the boundary conditions for stressbalance for MLHO: VxBase, VyBase, VxShear, VyShear
-%
-%   Usage:
-%      md=SetMLHOBC(md)
-%
-
-
-%node on Dirichlet
-if md.flowequation.isMLHO
-        md.stressbalance.spcvx_base=md.stressbalance.spcvx;
-        md.stressbalance.spcvy_base=md.stressbalance.spcvy;
-
-        md.stressbalance.spcvx_shear=NaN*ones(size(md.stressbalance.spcvx_base));
-        md.stressbalance.spcvy_shear=NaN*ones(size(md.stressbalance.spcvy_base));
-end

../trunk-jpl/src/m/classes/flowequation.m

@@ -8 +8 @@
                 isSIA                          = 0;
                 isSSA                          = 0;
                 isL1L2                         = 0;
-                isMLHO                         = 0;
+                isMOLHO                         = 0;
                 isHO                           = 0;
                 isFS                           = 0;
                 isNitscheBC                    = 0;
@@ -59 +59 @@
                                 if isfield(objstruct,'bordermacayeal'),self.borderSSA = objstruct.bordermacayeal; end; 
                                 if isfield(objstruct,'borderpattyn'),  self.borderHO  = objstruct.borderpattyn;   end; 
                                 if isfield(objstruct,'borderstokes'),  self.borderFS  = objstruct.borderstokes;   end; 
-                        end
 
-                        %Nov 6 2021
-                        if any(self.vertex_equation==4)
-                                disp(['Monolayer Higher-Order (MLHO) detected in md.flowequation, this is still under development. Please double check your settings.']);
+                                %May 31 2022
+                                if isfield(objstruct,'isMLHO')
+                                        self.isMOLHO = objstruct.isMLHO;
+                                end
                         end
 
+
                 end% }}}
         end
         methods
@@ -103 +104 @@
                         md = checkfield(md,'fieldname','flowequation.isSIA','numel',[1],'values',[0 1]);
                         md = checkfield(md,'fieldname','flowequation.isSSA','numel',[1],'values',[0 1]);
                         md = checkfield(md,'fieldname','flowequation.isL1L2','numel',[1],'values',[0 1]);
-                        md = checkfield(md,'fieldname','flowequation.isMLHO','numel',[1],'values',[0 1]);
+                        md = checkfield(md,'fieldname','flowequation.isMOLHO','numel',[1],'values',[0 1]);
                         md = checkfield(md,'fieldname','flowequation.isHO','numel',[1],'values',[0 1]);
                         md = checkfield(md,'fieldname','flowequation.isFS','numel',[1],'values',[0 1]);
                         md = checkfield(md,'fieldname','flowequation.isNitscheBC','numel',[1],'values',[0 1]);
@@ -134 +135 @@
                         else
                                 error('Case not supported yet');
                         end
-                        if ~(self.isSIA || self.isSSA || self.isL1L2 || self.isMLHO || self.isHO || self.isFS),
+                        if ~(self.isSIA || self.isSSA || self.isL1L2 || self.isMOLHO || self.isHO || self.isFS),
                                 md = checkmessage(md,['no element types set for this model']);
                         end
                         if ismember('StressbalanceSIAAnalysis',analyses),
@@ -152 +153 @@
                         fielddisplay(self,'isSIA','is the Shallow Ice Approximation (SIA) used?');
                         fielddisplay(self,'isSSA','is the Shelfy-Stream Approximation (SSA) used?');
                         fielddisplay(self,'isL1L2','is the L1L2 approximation used?');
-                        fielddisplay(self,'isMLHO','is the Mono-Layer Higher-Order approximation used?');
+                        fielddisplay(self,'isMOLHO','is the MOno-Layer Higher-Order (MOLHO) approximation used?');
                         fielddisplay(self,'isHO','is the Higher-Order (HO) approximation used?');
                         fielddisplay(self,'isFS','are the Full-FS (FS) equations used?');
                         fielddisplay(self,'isNitscheBC','is weakly imposed condition used?');
@@ -171 +172 @@
                         WriteData(fid,prefix,'object',self,'fieldname','isSIA','format','Boolean');
                         WriteData(fid,prefix,'object',self,'fieldname','isSSA','format','Boolean');
                         WriteData(fid,prefix,'object',self,'fieldname','isL1L2','format','Boolean');
-                        WriteData(fid,prefix,'object',self,'fieldname','isMLHO','format','Boolean');
+                        WriteData(fid,prefix,'object',self,'fieldname','isMOLHO','format','Boolean');
                         WriteData(fid,prefix,'object',self,'fieldname','isHO','format','Boolean');
                         WriteData(fid,prefix,'object',self,'fieldname','isFS','format','Boolean');
                         WriteData(fid,prefix,'object',self,'fieldname','isNitscheBC','format','Boolean');
@@ -196 +197 @@
                         writejsdouble(fid,[modelname '.flowequation.isSIA'],self.isSIA);
                         writejsdouble(fid,[modelname '.flowequation.isSSA'],self.isSSA);
                         writejsdouble(fid,[modelname '.flowequation.isL1L2'],self.isL1L2);
-                        writejsdouble(fid,[modelname '.flowequation.isMLHO'],self.isMLHO);
+                        writejsdouble(fid,[modelname '.flowequation.isMOLHO'],self.isMOLHO);
                         writejsdouble(fid,[modelname '.flowequation.isHO'],self.isHO);
                         writejsdouble(fid,[modelname '.flowequation.isFS'],self.isFS);
          writejsstring(fid,[modelname '.flowequation.isNitscheBC'],self.isNitscheBC);

../trunk-jpl/src/c/analyses/AdjointHorizAnalysis.cpp

@@ -51 +51 @@
                         return CreateKMatrixHO(element);
                 case FSApproximationEnum: 
                         return CreateKMatrixFS(element);
-      case MLHOApproximationEnum:
+      case MOLHOApproximationEnum:
                 // a more accurate option, but integrate in the vertical direction numerically.
-                //      return CreateKMatrixMLHOVerticalIntergrated(element);
-                        return CreateKMatrixMLHO(element);
+                //      return CreateKMatrixMOLHOVerticalIntergrated(element);
+                        return CreateKMatrixMOLHO(element);
                 case NoneApproximationEnum:
                         return NULL;
                 default:
@@ -219 +219 @@
         xDelete<IssmDouble>(xyz_list);
         return Ke;
 }/*}}}*/
-ElementMatrix* AdjointHorizAnalysis::CreateKMatrixMLHO(Element* element){/*{{{*/
+ElementMatrix* AdjointHorizAnalysis::CreateKMatrixMOLHO(Element* element){/*{{{*/
 
         /* Check if ice in element */
         if(!element->IsIceInElement()) return NULL;
@@ -291 +291 @@
       vxshear_input->GetInputValue(&vxshear,gauss);
       vyshear_input->GetInputValue(&vyshear,gauss);
 
-                element->material->ViscosityMLHOAdjoint(&viscosity[0],dim,xyz_list,gauss,vxbase_input,vybase_input,vxshear_input,vyshear_input,thickness_input,n_input);
+                element->material->ViscosityMOLHOAdjoint(&viscosity[0],dim,xyz_list,gauss,vxbase_input,vybase_input,vxshear_input,vyshear_input,thickness_input,n_input);
 
                 effmu = 2.0*(1-n)/2.0/n;
 
@@ -352 +352 @@
    xDelete<IssmDouble>(basis);
    return Ke;
 }/*}}}*/
-ElementMatrix* AdjointHorizAnalysis::CreateKMatrixMLHOVerticalIntergrated(Element* element){/*{{{*/
+ElementMatrix* AdjointHorizAnalysis::CreateKMatrixMOLHOVerticalIntergrated(Element* element){/*{{{*/
 
         /* Check if ice in element */
         if(!element->IsIceInElement()) return NULL;
@@ -432 +432 @@
               zeta=0.5*(gauss_seg->coord1+1);
 
                    /* eps_eff^2 = exx^2 + eyy^2 + exy^2 + exz^2 + eyz^2 + exx*eyy (for a given zeta)*/
-                        element->StrainRateMLHO(&epsilon[0],xyz_list,gauss,
+                        element->StrainRateMOLHO(&epsilon[0],xyz_list,gauss,
                   vxbase_input,vybase_input,vxshear_input,vyshear_input,thickness_input,n_input,zeta);
                         epsilon_eff=sqrt(epsilon[0]*epsilon[0] + epsilon[1]*epsilon[1] + epsilon[2]*epsilon[2]
                     +  epsilon[3]*epsilon[3] + epsilon[4]*epsilon[4] + epsilon[0]*epsilon[1]);
@@ -617 +617 @@
                         return CreatePVectorHO(element);
                 case FSApproximationEnum: 
                         return CreatePVectorFS(element);
-      case MLHOApproximationEnum:
-                        return CreatePVectorMLHO(element);
+      case MOLHOApproximationEnum:
+                        return CreatePVectorMOLHO(element);
                 case NoneApproximationEnum:
                         return NULL;
                 default:
@@ -1104 +1104 @@
         return pe;
 
 }/*}}}*/
-ElementVector* AdjointHorizAnalysis::CreatePVectorMLHO(Element* element){/*{{{*/
+ElementVector* AdjointHorizAnalysis::CreatePVectorMOLHO(Element* element){/*{{{*/
 
         /*Intermediaries*/
         int      domaintype;
@@ -1142 +1142 @@
         int numnodes = basalelement->GetNumberOfNodes();
 
         /*Initialize Element vector and vectors*/
-        ElementVector* pe    = basalelement->NewElementVector(MLHOApproximationEnum);
+        ElementVector* pe    = basalelement->NewElementVector(MOLHOApproximationEnum);
         IssmDouble*    basis = xNew<IssmDouble>(numnodes);
 
         /*Retrieve all inputs and parameters*/
@@ -1211 +1211 @@
                                                         pe->values[i*4+3]+=duy*weight*Jdet*gauss->weight*basis[i]; 
                                                 }
                                                 else {
-                                                        _error_("2D vertical is not implemented for MLHO");
+                                                        _error_("2D vertical is not implemented for MOLHO");
                                                 }
                                         }
                                         break;
@@ -1239 +1239 @@
                                                         pe->values[i*4+3]+=duy*weight*Jdet*gauss->weight*basis[i];
                                                 }
                                                 else{
-                                                        _error_("2D vertical is not implemented for MLHO");
+                                                        _error_("2D vertical is not implemented for MOLHO");
                                                 }
                                         }
                                         break;
@@ -1268 +1268 @@
                                                         pe->values[i*4+3]+=duy*weight*Jdet*gauss->weight*basis[i];
                                                 }
                                                 else{
-                                                        _error_("2D vertical is not implemented for MLHO");
+                                                        _error_("2D vertical is not implemented for MOLHO");
                                                 }
                                         }
                                         break;
@@ -1294 +1294 @@
                                                         pe->values[i*4+3]+=duy*weight*Jdet*gauss->weight*basis[i];
                                                 }
                                                 else{
-                                                        _error_("2D vertical is not implemented for MLHO");
+                                                        _error_("2D vertical is not implemented for MOLHO");
                                                 }
                                         }
                                         break;
@@ -1318 +1318 @@
                                                         pe->values[i*4+3]+=duy*weight*Jdet*gauss->weight*basis[i];
                                                 }
                                                 else{
-                                                        _error_("2D vertical is not implemented for MLHO");
+                                                        _error_("2D vertical is not implemented for MOLHO");
                                                 }
                                         }
                                         break;
@@ -1668 +1668 @@
                                 case SSAApproximationEnum: GradientJDragSSA(element,gradient,control_interp,control_index); break;
                                 case L1L2ApproximationEnum:GradientJDragL1L2(element,gradient,control_interp,control_index); break;
                                 case HOApproximationEnum:  GradientJDragHO( element,gradient,control_interp,control_index); break;
-                                case MLHOApproximationEnum:  GradientJDragMLHO( element,gradient,control_interp,control_index); break;
+                                case MOLHOApproximationEnum:  GradientJDragMOLHO( element,gradient,control_interp,control_index); break;
                                 case FSApproximationEnum:  GradientJDragFS( element,gradient,control_interp,control_index); break;
                                 case NoneApproximationEnum: /*Gradient is 0*/                    break;
                                 default: _error_("approximation " << EnumToStringx(approximation) << " not supported yet");
@@ -1689 +1689 @@
                                 case SSAApproximationEnum: GradientJBbarSSA(element,gradient,control_interp,control_index); break;
                                 case L1L2ApproximationEnum:GradientJBbarL1L2(element,gradient,control_interp,control_index); break;
                                 case HOApproximationEnum:  GradientJBbarHO( element,gradient,control_interp,control_index); break;
-                                case MLHOApproximationEnum:  GradientJBbarMLHO( element,gradient,control_interp,control_index); break;
+                                case MOLHOApproximationEnum:  GradientJBbarMOLHO( element,gradient,control_interp,control_index); break;
                                 case FSApproximationEnum:  GradientJBbarFS( element,gradient,control_interp,control_index); break;
                                 case NoneApproximationEnum: /*Gradient is 0*/                    break;
                                 default: _error_("approximation " << EnumToStringx(approximation) << " not supported yet");
@@ -1699 +1699 @@
                         switch(approximation){
                                 case SSAApproximationEnum: GradientJBSSA(element,gradient,control_interp,control_index); break;
                                 case HOApproximationEnum:  GradientJBHO( element,gradient,control_interp,control_index); break;
-                        //      case MLHOApproximationEnum:  GradientJBMLHO( element,gradient,control_interp,control_index); break;
+                        //      case MOLHOApproximationEnum:  GradientJBMOLHO( element,gradient,control_interp,control_index); break;
                                 case FSApproximationEnum:  GradientJBFS( element,gradient,control_interp,control_index); break;
                                 case NoneApproximationEnum: /*Gradient is 0*/                    break;
                                 default: _error_("approximation " << EnumToStringx(approximation) << " not supported yet");
@@ -1809 +1809 @@
         /*WARNING: We use SSA as an estimate for now*/
         this->GradientJBbarSSA(element,gradient,control_interp,control_index);
 }/*}}}*/
-void           AdjointHorizAnalysis::GradientJBbarMLHO(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index){/*{{{*/
+void           AdjointHorizAnalysis::GradientJBbarMOLHO(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index){/*{{{*/
 
    if(control_interp!=P1Enum) _error_("not implemented yet...");
         /*Intermediaries*/
@@ -1892 +1892 @@
       while(gauss_seg->next()){
                         zeta=0.5*(gauss_seg->coord1+1);
 
-                        basalelement->StrainRateMLHO(&epsilon[0],xyz_list,gauss,
+                        basalelement->StrainRateMOLHO(&epsilon[0],xyz_list,gauss,
                   vxbase_input,vybase_input,vxshear_input,vyshear_input,thickness_input,n_input,zeta);
 
-                        basalelement->dViscositydBMLHO(&dmudB,dim,xyz_list,gauss,vxbase_input,vybase_input,vxshear_input,vyshear_input,thickness_input,n_input,zeta);
+                        basalelement->dViscositydBMOLHO(&dmudB,dim,xyz_list,gauss,vxbase_input,vybase_input,vxshear_input,vyshear_input,thickness_input,n_input,zeta);
 
                         e1[0] = 2.0*epsilon[0]+epsilon[1];
          e1[1] = epsilon[2];
@@ -2250 +2250 @@
         xDelete<int>(vertexpidlist);
         delete gauss;
 }/*}}}*/
-void           AdjointHorizAnalysis::GradientJBMLHO(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index){/*{{{*/
+void           AdjointHorizAnalysis::GradientJBMOLHO(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index){/*{{{*/
         _error_("not implemented yet...");
 }/*}}}*/
 void           AdjointHorizAnalysis::GradientJBSSA(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index){/*{{{*/
@@ -2650 +2650 @@
         delete gauss;
         delete friction;
 }/*}}}*/
-void           AdjointHorizAnalysis::GradientJDragMLHO(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index){/*{{{*/
+void           AdjointHorizAnalysis::GradientJDragMOLHO(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index){/*{{{*/
 
         /*return if floating (gradient is 0)*/
         if(element->IsAllFloating()) return;
@@ -3271 +3271 @@
         if(approximation==FSApproximationEnum || approximation==NoneApproximationEnum){
                 InputUpdateFromSolutionFS(solution,element);
         }
-        else if (approximation==MLHOApproximationEnum) {
-                InputUpdateFromSolutionMLHO(solution, element);
+        else if (approximation==MOLHOApproximationEnum) {
+                InputUpdateFromSolutionMOLHO(solution, element);
         }
         else{
                 InputUpdateFromSolutionHoriz(solution,element);
@@ -3419 +3419 @@
         xDelete<IssmDouble>(lambday);
         xDelete<int>(doflist);
 }/*}}}*/
-void           AdjointHorizAnalysis::InputUpdateFromSolutionMLHO(IssmDouble* solution,Element* element){/*{{{*/
+void           AdjointHorizAnalysis::InputUpdateFromSolutionMOLHO(IssmDouble* solution,Element* element){/*{{{*/
         int  i;
         int* doflist=NULL;
 
@@ -3432 +3432 @@
         int numdof = numnodes * 4;
 
         /*Fetch dof list and allocate solution vectors*/
-        element->GetDofListLocal(&doflist,MLHOApproximationEnum,GsetEnum);
+        element->GetDofListLocal(&doflist,MOLHOApproximationEnum,GsetEnum);
         IssmDouble* values    = xNew<IssmDouble>(numdof);
         IssmDouble* lambdax   = xNew<IssmDouble>(numnodes);
         IssmDouble* lambday   = xNew<IssmDouble>(numnodes);

../trunk-jpl/src/c/analyses/AdjointHorizAnalysis.h

@@ -27 +27 @@
                 ElementMatrix* CreateKMatrix(Element* element);
                 ElementMatrix* CreateKMatrixFS(Element* element);
                 ElementMatrix* CreateKMatrixHO(Element* element);
-                ElementMatrix* CreateKMatrixMLHO(Element* element);
-                ElementMatrix* CreateKMatrixMLHOVerticalIntergrated(Element* element);
+                ElementMatrix* CreateKMatrixMOLHO(Element* element);
+                ElementMatrix* CreateKMatrixMOLHOVerticalIntergrated(Element* element);
                 ElementMatrix* CreateKMatrixL1L2(Element* element);
                 ElementMatrix* CreateKMatrixSSA(Element* element);
                 ElementVector* CreatePVector(Element* element);
@@ -35 +35 @@
                 ElementVector* CreatePVectorFS(Element* element);
                 ElementVector* CreatePVectorL1L2(Element* element);
                 ElementVector* CreatePVectorHO(Element* element);
-                ElementVector* CreatePVectorMLHO(Element* element);
+                ElementVector* CreatePVectorMOLHO(Element* element);
                 ElementVector* CreatePVectorSSA(Element* element);
                 void           GetSolutionFromInputs(Vector<IssmDouble>* solution,Element* element);
                 void           GradientJ(Vector<IssmDouble>* gradient,Element*  element,int control_type,int control_interp,int control_index);
@@ -44 +44 @@
                 void           GradientJBinitial(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index);
                 void           GradientJBbarL1L2(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index);
                 void           GradientJBbarHO(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index);
-                void           GradientJBbarMLHO(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index);
+                void           GradientJBbarMOLHO(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index);
                 void           GradientJBbarSSA(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index);
                 void           GradientJBFS(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index);
                 void           GradientJBGradient(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index);
                 void           GradientJBHO(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index);
-                void           GradientJBMLHO(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index);
+                void           GradientJBMOLHO(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index);
                 void           GradientJBSSA(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index);
                 void           GradientJDragFS(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index);
                 void           GradientJDragGradient(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index);
                 void           GradientJDragL1L2(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index);
                 void           GradientJDragHO(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index);
-                void           GradientJDragMLHO(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index);
+                void           GradientJDragMOLHO(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index);
                 void           GradientJDragSSA(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index);
                 void           GradientJDragHydroFS(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index);
                 void           GradientJDragHydroL1L2(Element* element,Vector<IssmDouble>* gradient,int control_interp,int control_index);
@@ -65 +65 @@
                 void           InputUpdateFromSolution(IssmDouble* solution,Element* element);
                 void           InputUpdateFromSolutionFS(IssmDouble* solution,Element* element);
                 void           InputUpdateFromSolutionHoriz(IssmDouble* solution,Element* element);
-                void           InputUpdateFromSolutionMLHO(IssmDouble* solution,Element* element);
+                void           InputUpdateFromSolutionMOLHO(IssmDouble* solution,Element* element);
                 void           UpdateConstraints(FemModel* femmodel);
 };
 #endif

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

@@ -30 +30 @@
         IssmDouble g;
         IssmDouble rho_ice;
         IssmDouble FSreconditioning;
-        bool       isSIA,isSSA,isL1L2,isMLHO,isHO,isFS,iscoupling;
+        bool       isSIA,isSSA,isL1L2,isMOLHO,isHO,isFS,iscoupling;
         bool       spcpresent = false;
         int        Mx,Nx;
         int        My,Ny;
@@ -58 +58 @@
         iomodel->FindConstant(&isSIA,"md.flowequation.isSIA");
         iomodel->FindConstant(&isSSA,"md.flowequation.isSSA");
         iomodel->FindConstant(&isL1L2,"md.flowequation.isL1L2");
-        iomodel->FindConstant(&isMLHO,"md.flowequation.isMLHO");
+        iomodel->FindConstant(&isMOLHO,"md.flowequation.isMOLHO");
         iomodel->FindConstant(&isHO,"md.flowequation.isHO");
         iomodel->FindConstant(&isFS,"md.flowequation.isFS");
 
         /*Is this model only SIA??*/
-        if(!isSSA && !isHO && !isFS && !isL1L2 && !isMLHO) return;
+        if(!isSSA && !isHO && !isFS && !isL1L2 && !isMOLHO) return;
 
         /*Do we have coupling*/
-        if((isSIA?1.:0.) + (isSSA?1.:0.) + (isL1L2?1.:0.) + (isMLHO?1.:0.) + (isHO?1.:0.) + (isFS?1.:0.) >1.)
+        if((isSIA?1.:0.) + (isSSA?1.:0.) + (isL1L2?1.:0.) + (isMOLHO?1.:0.) + (isHO?1.:0.) + (isFS?1.:0.) >1.)
          iscoupling = true;
         else
          iscoupling = false;
@@ -77 +77 @@
                 /*Get finite element type*/
                 if(isSSA)       iomodel->FindConstant(&finiteelement,"md.flowequation.fe_SSA");
                 else if(isL1L2) finiteelement = P1Enum;
-                else if(isMLHO) finiteelement = P1Enum;
+                else if(isMOLHO) finiteelement = P1Enum;
                 else if(isHO)   iomodel->FindConstant(&finiteelement,"md.flowequation.fe_HO");
                 else if(isFS){  iomodel->FindConstant(&finiteelement,"md.flowequation.fe_FS");
                         /*Deduce velocity interpolation from finite element*/
@@ -187 +187 @@
                         iomodel->DeleteData(z,"md.mesh.z");
                 }
                 else{
-                        if(!isMLHO){
+                        if(!isMOLHO){
                                 IoModelToConstraintsx(constraints,iomodel,"md.stressbalance.spcvx",StressbalanceAnalysisEnum,finiteelement,0);
                                 if(iomodel->domaintype!=Domain2DverticalEnum){
                                         IoModelToConstraintsx(constraints,iomodel,"md.stressbalance.spcvy",StressbalanceAnalysisEnum,finiteelement,1);
                                 }
                         }
-                        else{//MLHO 
+                        else{//MOLHO 
                                 IoModelToConstraintsx(constraints,iomodel,"md.stressbalance.spcvx_base",StressbalanceAnalysisEnum,finiteelement,0);
                                 IoModelToConstraintsx(constraints,iomodel,"md.stressbalance.spcvx_shear",StressbalanceAnalysisEnum,finiteelement,1);
                                 IoModelToConstraintsx(constraints,iomodel,"md.stressbalance.spcvy_base",StressbalanceAnalysisEnum,finiteelement,2);
@@ -461 +461 @@
         int         i;
         int         count;
         int         penpair_ids[2];
-        bool        isSSA,isL1L2,isMLHO,isHO,isFS;
+        bool        isSSA,isL1L2,isMOLHO,isHO,isFS;
         int         numpenalties,numrifts,numriftsegments;
         IssmDouble *riftinfo       = NULL;
         IssmDouble *penalties      = NULL;
@@ -469 +469 @@
 
         /*Fetch parameters: */
         iomodel->FindConstant(&isL1L2,"md.flowequation.isL1L2");
-        iomodel->FindConstant(&isMLHO,"md.flowequation.isMLHO");
+        iomodel->FindConstant(&isMOLHO,"md.flowequation.isMOLHO");
         iomodel->FindConstant(&isFS,"md.flowequation.isFS");
         iomodel->FindConstant(&isSSA,"md.flowequation.isSSA");
         iomodel->FindConstant(&isHO,"md.flowequation.isHO");
@@ -476 +476 @@
         iomodel->FindConstant(&numrifts,"md.rifts.numrifts");
 
         /*Is this SIA only?*/
-        if(!isSSA && !isHO && !isFS && !isL1L2 && !isMLHO) return;
+        if(!isSSA && !isHO && !isFS && !isL1L2 && !isMOLHO) return;
 
         /*Initialize counter: */
         count=0;
@@ -521 +521 @@
 void StressbalanceAnalysis::CreateNodes(Nodes* nodes,IoModel* iomodel,bool isamr){/*{{{*/
 
         /*Intermediary*/
-        bool isSSA,isL1L2,isMLHO,isHO,isFS,iscoupling;
+        bool isSSA,isL1L2,isMOLHO,isHO,isFS,iscoupling;
         int  finiteelement=-1,approximation=-1;
 
         /*Fetch parameters: */
         iomodel->FindConstant(&isSSA,"md.flowequation.isSSA");
         iomodel->FindConstant(&isL1L2,"md.flowequation.isL1L2");
-        iomodel->FindConstant(&isMLHO,"md.flowequation.isMLHO");
+        iomodel->FindConstant(&isMOLHO,"md.flowequation.isMOLHO");
         iomodel->FindConstant(&isHO,"md.flowequation.isHO");
         iomodel->FindConstant(&isFS,"md.flowequation.isFS");
 
         /*Now, check that we have non SIA elements */
-        if(!isSSA && !isL1L2 && !isMLHO && !isHO && !isFS) return;
+        if(!isSSA && !isL1L2 && !isMOLHO && !isHO && !isFS) return;
 
         /*Do we have coupling*/
-        if( (isSSA?1.:0.) + (isL1L2?1.:0.) + (isMLHO?1.:0.) + (isHO?1.:0.) + (isFS?1.:0.) >1.)
+        if( (isSSA?1.:0.) + (isL1L2?1.:0.) + (isMOLHO?1.:0.) + (isHO?1.:0.) + (isFS?1.:0.) >1.)
          iscoupling = true;
         else
          iscoupling = false;
@@ -552 +552 @@
                         approximation = L1L2ApproximationEnum;
                         finiteelement = P1Enum;
                 }
-                else if(isMLHO){
-                        approximation = MLHOApproximationEnum;
+                else if(isMOLHO){
+                        approximation = MOLHOApproximationEnum;
                         finiteelement = P1Enum;
                 }
                 else if(isHO){
@@ -629 +629 @@
                          }
                          break;
                 case L1L2ApproximationEnum: numdofs = 2; break;
-                case MLHOApproximationEnum: numdofs = 4; break;
+                case MOLHOApproximationEnum: numdofs = 4; break;
                 case HOApproximationEnum:
                          switch(domaintype){
                                  case Domain3DEnum:         numdofs=2; break;
@@ -694 +694 @@
         int    approximation,frictionlaw;
         int    FrictionCoupling;
         int*   finiteelement_list=NULL;
-        bool   isSSA,isL1L2,isMLHO,isHO,isFS,iscoupling;
+        bool   isSSA,isL1L2,isMOLHO,isHO,isFS,iscoupling;
         bool   control_analysis;
         bool   dakota_analysis;
         bool   ismovingfront;
@@ -702 +702 @@
         /*Fetch constants needed: */
         iomodel->FindConstant(&isSSA,"md.flowequation.isSSA");
         iomodel->FindConstant(&isL1L2,"md.flowequation.isL1L2");
-        iomodel->FindConstant(&isMLHO,"md.flowequation.isMLHO");
+        iomodel->FindConstant(&isMOLHO,"md.flowequation.isMOLHO");
         iomodel->FindConstant(&isHO,"md.flowequation.isHO");
         iomodel->FindConstant(&isFS,"md.flowequation.isFS");
         iomodel->FindConstant(&control_analysis,"md.inversion.iscontrol");
@@ -712 +712 @@
         iomodel->FindConstant(&frictionlaw,"md.friction.law");
 
         /*return if no processing required*/
-        if(!isSSA && !isL1L2 && !isMLHO && !isHO && !isFS) return;
+        if(!isSSA && !isL1L2 && !isMOLHO && !isHO && !isFS) return;
 
         /*Fetch data needed and allocate vectors: */
         iomodel->FetchData(1,"md.flowequation.element_equation");
@@ -719 +719 @@
         finiteelement_list=xNewZeroInit<int>(iomodel->numberofelements);
 
         /*Do we have coupling*/
-        if( (isSSA?1.:0.) + (isL1L2?1.:0.) + (isMLHO?1.:0.) + (isHO?1.:0.) + (isFS?1.:0.) >1.)
+        if( (isSSA?1.:0.) + (isL1L2?1.:0.) + (isMOLHO?1.:0.) + (isHO?1.:0.) + (isFS?1.:0.) >1.)
          iscoupling = true;
         else
          iscoupling = false;
@@ -728 +728 @@
         if(!iscoupling){
                 if(isSSA)       iomodel->FindConstant(&finiteelement,"md.flowequation.fe_SSA");
                 else if(isL1L2) finiteelement = P1Enum;
-                else if(isMLHO) finiteelement = P1Enum;
+                else if(isMOLHO) finiteelement = P1Enum;
                 else if(isHO)   iomodel->FindConstant(&finiteelement,"md.flowequation.fe_HO");
                 else if(isFS)   iomodel->FindConstant(&finiteelement,"md.flowequation.fe_FS");
                 for(int i=0;i<iomodel->numberofelements;i++){
@@ -780 +780 @@
         iomodel->FetchDataToInput(inputs,elements,"md.mask.ocean_levelset",MaskOceanLevelsetEnum);
         iomodel->FetchDataToInput(inputs,elements,"md.initialization.vx",VxEnum,0.);
         iomodel->FetchDataToInput(inputs,elements,"md.initialization.vy",VyEnum,0.);
-        /*MLHO*/
-        if(isMLHO){
+        /*MOLHO*/
+        if(isMOLHO){
                 iomodel->FetchDataToInput(inputs,elements,"md.initialization.vx",VxShearEnum,0.);
                 iomodel->FetchDataToInput(inputs,elements,"md.initialization.vy",VyShearEnum,0.);
-                /*3D MLHO also need to have VxBase and VyBase for reconstruting Vx and Vy*/
+                /*3D MOLHO also need to have VxBase and VyBase for reconstruting Vx and Vy*/
                 if (iomodel->domaintype==Domain3DEnum) {
                         iomodel->FetchDataToInput(inputs,elements,"md.initialization.vx",VxBaseEnum,0.);
                         iomodel->FetchDataToInput(inputs,elements,"md.initialization.vy",VyBaseEnum,0.);
@@ -989 +989 @@
         parameters->AddObject(iomodel->CopyConstantObject("md.flowequation.isSIA",FlowequationIsSIAEnum));
         parameters->AddObject(iomodel->CopyConstantObject("md.flowequation.isSSA",FlowequationIsSSAEnum));
         parameters->AddObject(iomodel->CopyConstantObject("md.flowequation.isL1L2",FlowequationIsL1L2Enum));
-        parameters->AddObject(iomodel->CopyConstantObject("md.flowequation.isMLHO",FlowequationIsMLHOEnum));
+        parameters->AddObject(iomodel->CopyConstantObject("md.flowequation.isMOLHO",FlowequationIsMOLHOEnum));
         parameters->AddObject(iomodel->CopyConstantObject("md.flowequation.isHO",FlowequationIsHOEnum));
         parameters->AddObject(iomodel->CopyConstantObject("md.flowequation.isFS",FlowequationIsFSEnum));
         parameters->AddObject(iomodel->CopyConstantObject("md.flowequation.fe_FS",FlowequationFeFSEnum));
@@ -1092 +1092 @@
 void           StressbalanceAnalysis::Core(FemModel* femmodel){/*{{{*/
 
         /*Intermediaries*/
-        bool isSSA,isL1L2,isMLHO,isHO,isFS;
+        bool isSSA,isL1L2,isMOLHO,isHO,isFS;
         bool conserve_loads = true;
         int  newton,domaintype,fe_FS;
 
@@ -1099 +1099 @@
         /* recover parameters:*/
         femmodel->parameters->FindParam(&isSSA,FlowequationIsSSAEnum);
         femmodel->parameters->FindParam(&isL1L2,FlowequationIsL1L2Enum);
-        femmodel->parameters->FindParam(&isMLHO,FlowequationIsMLHOEnum);
+        femmodel->parameters->FindParam(&isMOLHO,FlowequationIsMOLHOEnum);
         femmodel->parameters->FindParam(&isHO,FlowequationIsHOEnum);
         femmodel->parameters->FindParam(&isFS,FlowequationIsFSEnum);
         femmodel->parameters->FindParam(&fe_FS,FlowequationFeFSEnum);
@@ -1106 +1106 @@
         femmodel->parameters->FindParam(&domaintype,DomainTypeEnum);
         femmodel->parameters->FindParam(&newton,StressbalanceIsnewtonEnum);
 
-        if(isFS && !(isSSA || isHO || isL1L2 || isMLHO)){
+        if(isFS && !(isSSA || isHO || isL1L2 || isMOLHO)){
                 femmodel->SetCurrentConfiguration(StressbalanceAnalysisEnum);
 
                 bool is_schur_cg_solver = false;
@@ -1128 +1128 @@
                 else
                  solutionsequence_nonlinear(femmodel,conserve_loads);
         }
-        else if(!isFS && (isSSA || isHO || isL1L2 || isMLHO)){
+        else if(!isFS && (isSSA || isHO || isL1L2 || isMOLHO)){
                 femmodel->SetCurrentConfiguration(StressbalanceAnalysisEnum);
                 if(newton>0)
                  solutionsequence_newton(femmodel);
@@ -1142 +1142 @@
                         extrudefrombase_core(femmodel);
                 }
         }
-        else if ((isSSA || isL1L2 || isMLHO || isHO) && isFS){
+        else if ((isSSA || isL1L2 || isMOLHO || isHO) && isFS){
                 if(VerboseSolution()) _printf0_("   computing coupling between lower order models and FS\n");
                 solutionsequence_FScoupling_nonlinear(femmodel,conserve_loads);
         }
@@ -1194 +1194 @@
                         return CreateKMatrixSSA(element);
                 case L1L2ApproximationEnum:
                         return CreateKMatrixL1L2(element);
-                case MLHOApproximationEnum:
-                        return CreateKMatrixMLHO(element);
+                case MOLHOApproximationEnum:
+                        return CreateKMatrixMOLHO(element);
                 case HOApproximationEnum:
                         return CreateKMatrixHO(element);
                 case FSApproximationEnum:
@@ -1223 +1223 @@
                         return CreatePVectorSSA(element);
                 case L1L2ApproximationEnum:
                         return CreatePVectorL1L2(element);
-                case MLHOApproximationEnum:
-                        return CreatePVectorMLHO(element);
+                case MOLHOApproximationEnum:
+                        return CreatePVectorMOLHO(element);
                 case HOApproximationEnum:
                         return CreatePVectorHO(element);
                 case FSApproximationEnum:
@@ -1252 +1252 @@
                 case SSAApproximationEnum: case HOApproximationEnum: case L1L2ApproximationEnum: case SIAApproximationEnum:
                         GetSolutionFromInputsHoriz(solution,element);
                         return;
-                case MLHOApproximationEnum:
-                        GetSolutionFromInputsMLHO(solution,element);
+                case MOLHOApproximationEnum:
+                        GetSolutionFromInputsMOLHO(solution,element);
                         return;
                 case SSAHOApproximationEnum: case HOFSApproximationEnum: case SSAFSApproximationEnum:
                         /*the elements around will create the solution*/
@@ -1334 +1334 @@
                 case L1L2ApproximationEnum:
                         InputUpdateFromSolutionL1L2(solution,element);
                         return;
-                case MLHOApproximationEnum:
-                        InputUpdateFromSolutionMLHO(solution,element);
+                case MOLHOApproximationEnum:
+                        InputUpdateFromSolutionMOLHO(solution,element);
                         return;
                 case SSAHOApproximationEnum:
                         InputUpdateFromSolutionSSAHO(solution,element);
@@ -2759 +2759 @@
         if(basalelement->IsSpawnedElement()){basalelement->DeleteMaterials(); delete basalelement;};
 }/*}}}*/
 
-/*MLHO*/
-ElementMatrix* StressbalanceAnalysis::CreateKMatrixMLHO(Element* element){/*{{{*/
+/*MOLHO*/
+ElementMatrix* StressbalanceAnalysis::CreateKMatrixMOLHO(Element* element){/*{{{*/
 
         /* Check if ice in element */
         if(!element->IsIceInElement()) return NULL;
@@ -2783 +2783 @@
         }
 
         /*compute all stiffness matrices for this element*/
-        ElementMatrix* Ke1=CreateKMatrixMLHOViscous(basalelement);
-        ElementMatrix* Ke2=CreateKMatrixMLHOFriction(basalelement);
+        ElementMatrix* Ke1=CreateKMatrixMOLHOViscous(basalelement);
+        ElementMatrix* Ke2=CreateKMatrixMOLHOFriction(basalelement);
         ElementMatrix* Ke =new ElementMatrix(Ke1,Ke2);
 
         /*clean-up and return*/
@@ -2793 +2793 @@
         delete Ke2;
         return Ke;
 }/*}}}*/
-ElementMatrix* StressbalanceAnalysis::CreateKMatrixMLHOFriction(Element* element){/*{{{*/
+ElementMatrix* StressbalanceAnalysis::CreateKMatrixMOLHOFriction(Element* element){/*{{{*/
 
         if(element->IsAllFloating()) return NULL;
 
@@ -2820 +2820 @@
         int numdof   = numnodes*dim;
 
         /*Initialize Element matrix and vectors*/
-        ElementMatrix* Ke = element->NewElementMatrix(MLHOApproximationEnum);
+        ElementMatrix* Ke = element->NewElementMatrix(MOLHOApproximationEnum);
         IssmDouble*    basis  = xNew<IssmDouble>(numnodes);
 
         /*Retrieve all inputs and parameters*/
@@ -2875 +2875 @@
         xDelete<IssmDouble>(basis);
         return Ke;
 }/*}}}*/
-ElementMatrix* StressbalanceAnalysis::CreateKMatrixMLHOViscous(Element* element){/*{{{*/
+ElementMatrix* StressbalanceAnalysis::CreateKMatrixMOLHOViscous(Element* element){/*{{{*/
 
         /* Check if ice in element */
         if(!element->IsIceInElement()) return NULL;
@@ -2890 +2890 @@
         int numnodes = element->GetNumberOfNodes();
 
         /*Initialize Element matrix and vectors*/
-        ElementMatrix* Ke     = element->NewElementMatrix(MLHOApproximationEnum);
+        ElementMatrix* Ke     = element->NewElementMatrix(MOLHOApproximationEnum);
         IssmDouble*    dbasis = xNew<IssmDouble>(2*numnodes); // like SSA
         IssmDouble*    basis  = xNew<IssmDouble>(numnodes); // like SSA
 
@@ -2914 +2914 @@
 
                 thickness_input->GetInputValue(&thickness, gauss);
                 n_input->GetInputValue(&n,gauss);
-                element->material->ViscosityMLHO(&viscosity[0],dim,xyz_list,gauss,vxbase_input,vybase_input,vxshear_input,vyshear_input,thickness_input,n_input);
+                element->material->ViscosityMOLHO(&viscosity[0],dim,xyz_list,gauss,vxbase_input,vybase_input,vxshear_input,vyshear_input,thickness_input,n_input);
 
                 for(int i=0;i<numnodes;i++){//shape functions on tria element
                         for(int j=0;j<numnodes;j++){
@@ -2980 +2980 @@
         }
 
         /*Transform Coordinate System*/
-        //basalelement->TransformStiffnessMatrixCoord(Ke,XYMLHOEnum);
+        //basalelement->TransformStiffnessMatrixCoord(Ke,XYMOLHOEnum);
 
         /*Clean up and return*/
         delete gauss;
@@ -2989 +2989 @@
         xDelete<IssmDouble>(basis);
         return Ke;
 }/*}}}*/
-ElementVector* StressbalanceAnalysis::CreatePVectorMLHO(Element* element){/*{{{*/
+ElementVector* StressbalanceAnalysis::CreatePVectorMOLHO(Element* element){/*{{{*/
 
         /* Check if ice in element */
         if(!element->IsIceInElement()) return NULL;
@@ -3012 +3012 @@
         }
 
         /*compute all load vectors for this element*/
-        ElementVector* pe1=CreatePVectorMLHODrivingStress(basalelement);
-        ElementVector* pe2=CreatePVectorMLHOFront(basalelement);
+        ElementVector* pe1=CreatePVectorMOLHODrivingStress(basalelement);
+        ElementVector* pe2=CreatePVectorMOLHOFront(basalelement);
         ElementVector* pe =new ElementVector(pe1,pe2);
 
         /*clean-up and return*/
@@ -3022 +3022 @@
         delete pe2;
         return pe;
 }/*}}}*/
-ElementVector* StressbalanceAnalysis::CreatePVectorMLHODrivingStress(Element* element){/*{{{*/
+ElementVector* StressbalanceAnalysis::CreatePVectorMOLHODrivingStress(Element* element){/*{{{*/
 
         /*Intermediaries */
         IssmDouble  thickness,Jdet,slope[2],n;
@@ -3032 +3032 @@
         int numnodes = element->GetNumberOfNodes();
 
         /*Initialize Element vector and vectors*/
-        ElementVector* pe    = element->NewElementVector(MLHOApproximationEnum);
+        ElementVector* pe    = element->NewElementVector(MOLHOApproximationEnum);
         IssmDouble*    basis = xNew<IssmDouble>(numnodes);
 
         /*Retrieve all inputs and parameters*/
@@ -3061 +3061 @@
         }
 
         /*Transform coordinate system*/
-        //element->TransformLoadVectorCoord(pe,XYMLHOEnum);
+        //element->TransformLoadVectorCoord(pe,XYMOLHOEnum);
 
         /*Clean up and return*/
         xDelete<IssmDouble>(xyz_list);
@@ -3069 +3069 @@
         delete gauss;
         return pe;
 }/*}}}*/
-ElementVector* StressbalanceAnalysis::CreatePVectorMLHOFront(Element* element){/*{{{*/
+ElementVector* StressbalanceAnalysis::CreatePVectorMOLHOFront(Element* element){/*{{{*/
 
         /*If no front, return NULL*/
         if(!element->IsIcefront()) return NULL;
@@ -3085 +3085 @@
         int numnodes = element->GetNumberOfNodes();
 
         /*Initialize Element vector and other vectors*/
-        ElementVector* pe    = element->NewElementVector(MLHOApproximationEnum);
+        ElementVector* pe    = element->NewElementVector(MOLHOApproximationEnum);
         IssmDouble*    basis = xNew<IssmDouble>(numnodes);
 
         /*Retrieve all inputs and parameters*/
@@ -3134 +3134 @@
         }
 
         /*Transform coordinate system*/
-        //element->TransformLoadVectorCoord(pe,XYMLHOEnum);
+        //element->TransformLoadVectorCoord(pe,XYMOLHOEnum);
 
         /*Clean up and return*/
         xDelete<IssmDouble>(xyz_list);
@@ -3143 +3143 @@
         delete gauss;
         return pe;
 }/*}}}*/
-void           StressbalanceAnalysis::InputUpdateFromSolutionMLHO(IssmDouble* solution,Element* element){/*{{{*/
+void           StressbalanceAnalysis::InputUpdateFromSolutionMOLHO(IssmDouble* solution,Element* element){/*{{{*/
 
         int         i,dim,domaintype;
         IssmDouble  rho_ice,g;
@@ -3196 +3196 @@
         int numdof   = numnodes*dim*2; //2xdim DOFs per node
 
         /*Fetch dof list and allocate solution vectors*/
-        basalelement->GetDofListLocal(&doflist,MLHOApproximationEnum,GsetEnum); 
+        basalelement->GetDofListLocal(&doflist,MOLHOApproximationEnum,GsetEnum); 
         IssmDouble* values    = xNew<IssmDouble>(numdof);
         IssmDouble* vbx       = xNew<IssmDouble>(numnodes);
         IssmDouble* vby       = xNew<IssmDouble>(numnodes);
@@ -3219 +3219 @@
         if(dim==2) basalelement->TransformSolutionCoord(&values[0],XYEnum);
 
    /*Ok, we have vx and vy in values, fill in vx and vy arrays: */
-   for(i=0;i<numnodes;i++){ //numnodes of the 2D mesh in which the MLHO is written
+   for(i=0;i<numnodes;i++){ //numnodes of the 2D mesh in which the MOLHO is written
       vbx[i] =values[i*4+0]; //base vx
       vshx[i]=values[i*4+1]; //shear vx
                 vsx[i] =vbx[i]+vshx[i]; //surface vx
@@ -3256 +3256 @@
 
    switch(domaintype){
       case Domain2DhorizontalEnum:
-                        for(i=0;i<numnodes;i++){ //numnodes of the 2D mesh in which the MLHO is written
+                        for(i=0;i<numnodes;i++){ //numnodes of the 2D mesh in which the MOLHO is written
                                 vx[i]=vbx[i]+vshx[i]*(n[i]+1)/(n[i]+2);
                                 vy[i]=vby[i]+vshy[i]*(n[i]+1)/(n[i]+2);
                                 vel[i]=sqrt(vx[i]*vx[i] + vy[i]*vy[i]); 
@@ -3268 +3268 @@
       case Domain3DEnum:
                    basalelement->GetInputListOnNodes(&H[0],ThicknessEnum,0.);
                    basalelement->GetInputListOnNodes(&s[0],SurfaceEnum,0.);
-                        element->Recover3DMLHOInput(VxEnum, numnodes, vbx, vshx, n, H, s);
-                        element->Recover3DMLHOInput(VyEnum, numnodes, vby, vshy, n, H, s);
+                        element->Recover3DMOLHOInput(VxEnum, numnodes, vbx, vshx, n, H, s);
+                        element->Recover3DMOLHOInput(VyEnum, numnodes, vby, vshy, n, H, s);
                         break;
                 default: _error_("mesh "<<EnumToStringx(domaintype)<<" not supported yet");
         }
@@ -3293 +3293 @@
         xDelete<int>(doflist);
         if(basalelement->IsSpawnedElement()){basalelement->DeleteMaterials(); delete basalelement;};
 }/*}}}*/
-void           StressbalanceAnalysis::GetSolutionFromInputsMLHO(Vector<IssmDouble>* solution,Element* element){/*{{{*/
+void           StressbalanceAnalysis::GetSolutionFromInputsMOLHO(Vector<IssmDouble>* solution,Element* element){/*{{{*/
 
         IssmDouble   vbx,vby,vshx,vshy;
         int          domaintype,dim,approximation,dofpernode;
@@ -3312 +3312 @@
         int numnodes = element->GetNumberOfNodes();
         int numdof   = numnodes*dofpernode;
         element->GetInputValue(&approximation,ApproximationEnum);
-        if(approximation!=MLHOApproximationEnum) _error_("mesh "<<EnumToStringx(approximation)<<" not supported here");
+        if(approximation!=MOLHOApproximationEnum) _error_("mesh "<<EnumToStringx(approximation)<<" not supported here");
 
         /*Fetch dof list and allocate solution vector*/
         element->GetDofList(&doflist,approximation,GsetEnum); 

../trunk-jpl/src/c/analyses/StressbalanceAnalysis.h

@@ -56 +56 @@
                 ElementVector* CreatePVectorL1L2Front(Element* element);
                 ElementVector* CreatePVectorL1L2DrivingStress(Element* element);
                 void           InputUpdateFromSolutionL1L2(IssmDouble* solution,Element* element);
-                /*MLHO*/
-                ElementMatrix* CreateKMatrixMLHO(Element* element);
-                ElementMatrix* CreateKMatrixMLHOFriction(Element* element);
-                ElementMatrix* CreateKMatrixMLHOViscous(Element* element);
-                ElementVector* CreatePVectorMLHO(Element* element);
-                ElementVector* CreatePVectorMLHOFront(Element* element);
-                ElementVector* CreatePVectorMLHODrivingStress(Element* element);
-                void           InputUpdateFromSolutionMLHO(IssmDouble* solution,Element* element);
-                void           GetSolutionFromInputsMLHO(Vector<IssmDouble>* solution,Element* element);
+                /*MOLHO*/
+                ElementMatrix* CreateKMatrixMOLHO(Element* element);
+                ElementMatrix* CreateKMatrixMOLHOFriction(Element* element);
+                ElementMatrix* CreateKMatrixMOLHOViscous(Element* element);
+                ElementVector* CreatePVectorMOLHO(Element* element);
+                ElementVector* CreatePVectorMOLHOFront(Element* element);
+                ElementVector* CreatePVectorMOLHODrivingStress(Element* element);
+                void           InputUpdateFromSolutionMOLHO(IssmDouble* solution,Element* element);
+                void           GetSolutionFromInputsMOLHO(Vector<IssmDouble>* solution,Element* element);
                 /*HO*/
                 ElementMatrix* CreateJacobianMatrixHO(Element* element);
                 ElementMatrix* CreateKMatrixHO(Element* element);

../trunk-jpl/src/c/analyses/StressbalanceSIAAnalysis.cpp

@@ -9 +9 @@
 void StressbalanceSIAAnalysis::CreateConstraints(Constraints* constraints,IoModel* iomodel){/*{{{*/
 
         /*Intermediaries*/
-        bool       isSIA,isSSA,isL1L2,isMLHO,isHO,isFS,iscoupling;
+        bool       isSIA,isSSA,isL1L2,isMOLHO,isHO,isFS,iscoupling;
 
         /*Fetch parameters: */
         iomodel->FindConstant(&isSIA,"md.flowequation.isSIA");
         iomodel->FindConstant(&isSSA,"md.flowequation.isSSA");
         iomodel->FindConstant(&isL1L2,"md.flowequation.isL1L2");
-        iomodel->FindConstant(&isMLHO,"md.flowequation.isMLHO");
+        iomodel->FindConstant(&isMOLHO,"md.flowequation.isMOLHO");
         iomodel->FindConstant(&isHO,"md.flowequation.isHO");
         iomodel->FindConstant(&isFS,"md.flowequation.isFS");
 
@@ -23 +23 @@
         if (!isSIA) return;
 
         /*Do we have coupling*/
-        if((isSIA?1.:0.) + (isSSA?1.:0.) + (isL1L2?1.:0.) + (isMLHO?1.:0.) + (isHO?1.:0.) + (isFS?1.:0.) >1.)
+        if((isSIA?1.:0.) + (isSSA?1.:0.) + (isL1L2?1.:0.) + (isMOLHO?1.:0.) + (isHO?1.:0.) + (isFS?1.:0.) >1.)
          iscoupling = true;
         else
          iscoupling = false;

../trunk-jpl/src/c/analyses/StressbalanceVerticalAnalysis.cpp

@@ -11 +11 @@
 void StressbalanceVerticalAnalysis::CreateConstraints(Constraints* constraints,IoModel* iomodel){/*{{{*/
 
         /*Intermediary*/
-        bool        isSIA,isSSA,isL1L2,isMLHO,isHO,isFS,iscoupling;
+        bool        isSIA,isSSA,isL1L2,isMOLHO,isHO,isFS,iscoupling;
         int         Mz,Nz;
         IssmDouble *spcvz = NULL;
 
@@ -22 +22 @@
         iomodel->FindConstant(&isSIA,"md.flowequation.isSIA");
         iomodel->FindConstant(&isSSA,"md.flowequation.isSSA");
         iomodel->FindConstant(&isL1L2,"md.flowequation.isL1L2");
-        iomodel->FindConstant(&isMLHO,"md.flowequation.isMLHO");
+        iomodel->FindConstant(&isMOLHO,"md.flowequation.isMOLHO");
         iomodel->FindConstant(&isHO,"md.flowequation.isHO");
         iomodel->FindConstant(&isFS,"md.flowequation.isFS");
 
         /*Do we have coupling*/
-        if((isSIA?1.:0.) + (isSSA?1.:0.) + (isL1L2?1.:0.) + (isMLHO?1.:0.) + (isHO?1.:0.) + (isFS?1.:0.) >1.)
+        if((isSIA?1.:0.) + (isSSA?1.:0.) + (isL1L2?1.:0.) + (isMOLHO?1.:0.) + (isHO?1.:0.) + (isFS?1.:0.) >1.)
          iscoupling = true;
         else
          iscoupling = false;

../trunk-jpl/src/c/classes/Elements/Element.cpp

@@ -1040 +1040 @@
 
 }
 /*}}}*/
-void       Element::dViscositydBMLHO(IssmDouble* pdmudB,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vxbase_input,Input* vybase_input, Input* vxshear_input ,Input* vyshear_input,Input* thickness_input,Input* n_input, IssmDouble zeta){/*{{{*/
+void       Element::dViscositydBMOLHO(IssmDouble* pdmudB,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vxbase_input,Input* vybase_input, Input* vxshear_input ,Input* vyshear_input,Input* thickness_input,Input* n_input, IssmDouble zeta){/*{{{*/
 
         /*Intermediaries*/
         int materialstype;
@@ -1050 +1050 @@
         IssmDouble eps0=1.e-27;
 
         /* eps_eff^2 = exx^2 + eyy^2 + exy^2 + exz^2 + eyz^2 + exx*eyy */
-   this->StrainRateMLHO(&epsilon[0],xyz_list,gauss,
+   this->StrainRateMOLHO(&epsilon[0],xyz_list,gauss,
                                 vxbase_input,vybase_input,vxshear_input,vyshear_input,thickness_input,n_input,zeta);
    eps_eff=sqrt(epsilon[0]*epsilon[0] + epsilon[1]*epsilon[1] + epsilon[2]*epsilon[2]
                                 + epsilon[3]*epsilon[3] + epsilon[4]*epsilon[4] + epsilon[0]*epsilon[1] + eps0*eps0);
@@ -4618 +4618 @@
         epsilon[1] = 0.5*dvx[1];
 
 }/*}}}*/
-void       Element::StrainRateMLHO(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vxbase_input,Input* vybase_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input,IssmDouble zeta){/*{{{*/
-        /*Compute the 2d Blatter/MLHO Strain Rate (5 components) for a given vertical coordinate (zeta):
+void       Element::StrainRateMOLHO(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vxbase_input,Input* vybase_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input,IssmDouble zeta){/*{{{*/
+        /*Compute the 2d Blatter/MOLHO Strain Rate (5 components) for a given vertical coordinate (zeta):
          *
          * epsilon=[exx eyy exy exz eyz]
          *

../trunk-jpl/src/c/classes/Elements/Element.h

@@ -85 +85 @@
                 void               dViscositydBFS(IssmDouble* pdmudB,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input);
                 void               dViscositydBHO(IssmDouble* pdmudB,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input);
                 void               dViscositydBSSA(IssmDouble* pdmudB,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input);
-                void               dViscositydBMLHO(IssmDouble* pdmudB,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vxbase_input,Input* vybase_input, Input* vxshear_input ,Input* vyshear_input,Input* thickness_input,Input* n_input, IssmDouble zeta);
+                void               dViscositydBMOLHO(IssmDouble* pdmudB,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vxbase_input,Input* vybase_input, Input* vxshear_input ,Input* vyshear_input,Input* thickness_input,Input* n_input, IssmDouble zeta);
                 void               dViscositydDSSA(IssmDouble* pdmudB,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input);
                 void               Echo();
                 void               FindParam(bool* pvalue,int paramenum);
@@ -189 +189 @@
                 void               StrainRateFS(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input);
                 void               StrainRateHO(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input);
                 void               StrainRateHO2dvertical(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input);
-                void               StrainRateMLHO(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input,IssmDouble zeta);
+                void               StrainRateMOLHO(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input,IssmDouble zeta);
                 void               StrainRateSSA(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input);
                 void               StrainRateSSA1d(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input);
                 void               StressMaxPrincipalCreateInput(void);
@@ -350 +350 @@
                 virtual int        NumberofNodesVelocity(void)=0;
                 virtual void       PotentialUngrounding(Vector<IssmDouble>* potential_sheet_ungrounding)=0;
                 virtual int        PressureInterpolation()=0;
-      virtual void       Recover3DMLHOInput(int targetVel_enum, int numnodes, IssmDouble* vb,  IssmDouble* vsh, IssmDouble* n, IssmDouble* H, IssmDouble* s){_error_("not implemented yet");};
+      virtual void       Recover3DMOLHOInput(int targetVel_enum, int numnodes, IssmDouble* vb,  IssmDouble* vsh, IssmDouble* n, IssmDouble* H, IssmDouble* s){_error_("not implemented yet");};
                 virtual void       ReduceMatrices(ElementMatrix* Ke,ElementVector* pe)=0;
                 virtual void       ResetFSBasalBoundaryCondition()=0;
                 virtual void       ResetHooks()=0;

../trunk-jpl/src/c/classes/Elements/Penta.cpp

@@ -3281 +3281 @@
         return PentaRef::PressureInterpolation(this->element_type);
 }
 /*}}}*/
-void       Penta::Recover3DMLHOInput(int targetVel_enum, int numnodes, IssmDouble* vb,  IssmDouble* vsh, IssmDouble* n, IssmDouble* H, IssmDouble* s){/*{{{*/
+void       Penta::Recover3DMOLHOInput(int targetVel_enum, int numnodes, IssmDouble* vb,  IssmDouble* vsh, IssmDouble* n, IssmDouble* H, IssmDouble* s){/*{{{*/
    /* Recover the velocity acording to v=vb+(1-\zeta^{n+1})vsh, where \zeta=(s-z)/H
     * The variables vb, vsh, n, H and s are all from the 2D horizontal mesh(Tria), with "numnodes" DOFs
     * To project to penta the DOFs are doubled in size

../trunk-jpl/src/c/classes/Elements/Penta.h

@@ -169 +169 @@
                 int            NumberofNodesVelocity(void);
                 void           PotentialUngrounding(Vector<IssmDouble>* potential_sheet_ungrounding);
                 int            PressureInterpolation();
-      void           Recover3DMLHOInput(int targetVel_enum, int numnodes, IssmDouble* vb,  IssmDouble* vsh, IssmDouble* n, IssmDouble* H, IssmDouble* s);
+      void           Recover3DMOLHOInput(int targetVel_enum, int numnodes, IssmDouble* vb,  IssmDouble* vsh, IssmDouble* n, IssmDouble* H, IssmDouble* s);
                 void           ReduceMatrices(ElementMatrix* Ke,ElementVector* pe);
                 void           ResetFSBasalBoundaryCondition(void);
                 void           ResetHooks();

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

@@ -3218 +3218 @@
         this->parameters->FindParam(&temp,FlowequationIsSIAEnum); iomodel->AddConstant(new IoConstant(temp,"md.flowequation.isSIA"));
         this->parameters->FindParam(&temp,FlowequationIsSSAEnum); iomodel->AddConstant(new IoConstant(temp,"md.flowequation.isSSA"));
         this->parameters->FindParam(&temp,FlowequationIsL1L2Enum); iomodel->AddConstant(new IoConstant(temp,"md.flowequation.isL1L2"));
-        this->parameters->FindParam(&temp,FlowequationIsMLHOEnum); iomodel->AddConstant(new IoConstant(temp,"md.flowequation.isMLHO"));
+        this->parameters->FindParam(&temp,FlowequationIsMOLHOEnum); iomodel->AddConstant(new IoConstant(temp,"md.flowequation.isMOLHO"));
         this->parameters->FindParam(&temp,FlowequationIsHOEnum); iomodel->AddConstant(new IoConstant(temp,"md.flowequation.isHO"));
         this->parameters->FindParam(&temp,FlowequationIsFSEnum); iomodel->AddConstant(new IoConstant(temp,"md.flowequation.isFS"));
         this->parameters->FindParam(&tempint,MasstransportStabilizationEnum); iomodel->AddConstant(new IoConstant(tempint,"md.masstransport.stabilization"));

../trunk-jpl/src/c/classes/Loads/Penpair.cpp

@@ -354 +354 @@
                                 case HOApproximationEnum:   return PenaltyCreateKMatrixStressbalanceSSAHO(kmax); 
                                 default: _error_("Approximation "<<EnumToStringx(approximation1)<<" not supported yet");
                         }
-                case MLHOApproximationEnum:
+                case MOLHOApproximationEnum:
                         switch(approximation1){
-                                case MLHOApproximationEnum:   return PenaltyCreateKMatrixStressbalanceSSAHO(kmax); 
+                                case MOLHOApproximationEnum:   return PenaltyCreateKMatrixStressbalanceSSAHO(kmax); 
                                 default: _error_("Approximation "<<EnumToStringx(approximation1)<<" not supported yet");
                         }
                 case L1L2ApproximationEnum:

../trunk-jpl/src/c/classes/Materials/Material.h

@@ -50 +50 @@
                 virtual void       ViscositySSADerivativeEpsSquare(IssmDouble* pmu_prime,IssmDouble* epsilon,Gauss* gauss)=0;
                 virtual void       ViscosityFS(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input)=0;
                 virtual void       ViscosityHO(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input)=0;
-                virtual void       ViscosityMLHO(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input)=0;
-                virtual void       ViscosityMLHOAdjoint(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input)=0;
+                virtual void       ViscosityMOLHO(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input)=0;
+                virtual void       ViscosityMOLHOAdjoint(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input)=0;
                 virtual void       ViscosityL1L2(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* surf)=0;
                 virtual void       ViscositySSA(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input)=0;
                 virtual void       ViscosityBFS(IssmDouble* pmudB,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input,IssmDouble epseff)=0;

../trunk-jpl/src/c/classes/Materials/Matestar.cpp

@@ -566 +566 @@
         /*Compute viscosity*/
         *pviscosity=GetViscosityGeneral(vx,vy,vz,&dvx[0],&dvy[0],&dvz[0],eps_eff,isdepthaveraged,gauss);
 }/*}}}*/
-void  Matestar::ViscosityMLHO(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input){/*{{{*/
+void  Matestar::ViscosityMOLHO(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input){/*{{{*/
         _error_("not implemented yet");
 }/*}}}*/
-void  Matestar::ViscosityMLHOAdjoint(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input){/*{{{*/
+void  Matestar::ViscosityMOLHOAdjoint(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input){/*{{{*/
         _error_("not implemented yet");
 }/*}}}*/
 void  Matestar::ViscosityL1L2(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* surface_input){/*{{{*/

../trunk-jpl/src/m/classes/inversion.py

@@ -150 +150 @@
 
         # Only SSA, HO and FS are supported right now
         if solution == 'StressbalanceSolution':
-            if not (md.flowequation.isSSA or md.flowequation.isMLHO or md.flowequation.isHO or md.flowequation.isFS or md.flowequation.isL1L2):
-                md.checkmessage("'inversion can only be performed for SSA, MLHO, HO or FS ice flow models")
+            if not (md.flowequation.isSSA or md.flowequation.isMOLHO or md.flowequation.isHO or md.flowequation.isFS or md.flowequation.isL1L2):
+                md.checkmessage("'inversion can only be performed for SSA, MOLHO, HO or FS ice flow models")
         if solution == 'BalancethicknessSolution':
             md = checkfield(md, 'fieldname', 'inversion.thickness_obs', 'size', [md.mesh.numberofvertices], 'NaN', 1, 'Inf', 1)
         elif solution == 'BalancethicknessSoftSolution':

../trunk-jpl/src/m/parameterization/setflowequation.js

@@ -2 +2 @@
 //SETFLOWEQUATION - associate a solution type to each element
 //
 //   This routine works like plotmodel: it works with an even number of inputs
-//   'SIA','SSA','L1L2','MLHO','HO','FS' and 'fill' are the possible options
+//   'SIA','SSA','L1L2','MOLHO','HO','FS' and 'fill' are the possible options
 //   that must be followed by the corresponding exp file or flags list
 //   It can either be a domain file (argus type, .exp extension), or an array of element flags. 
 //   If user wants every element outside the domain to be 
@@ -10 +10 @@
 //   an empty string '' will be considered as an empty domain
 //   a string 'all' will be considered as the entire domain
 //   You can specify the type of coupling, 'penalties' or 'tiling', to use with the input 'coupling'
-//   NB: L1L2 and MLHO cannot currently be coupled to any other ice flow model
+//   NB: L1L2 and MOLHO cannot currently be coupled to any other ice flow model
 //
 //   Usage:
 //      setflowequation(md,varargin)
@@ -37 +37 @@
         SSAflag  = FlagElements(md,options.getfieldvalue('SSA',''));
         HOflag   = FlagElements(md,options.getfieldvalue('HO',''));
         L1L2flag = FlagElements(md,options.getfieldvalue('L1L2',''));
-        MLHOflag = FlagElements(md,options.getfieldvalue('MLHO',''));
+        MOLHOflag = FlagElements(md,options.getfieldvalue('MOLHO',''));
         FSflag   = FlagElements(md,options.getfieldvalue('FS',''));
         filltype = options.getfieldvalue('fill','none');
         options.displayunused();
@@ -54 +54 @@
         }
 
         //check that each element has at least one flag
-        for(var i=0;i<md.mesh.numberofelements;i++)if((SIAflag[i] + SSAflag[i] + HOflag[i] + L1L2flag[i] + MLHOflag[i] + FSflag[i])==0)
+        for(var i=0;i<md.mesh.numberofelements;i++)if((SIAflag[i] + SSAflag[i] + HOflag[i] + L1L2flag[i] + MOLHOflag[i] + FSflag[i])==0)
         throw Error("elements type not assigned, supported models are 'SIA','SSA','HO' and 'FS'");
 
         //check that each element has only one flag
-        if (ArrayAnyAboveStrict(ArrayXPY(SIAflag,SSAflag,HOflag,L1L2flag,MLHOflag),1)){
+        if (ArrayAnyAboveStrict(ArrayXPY(SIAflag,SSAflag,HOflag,L1L2flag,MOLHOflag),1)){
                 console.log('setflowequation warning message: some elements have several types, higher order type is used for them')
 
                 for(var i=0;i<md.mesh.numberofelements;i++){
@@ -70 +70 @@
 
         //check that L1L2 is not coupled to any other model for now
         if (ArrayAnyEqual(L1L2flag,1) & ArrayAnyEqual(ArrayOr(SIAflag,SSAflag,HOflag,FSflag),1)) throw Error('L1L2 cannot be coupled to any other model');
-        if (ArrayAnyEqual(MLHOflag,1) & ArrayAnyEqual(ArrayOr(SIAflag,SSAflag,HOflag,FSflag),1)) throw Error('MLHO cannot be coupled to any other model');
+        if (ArrayAnyEqual(MOLHOflag,1) & ArrayAnyEqual(ArrayOr(SIAflag,SSAflag,HOflag,FSflag),1)) throw Error('MOLHO cannot be coupled to any other model');
 
         //Check that no HO or FS for 2d mesh
         if (md.mesh.domaintype() == '2Dhorizontal'){
@@ -99 +99 @@
         pos=ArrayFind(L1L2flag,1);
         for(var i=0;i<pos.length;i++) for(var j=0;j<md.mesh.elements[0].length;j++) nodeonL1L2[md.mesh.elements[pos[i]][j]-1]=1;
 
-        nodeonMLHO=NewArrayFill(md.mesh.numberofvertices,0);
-        pos=ArrayFind(MLHOflag,1);
-        for(var i=0;i<pos.length;i++) for(var j=0;j<md.mesh.elements[0].length;j++) nodeonMLHO[md.mesh.elements[pos[i]][j]-1]=1;
+        nodeonMOLHO=NewArrayFill(md.mesh.numberofvertices,0);
+        pos=ArrayFind(MOLHOflag,1);
+        for(var i=0;i<pos.length;i++) for(var j=0;j<md.mesh.elements[0].length;j++) nodeonMOLHO[md.mesh.elements[pos[i]][j]-1]=1;
 
         nodeonFS=NewArrayFill(md.mesh.numberofvertices,0);
         noneflag=NewArrayFill(md.mesh.numberofvertices,0);
@@ -261 +261 @@
         pos=ArrayFind(SIAflag,1);for(var i=0;i<pos.length;i++)md.flowequation.element_equation[pos[i]]=1;
         pos=ArrayFind(SSAflag,1);for(var i=0;i<pos.length;i++)md.flowequation.element_equation[pos[i]]=2;
         pos=ArrayFind(L1L2flag,1);for(var i=0;i<pos.length;i++)md.flowequation.element_equation[pos[i]]=3;
-        pos=ArrayFind(MLHOflag,1);for(var i=0;i<pos.length;i++)md.flowequation.element_equation[pos[i]]=4;
+        pos=ArrayFind(MOLHOflag,1);for(var i=0;i<pos.length;i++)md.flowequation.element_equation[pos[i]]=4;
         pos=ArrayFind(HOflag,1);for(var i=0;i<pos.length;i++)md.flowequation.element_equation[pos[i]]=5;
         pos=ArrayFind(FSflag,1);for(var i=0;i<pos.length;i++)md.flowequation.element_equation[pos[i]]=6;
         pos=ArrayFind(SSAHOflag,1);for(var i=0;i<pos.length;i++)md.flowequation.element_equation[pos[i]]=7;
@@ -280 +280 @@
 
         pos=ArrayFind(nodeonSSA,1);for(var i=0;i<pos.length;i++)md.flowequation.vertex_equation[pos[i]]=2;
         pos=ArrayFind(nodeonL1L2,1);for(var i=0;i<pos.length;i++)md.flowequation.vertex_equation[pos[i]]=3;
-        pos=ArrayFind(nodeonMLHO,1);for(var i=0;i<pos.length;i++)md.flowequation.vertex_equation[pos[i]]=4;
+        pos=ArrayFind(nodeonMOLHO,1);for(var i=0;i<pos.length;i++)md.flowequation.vertex_equation[pos[i]]=4;
         pos=ArrayFind(nodeonHO,1);for(var i=0;i<pos.length;i++)md.flowequation.vertex_equation[pos[i]]=5;
         pos=ArrayFind(nodeonFS,1);for(var i=0;i<pos.length;i++)md.flowequation.vertex_equation[pos[i]]=6;
         //DO SIA LAST! Otherwise spcs might not be set up correctly (SIA should have priority)
@@ -300 +300 @@
         md.flowequation.isSIA  = ArrayAnyEqual(md.flowequation.element_equation,1);
         md.flowequation.isSSA  = ArrayAnyEqual(md.flowequation.element_equation,2);
         md.flowequation.isL1L2 = ArrayAnyEqual(md.flowequation.element_equation,3);
-        md.flowequation.isMLHO = ArrayAnyEqual(md.flowequation.element_equation,4);
+        md.flowequation.isMOLHO = ArrayAnyEqual(md.flowequation.element_equation,4);
         md.flowequation.isHO   = ArrayAnyEqual(md.flowequation.element_equation,5);
         md.flowequation.isFS   = ArrayAnyEqual(md.flowequation.element_equation,6);
         return

../trunk-jpl/src/c/shared/Enum/Enum.vim

@@ -112 +112 @@
 syn keyword cConstant CalvingLawEnum
 syn keyword cConstant CalvingMinthicknessEnum
 syn keyword cConstant CalvingTestSpeedfactorEnum
+syn keyword cConstant CalvingTestIndependentRateEnum
 syn keyword cConstant CalvingUseParamEnum
 syn keyword cConstant CalvingScaleThetaEnum
 syn keyword cConstant CalvingAmpAlphaEnum
@@ -174 +175 @@
 syn keyword cConstant FlowequationIsFSEnum
 syn keyword cConstant FlowequationIsHOEnum
 syn keyword cConstant FlowequationIsL1L2Enum
-syn keyword cConstant FlowequationIsMLHOEnum
+syn keyword cConstant FlowequationIsMOLHOEnum
 syn keyword cConstant FlowequationIsSIAEnum
 syn keyword cConstant FlowequationIsSSAEnum
 syn keyword cConstant FlowequationIsNitscheEnum
@@ -657 +658 @@
 syn keyword cConstant BottomPressureOldEnum
 syn keyword cConstant CalvingCalvingrateEnum
 syn keyword cConstant CalvingHabFractionEnum
+syn keyword cConstant CalvingAblationrateEnum
 syn keyword cConstant CalvingMeltingrateEnum
 syn keyword cConstant CalvingStressThresholdFloatingiceEnum
 syn keyword cConstant CalvingStressThresholdGroundediceEnum
@@ -1382 +1384 @@
 syn keyword cConstant InversionVzObsEnum
 syn keyword cConstant JEnum
 syn keyword cConstant L1L2ApproximationEnum
-syn keyword cConstant MLHOApproximationEnum
+syn keyword cConstant MOLHOApproximationEnum
 syn keyword cConstant L2ProjectionBaseAnalysisEnum
 syn keyword cConstant L2ProjectionEPLAnalysisEnum
 syn keyword cConstant LACrouzeixRaviartEnum
@@ -1621 +1623 @@
 syn keyword cType Cfsurfacelogvel
 syn keyword cType Cfsurfacesquare
 syn keyword cType Channel
-syn keyword cType classes
 syn keyword cType Constraint
 syn keyword cType Constraints
 syn keyword cType Contour
@@ -1628 +1629 @@
 syn keyword cType Contours
 syn keyword cType ControlInput
 syn keyword cType Covertree
+syn keyword cType DataSetParam
 syn keyword cType DatasetInput
-syn keyword cType DataSetParam
 syn keyword cType Definition
 syn keyword cType DependentObject
 syn keyword cType DoubleInput
@@ -1642 +1643 @@
 syn keyword cType ElementHook
 syn keyword cType ElementInput
 syn keyword cType ElementMatrix
+syn keyword cType ElementVector
 syn keyword cType Elements
-syn keyword cType ElementVector
 syn keyword cType ExponentialVariogram
 syn keyword cType ExternalResult
 syn keyword cType FemModel
@@ -1650 +1651 @@
 syn keyword cType FileParam
 syn keyword cType Friction
 syn keyword cType Gauss
-syn keyword cType GaussianVariogram
-syn keyword cType gaussobjects
 syn keyword cType GaussPenta
 syn keyword cType GaussSeg
 syn keyword cType GaussTetra
 syn keyword cType GaussTria
+syn keyword cType GaussianVariogram
 syn keyword cType GenericExternalResult
 syn keyword cType GenericOption
 syn keyword cType GenericParam
@@ -1671 +1671 @@
 syn keyword cType IoModel
 syn keyword cType IssmDirectApplicInterface
 syn keyword cType IssmParallelDirectApplicInterface
-syn keyword cType krigingobjects
 syn keyword cType Load
 syn keyword cType Loads
 syn keyword cType Masscon
@@ -1682 +1681 @@
 syn keyword cType Matestar
 syn keyword cType Matice
 syn keyword cType Matlitho
-syn keyword cType matrixobjects
 syn keyword cType MatrixParam
 syn keyword cType Misfit
 syn keyword cType Moulin
@@ -1695 +1693 @@
 syn keyword cType Observation
 syn keyword cType Observations
 syn keyword cType Option
+syn keyword cType OptionUtilities
 syn keyword cType Options
-syn keyword cType OptionUtilities
 syn keyword cType Param
 syn keyword cType Parameters
 syn keyword cType Pengrid
@@ -1710 +1708 @@
 syn keyword cType Radar
 syn keyword cType Regionaloutput
 syn keyword cType Results
+syn keyword cType RiftStruct
 syn keyword cType Riftfront
-syn keyword cType RiftStruct
 syn keyword cType SealevelGeometry
 syn keyword cType Seg
 syn keyword cType SegInput
+syn keyword cType SegRef
 syn keyword cType Segment
-syn keyword cType SegRef
 syn keyword cType SpcDynamic
 syn keyword cType SpcStatic
 syn keyword cType SpcTransient
@@ -1735 +1733 @@
 syn keyword cType VectorParam
 syn keyword cType Vertex
 syn keyword cType Vertices
+syn keyword cType classes
+syn keyword cType gaussobjects
+syn keyword cType krigingobjects
+syn keyword cType matrixobjects
 syn keyword cType AdjointBalancethickness2Analysis
 syn keyword cType AdjointBalancethicknessAnalysis
 syn keyword cType AdjointHorizAnalysis

../trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp

@@ -114 +114 @@
                 case CalvingLawEnum : return "CalvingLaw";
                 case CalvingMinthicknessEnum : return "CalvingMinthickness";
                 case CalvingTestSpeedfactorEnum : return "CalvingTestSpeedfactor";
+                case CalvingTestIndependentRateEnum : return "CalvingTestIndependentRate";
                 case CalvingUseParamEnum : return "CalvingUseParam";
                 case CalvingScaleThetaEnum : return "CalvingScaleTheta";
                 case CalvingAmpAlphaEnum : return "CalvingAmpAlpha";
@@ -176 +177 @@
                 case FlowequationIsFSEnum : return "FlowequationIsFS";
                 case FlowequationIsHOEnum : return "FlowequationIsHO";
                 case FlowequationIsL1L2Enum : return "FlowequationIsL1L2";
-                case FlowequationIsMLHOEnum : return "FlowequationIsMLHO";
+                case FlowequationIsMOLHOEnum : return "FlowequationIsMOLHO";
                 case FlowequationIsSIAEnum : return "FlowequationIsSIA";
                 case FlowequationIsSSAEnum : return "FlowequationIsSSA";
                 case FlowequationIsNitscheEnum : return "FlowequationIsNitsche";
@@ -659 +660 @@
                 case BottomPressureOldEnum : return "BottomPressureOld";
                 case CalvingCalvingrateEnum : return "CalvingCalvingrate";
                 case CalvingHabFractionEnum : return "CalvingHabFraction";
+                case CalvingAblationrateEnum : return "CalvingAblationrate";
                 case CalvingMeltingrateEnum : return "CalvingMeltingrate";
                 case CalvingStressThresholdFloatingiceEnum : return "CalvingStressThresholdFloatingice";
                 case CalvingStressThresholdGroundediceEnum : return "CalvingStressThresholdGroundedice";
@@ -1384 +1386 @@
                 case InversionVzObsEnum : return "InversionVzObs";
                 case JEnum : return "J";
                 case L1L2ApproximationEnum : return "L1L2Approximation";
-                case MLHOApproximationEnum : return "MLHOApproximation";
+                case MOLHOApproximationEnum : return "MOLHOApproximation";
                 case L2ProjectionBaseAnalysisEnum : return "L2ProjectionBaseAnalysis";
                 case L2ProjectionEPLAnalysisEnum : return "L2ProjectionEPLAnalysis";
                 case LACrouzeixRaviartEnum : return "LACrouzeixRaviart";

../trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp

@@ -114 +114 @@
               else if (strcmp(name,"CalvingLaw")==0) return CalvingLawEnum;
               else if (strcmp(name,"CalvingMinthickness")==0) return CalvingMinthicknessEnum;
               else if (strcmp(name,"CalvingTestSpeedfactor")==0) return CalvingTestSpeedfactorEnum;
+              else if (strcmp(name,"CalvingTestIndependentRate")==0) return CalvingTestIndependentRateEnum;
               else if (strcmp(name,"CalvingUseParam")==0) return CalvingUseParamEnum;
               else if (strcmp(name,"CalvingScaleTheta")==0) return CalvingScaleThetaEnum;
               else if (strcmp(name,"CalvingAmpAlpha")==0) return CalvingAmpAlphaEnum;
@@ -135 +136 @@
               else if (strcmp(name,"CumBslcHydroPartition")==0) return CumBslcHydroPartitionEnum;
               else if (strcmp(name,"CumBslcOceanPartition")==0) return CumBslcOceanPartitionEnum;
               else if (strcmp(name,"CumGmtslc")==0) return CumGmtslcEnum;
-              else if (strcmp(name,"CumGmslc")==0) return CumGmslcEnum;
          else stage=2;
    }
    if(stage==2){
-              if (strcmp(name,"DamageC1")==0) return DamageC1Enum;
+              if (strcmp(name,"CumGmslc")==0) return CumGmslcEnum;
+              else if (strcmp(name,"DamageC1")==0) return DamageC1Enum;
               else if (strcmp(name,"DamageC2")==0) return DamageC2Enum;
               else if (strcmp(name,"DamageC3")==0) return DamageC3Enum;
               else if (strcmp(name,"DamageC4")==0) return DamageC4Enum;
@@ -179 +180 @@
               else if (strcmp(name,"FlowequationIsFS")==0) return FlowequationIsFSEnum;
               else if (strcmp(name,"FlowequationIsHO")==0) return FlowequationIsHOEnum;
               else if (strcmp(name,"FlowequationIsL1L2")==0) return FlowequationIsL1L2Enum;
-              else if (strcmp(name,"FlowequationIsMLHO")==0) return FlowequationIsMLHOEnum;
+              else if (strcmp(name,"FlowequationIsMOLHO")==0) return FlowequationIsMOLHOEnum;
               else if (strcmp(name,"FlowequationIsSIA")==0) return FlowequationIsSIAEnum;
               else if (strcmp(name,"FlowequationIsSSA")==0) return FlowequationIsSSAEnum;
               else if (strcmp(name,"FlowequationIsNitsche")==0) return FlowequationIsNitscheEnum;
@@ -258 +259 @@
               else if (strcmp(name,"InputToSmooth")==0) return InputToSmoothEnum;
               else if (strcmp(name,"InversionAlgorithm")==0) return InversionAlgorithmEnum;
               else if (strcmp(name,"InversionControlParameters")==0) return InversionControlParametersEnum;
-              else if (strcmp(name,"InversionControlScalingFactors")==0) return InversionControlScalingFactorsEnum;
          else stage=3;
    }
    if(stage==3){
-              if (strcmp(name,"InversionCostFunctions")==0) return InversionCostFunctionsEnum;
+              if (strcmp(name,"InversionControlScalingFactors")==0) return InversionControlScalingFactorsEnum;
+              else if (strcmp(name,"InversionCostFunctions")==0) return InversionCostFunctionsEnum;
               else if (strcmp(name,"InversionDxmin")==0) return InversionDxminEnum;
               else if (strcmp(name,"InversionGatol")==0) return InversionGatolEnum;
               else if (strcmp(name,"InversionGradientScaling")==0) return InversionGradientScalingEnum;
@@ -381 +382 @@
               else if (strcmp(name,"SamplingRobin")==0) return SamplingRobinEnum;
               else if (strcmp(name,"SamplingSeed")==0) return SamplingSeedEnum;
               else if (strcmp(name,"SaveResults")==0) return SaveResultsEnum;
-              else if (strcmp(name,"SolidearthPartitionIce")==0) return SolidearthPartitionIceEnum;
          else stage=4;
    }
    if(stage==4){
-              if (strcmp(name,"SolidearthPartitionHydro")==0) return SolidearthPartitionHydroEnum;
+              if (strcmp(name,"SolidearthPartitionIce")==0) return SolidearthPartitionIceEnum;
+              else if (strcmp(name,"SolidearthPartitionHydro")==0) return SolidearthPartitionHydroEnum;
               else if (strcmp(name,"SolidearthPartitionOcean")==0) return SolidearthPartitionOceanEnum;
               else if (strcmp(name,"SolidearthNpartIce")==0) return SolidearthNpartIceEnum;
               else if (strcmp(name,"SolidearthNpartOcean")==0) return SolidearthNpartOceanEnum;
@@ -504 +505 @@
               else if (strcmp(name,"SmbIsturbulentflux")==0) return SmbIsturbulentfluxEnum;
               else if (strcmp(name,"SmbK")==0) return SmbKEnum;
               else if (strcmp(name,"SmbLapseRates")==0) return SmbLapseRatesEnum;
-              else if (strcmp(name,"SmbNumBasins")==0) return SmbNumBasinsEnum;
          else stage=5;
    }
    if(stage==5){
-              if (strcmp(name,"SmbNumElevationBins")==0) return SmbNumElevationBinsEnum;
+              if (strcmp(name,"SmbNumBasins")==0) return SmbNumBasinsEnum;
+              else if (strcmp(name,"SmbNumElevationBins")==0) return SmbNumElevationBinsEnum;
               else if (strcmp(name,"SmbNumRequestedOutputs")==0) return SmbNumRequestedOutputsEnum;
               else if (strcmp(name,"SmbPfac")==0) return SmbPfacEnum;
               else if (strcmp(name,"SmbPhi")==0) return SmbPhiEnum;
@@ -627 +628 @@
               else if (strcmp(name,"BasalforcingsDeepwaterMeltingRateValuesAutoregression")==0) return BasalforcingsDeepwaterMeltingRateValuesAutoregressionEnum;
               else if (strcmp(name,"BasalforcingsFloatingiceMeltingRate")==0) return BasalforcingsFloatingiceMeltingRateEnum;
               else if (strcmp(name,"BasalforcingsGeothermalflux")==0) return BasalforcingsGeothermalfluxEnum;
-              else if (strcmp(name,"BasalforcingsGroundediceMeltingRate")==0) return BasalforcingsGroundediceMeltingRateEnum;
          else stage=6;
    }
    if(stage==6){
-              if (strcmp(name,"BasalforcingsLinearBasinId")==0) return BasalforcingsLinearBasinIdEnum;
+              if (strcmp(name,"BasalforcingsGroundediceMeltingRate")==0) return BasalforcingsGroundediceMeltingRateEnum;
+              else if (strcmp(name,"BasalforcingsLinearBasinId")==0) return BasalforcingsLinearBasinIdEnum;
               else if (strcmp(name,"BasalforcingsPerturbationMeltingRate")==0) return BasalforcingsPerturbationMeltingRateEnum;
               else if (strcmp(name,"BasalforcingsSpatialDeepwaterElevation")==0) return BasalforcingsSpatialDeepwaterElevationEnum;
               else if (strcmp(name,"BasalforcingsSpatialDeepwaterMeltingRate")==0) return BasalforcingsSpatialDeepwaterMeltingRateEnum;
@@ -674 +675 @@
               else if (strcmp(name,"BottomPressureOld")==0) return BottomPressureOldEnum;
               else if (strcmp(name,"CalvingCalvingrate")==0) return CalvingCalvingrateEnum;
               else if (strcmp(name,"CalvingHabFraction")==0) return CalvingHabFractionEnum;
+              else if (strcmp(name,"CalvingAblationrate")==0) return CalvingAblationrateEnum;
               else if (strcmp(name,"CalvingMeltingrate")==0) return CalvingMeltingrateEnum;
               else if (strcmp(name,"CalvingStressThresholdFloatingice")==0) return CalvingStressThresholdFloatingiceEnum;
               else if (strcmp(name,"CalvingStressThresholdGroundedice")==0) return CalvingStressThresholdGroundediceEnum;
@@ -749 +751 @@
               else if (strcmp(name,"FrictionM")==0) return FrictionMEnum;
               else if (strcmp(name,"FrictionP")==0) return FrictionPEnum;
               else if (strcmp(name,"FrictionPressureAdjustedTemperature")==0) return FrictionPressureAdjustedTemperatureEnum;
-              else if (strcmp(name,"FrictionQ")==0) return FrictionQEnum;
-              else if (strcmp(name,"FrictionSedimentCompressibilityCoefficient")==0) return FrictionSedimentCompressibilityCoefficientEnum;
          else stage=7;
    }
    if(stage==7){
-              if (strcmp(name,"FrictionTillFrictionAngle")==0) return FrictionTillFrictionAngleEnum;
+              if (strcmp(name,"FrictionQ")==0) return FrictionQEnum;
+              else if (strcmp(name,"FrictionSedimentCompressibilityCoefficient")==0) return FrictionSedimentCompressibilityCoefficientEnum;
+              else if (strcmp(name,"FrictionTillFrictionAngle")==0) return FrictionTillFrictionAngleEnum;
               else if (strcmp(name,"FrictionWaterLayer")==0) return FrictionWaterLayerEnum;
               else if (strcmp(name,"Frictionf")==0) return FrictionfEnum;
               else if (strcmp(name,"FrontalForcingsBasinId")==0) return FrontalForcingsBasinIdEnum;
@@ -872 +874 @@
               else if (strcmp(name,"SealevelBarystaticIceMask")==0) return SealevelBarystaticIceMaskEnum;
               else if (strcmp(name,"SealevelBarystaticIceWeights")==0) return SealevelBarystaticIceWeightsEnum;
               else if (strcmp(name,"SealevelBarystaticIceArea")==0) return SealevelBarystaticIceAreaEnum;
-              else if (strcmp(name,"SealevelBarystaticIceLatbar")==0) return SealevelBarystaticIceLatbarEnum;
-              else if (strcmp(name,"SealevelBarystaticIceLongbar")==0) return SealevelBarystaticIceLongbarEnum;
          else stage=8;
    }
    if(stage==8){
-              if (strcmp(name,"SealevelBarystaticIceLoad")==0) return SealevelBarystaticIceLoadEnum;
+              if (strcmp(name,"SealevelBarystaticIceLatbar")==0) return SealevelBarystaticIceLatbarEnum;
+              else if (strcmp(name,"SealevelBarystaticIceLongbar")==0) return SealevelBarystaticIceLongbarEnum;
+              else if (strcmp(name,"SealevelBarystaticIceLoad")==0) return SealevelBarystaticIceLoadEnum;
               else if (strcmp(name,"SealevelBarystaticHydroMask")==0) return SealevelBarystaticHydroMaskEnum;
               else if (strcmp(name,"SealevelBarystaticHydroWeights")==0) return SealevelBarystaticHydroWeightsEnum;
               else if (strcmp(name,"SealevelBarystaticHydroArea")==0) return SealevelBarystaticHydroAreaEnum;
@@ -995 +997 @@
               else if (strcmp(name,"SmbGdn")==0) return SmbGdnEnum;
               else if (strcmp(name,"SmbGdnini")==0) return SmbGdniniEnum;
               else if (strcmp(name,"SmbGsp")==0) return SmbGspEnum;
-              else if (strcmp(name,"SmbGspini")==0) return SmbGspiniEnum;
-              else if (strcmp(name,"SmbHref")==0) return SmbHrefEnum;
          else stage=9;
    }
    if(stage==9){
-              if (strcmp(name,"SmbIsInitialized")==0) return SmbIsInitializedEnum;
+              if (strcmp(name,"SmbGspini")==0) return SmbGspiniEnum;
+              else if (strcmp(name,"SmbHref")==0) return SmbHrefEnum;
+              else if (strcmp(name,"SmbIsInitialized")==0) return SmbIsInitializedEnum;
               else if (strcmp(name,"SmbMAdd")==0) return SmbMAddEnum;
               else if (strcmp(name,"SmbMassBalance")==0) return SmbMassBalanceEnum;
               else if (strcmp(name,"SmbMassBalanceSubstep")==0) return SmbMassBalanceSubstepEnum;
@@ -1118 +1120 @@
               else if (strcmp(name,"VxShear")==0) return VxShearEnum;
               else if (strcmp(name,"VxSurface")==0) return VxSurfaceEnum;
               else if (strcmp(name,"VyAverage")==0) return VyAverageEnum;
-              else if (strcmp(name,"VyBase")==0) return VyBaseEnum;
-              else if (strcmp(name,"Vy")==0) return VyEnum;
          else stage=10;
    }
    if(stage==10){
-              if (strcmp(name,"VyMesh")==0) return VyMeshEnum;
+              if (strcmp(name,"VyBase")==0) return VyBaseEnum;
+              else if (strcmp(name,"Vy")==0) return VyEnum;
+              else if (strcmp(name,"VyMesh")==0) return VyMeshEnum;
               else if (strcmp(name,"VyObs")==0) return VyObsEnum;
               else if (strcmp(name,"VyShear")==0) return VyShearEnum;
               else if (strcmp(name,"VySurface")==0) return VySurfaceEnum;
@@ -1241 +1243 @@
               else if (strcmp(name,"Outputdefinition95")==0) return Outputdefinition95Enum;
               else if (strcmp(name,"Outputdefinition96")==0) return Outputdefinition96Enum;
               else if (strcmp(name,"Outputdefinition97")==0) return Outputdefinition97Enum;
-              else if (strcmp(name,"Outputdefinition98")==0) return Outputdefinition98Enum;
-              else if (strcmp(name,"Outputdefinition99")==0) return Outputdefinition99Enum;
          else stage=11;
    }
    if(stage==11){
-              if (strcmp(name,"Outputdefinition9")==0) return Outputdefinition9Enum;
+              if (strcmp(name,"Outputdefinition98")==0) return Outputdefinition98Enum;
+              else if (strcmp(name,"Outputdefinition99")==0) return Outputdefinition99Enum;
+              else if (strcmp(name,"Outputdefinition9")==0) return Outputdefinition9Enum;
               else if (strcmp(name,"Outputdefinition100")==0) return Outputdefinition100Enum;
               else if (strcmp(name,"InputsEND")==0) return InputsENDEnum;
               else if (strcmp(name,"Absolute")==0) return AbsoluteEnum;
@@ -1364 +1366 @@
               else if (strcmp(name,"GaussSeg")==0) return GaussSegEnum;
               else if (strcmp(name,"GaussTetra")==0) return GaussTetraEnum;
               else if (strcmp(name,"GaussTria")==0) return GaussTriaEnum;
-              else if (strcmp(name,"GenericOption")==0) return GenericOptionEnum;
-              else if (strcmp(name,"GenericParam")==0) return GenericParamEnum;
          else stage=12;
    }
    if(stage==12){
-              if (strcmp(name,"GenericExternalResult")==0) return GenericExternalResultEnum;
+              if (strcmp(name,"GenericOption")==0) return GenericOptionEnum;
+              else if (strcmp(name,"GenericParam")==0) return GenericParamEnum;
+              else if (strcmp(name,"GenericExternalResult")==0) return GenericExternalResultEnum;
               else if (strcmp(name,"Gradient1")==0) return Gradient1Enum;
               else if (strcmp(name,"Gradient2")==0) return Gradient2Enum;
               else if (strcmp(name,"Gradient3")==0) return Gradient3Enum;
@@ -1417 +1419 @@
               else if (strcmp(name,"InversionVzObs")==0) return InversionVzObsEnum;
               else if (strcmp(name,"J")==0) return JEnum;
               else if (strcmp(name,"L1L2Approximation")==0) return L1L2ApproximationEnum;
-              else if (strcmp(name,"MLHOApproximation")==0) return MLHOApproximationEnum;
+              else if (strcmp(name,"MOLHOApproximation")==0) return MOLHOApproximationEnum;
               else if (strcmp(name,"L2ProjectionBaseAnalysis")==0) return L2ProjectionBaseAnalysisEnum;
               else if (strcmp(name,"L2ProjectionEPLAnalysis")==0) return L2ProjectionEPLAnalysisEnum;
               else if (strcmp(name,"LACrouzeixRaviart")==0) return LACrouzeixRaviartEnum;
@@ -1487 +1489 @@
               else if (strcmp(name,"NoMeltOnPartiallyFloating")==0) return NoMeltOnPartiallyFloatingEnum;
               else if (strcmp(name,"Nodal")==0) return NodalEnum;
               else if (strcmp(name,"Nodalvalue")==0) return NodalvalueEnum;
-              else if (strcmp(name,"NodeSId")==0) return NodeSIdEnum;
-              else if (strcmp(name,"NoneApproximation")==0) return NoneApproximationEnum;
          else stage=13;
    }
    if(stage==13){
-              if (strcmp(name,"None")==0) return NoneEnum;
+              if (strcmp(name,"NodeSId")==0) return NodeSIdEnum;
+              else if (strcmp(name,"NoneApproximation")==0) return NoneApproximationEnum;
+              else if (strcmp(name,"None")==0) return NoneEnum;
               else if (strcmp(name,"Numberedcostfunction")==0) return NumberedcostfunctionEnum;
               else if (strcmp(name,"NyeCO2")==0) return NyeCO2Enum;
               else if (strcmp(name,"NyeH2O")==0) return NyeH2OEnum;
@@ -1610 +1612 @@
               else if (strcmp(name,"TotalSmb")==0) return TotalSmbEnum;
               else if (strcmp(name,"TotalSmbScaled")==0) return TotalSmbScaledEnum;
               else if (strcmp(name,"TransientArrayParam")==0) return TransientArrayParamEnum;
-              else if (strcmp(name,"TransientInput")==0) return TransientInputEnum;
-              else if (strcmp(name,"TransientParam")==0) return TransientParamEnum;
          else stage=14;
    }
    if(stage==14){
-              if (strcmp(name,"TransientSolution")==0) return TransientSolutionEnum;
+              if (strcmp(name,"TransientInput")==0) return TransientInputEnum;
+              else if (strcmp(name,"TransientParam")==0) return TransientParamEnum;
+              else if (strcmp(name,"TransientSolution")==0) return TransientSolutionEnum;
               else if (strcmp(name,"Tria")==0) return TriaEnum;
               else if (strcmp(name,"TriaInput")==0) return TriaInputEnum;
               else if (strcmp(name,"UzawaPressureAnalysis")==0) return UzawaPressureAnalysisEnum;

../trunk-jpl/src/m/boundaryconditions/SetMOLHOBC.m

@@ +1 @@
+function md=SetMOLHOBC(md)
+%SETMOLHOBC - Create the boundary conditions for stressbalance for MOLHO: VxBase, VyBase, VxShear, VyShear
+%
+%   Usage:
+%      md=SetMOLHOBC(md)
+%
+
+
+%node on Dirichlet
+if md.flowequation.isMOLHO
+        md.stressbalance.spcvx_base=md.stressbalance.spcvx;
+        md.stressbalance.spcvy_base=md.stressbalance.spcvy;
+
+        md.stressbalance.spcvx_shear=NaN*ones(size(md.stressbalance.spcvx_base));
+        md.stressbalance.spcvy_shear=NaN*ones(size(md.stressbalance.spcvy_base));
+end

../trunk-jpl/src/m/classes/flowequation.js

@@ -21 +21 @@
                 fielddisplay(this,'isSIA','is the Shallow Ice Approximation (SIA) used ?');
                 fielddisplay(this,'isSSA','is the Shelfy-Stream Approximation (SSA) used ?');
                 fielddisplay(this,'isL1L2','is the L1L2 approximation used ?');
-                fielddisplay(this,'isMLHO','is the Mono-Layer Higher-Order approximation used?');
+                fielddisplay(this,'isMOLHO','is the MOno-Layer Higher-Order (MOLHO) approximation used?');
                 fielddisplay(this,'isHO','is the Higher-Order (HO) approximation used ?');
                 fielddisplay(this,'isFS','are the Full-FS (FS) equations used ?');
                 fielddisplay(this,'isNitscheBC','is weakly imposed condition used?');
@@ -40 +40 @@
                 return "flowequation";
 
         }// }}}
-    this.extrude = function(md) {//{{{
-        this.element_equation=project3d(md,'vector',this.element_equation,'type','element');
-        this.vertex_equation=project3d(md,'vector',this.vertex_equation,'type','node');
-        this.borderSSA=project3d(md,'vector',this.borderSSA,'type','node');
-        this.borderHO=project3d(md,'vector',this.borderHO,'type','node');
-        this.borderFS=project3d(md,'vector',this.borderFS,'type','node');
-        return this;
+        this.extrude = function(md) {//{{{
+                this.element_equation=project3d(md,'vector',this.element_equation,'type','element');
+                this.vertex_equation=project3d(md,'vector',this.vertex_equation,'type','node');
+                this.borderSSA=project3d(md,'vector',this.borderSSA,'type','node');
+                this.borderHO=project3d(md,'vector',this.borderHO,'type','node');
+                this.borderFS=project3d(md,'vector',this.borderFS,'type','node');
+                return this;
     }//}}}
                 this.checkconsistency = function(md,solution,analyses) {//{{{
 
@@ -58 +58 @@
                         checkfield(md,'fieldname','flowequation.isSIA','numel',[1],'values',[0, 1]);
                         checkfield(md,'fieldname','flowequation.isSSA','numel',[1],'values',[0, 1]);
                         checkfield(md,'fieldname','flowequation.isL1L2','numel',[1],'values',[0, 1]);
-                        checkfield(md,'fieldname','flowequation.isMLHO','numel',[1],'values',[0, 1]);
+                        checkfield(md,'fieldname','flowequation.isMOLHO','numel',[1],'values',[0, 1]);
                         checkfield(md,'fieldname','flowequation.isHO','numel',[1],'values',[0, 1]);
                         checkfield(md,'fieldname','flowequation.isFS','numel',[1],'values',[0, 1]);
                         checkfield(md,'fieldname','flowequation.isNitscheBC','numel',[1],'values',[0, 1]);
@@ -92 +92 @@
                         }
                         else throw Error('Case not supported yet');
                         
-                        if (!(this.isSIA | this.isSSA | this.isL1L2 | this.isMLHO | this.isHO | this.isFS)){
+                        if (!(this.isSIA | this.isSSA | this.isL1L2 | this.isMOLHO | this.isHO | this.isFS)){
                                 checkmessage(md,['no element types set for this model']);
                         }
                         if(ArrayAnyEqual(ArrayIsMember('StressbalanceSIAAnalysis', analyses),1)){
@@ -107 +107 @@
                         WriteData(fid,prefix,'object',this,'fieldname','isSIA','format','Boolean');
                         WriteData(fid,prefix,'object',this,'fieldname','isSSA','format','Boolean');
                         WriteData(fid,prefix,'object',this,'fieldname','isL1L2','format','Boolean');
-                        WriteData(fid,prefix,'object',this,'fieldname','isMLHO','format','Boolean');
+                        WriteData(fid,prefix,'object',this,'fieldname','isMOLHO','format','Boolean');
                         WriteData(fid,prefix,'object',this,'fieldname','isHO','format','Boolean');
                         WriteData(fid,prefix,'object',this,'fieldname','isFS','format','Boolean');
                         WriteData(fid,prefix,'object',this,'fieldname','isNitscheBC','format','Boolean');
@@ -137 +137 @@
         this.isSIA                          = 0;
         this.isSSA                          = 0;
         this.isL1L2                         = 0;
-        this.isMLHO                         = 0;
+        this.isMOLHO                         = 0;
         this.isHO                           = 0;
         this.isFS                           = 0;
         this.isNitscheBC                    = 0;

../trunk-jpl/src/m/classes/inversion.m

@@ -103 +103 @@
 
                         %Only SSA, HO and FS are supported right now
                         if strcmp(solution,'StressbalanceSolution')
-                                if ~(md.flowequation.isSSA || md.flowequation.isMLHO || md.flowequation.isHO || md.flowequation.isFS || md.flowequation.isL1L2),
-                                        md = checkmessage(md,['inversion can only be performed for SSA, MLHO, HO or FS ice flow models']);
+                                if ~(md.flowequation.isSSA || md.flowequation.isMOLHO || md.flowequation.isHO || md.flowequation.isFS || md.flowequation.isL1L2),
+                                        md = checkmessage(md,['inversion can only be performed for SSA, MOLHO, HO or FS ice flow models']);
                                 end
                         end
                         if strcmp(solution,'BalancethicknessSolution')

../trunk-jpl/src/m/classes/stressbalance.py

@@ -58 +58 @@
         s += '{}\n'.format(fielddisplay(self, 'spcvy', 'y-axis velocity constraint (NaN means no constraint) [m / yr]'))
         s += '{}\n'.format(fielddisplay(self, 'spcvz', 'z-axis velocity constraint (NaN means no constraint) [m / yr]'))
         s += '{}\n'.format(fielddisplay(self, 'icefront', 'segments on ice front list (last column 0: Air, 1: Water, 2: Ice'))
-        s += '      MLHO boundary conditions:\n'
+        s += '      MOLHO boundary conditions:\n'
         s += '{}\n'.format(fielddisplay(self, 'spcvx_base', 'x-axis basal velocity constraint (NaN means no constraint) [m / yr]'))
         s += '{}\n'.format(fielddisplay(self, 'spcvy_base', 'y-axis basal velocity constraint (NaN means no constraint) [m / yr]'))
         s += '{}\n'.format(fielddisplay(self, 'spcvx_shear', 'x-axis shear velocity constraint (NaN means no constraint) [m / yr]'))
@@ -85 +85 @@
         self.referential = project3d(md, 'vector', self.referential, 'type', 'node')
         self.loadingforce = project3d(md, 'vector', self.loadingforce, 'type', 'node')
 
-        if md.flowequation.isMLHO:
+        if md.flowequation.isMOLHO:
             self.spcvx_base = project3d(md, 'vector', self.spcvx_base, 'type', 'node')
             self.spcvy_base = project3d(md, 'vector', self.spcvy_base, 'type', 'node')
             self.spcvx_shear = project3d(md, 'vector', self.spcvx_shear, 'type', 'poly', 'degree', 4)
@@ -175 +175 @@
             pos = np.nonzero(np.logical_and(md.mask.ocean_levelset, md.mesh.vertexonbase))
             if np.any(np.logical_not(np.isnan(md.stressbalance.referential[pos, :]))):
                 md.checkmessage("no referential should be specified for basal vertices of grounded ice")
-        if md.flowequation.isMLHO:
+        if md.flowequation.isMOLHO:
             md = checkfield(md, 'fieldname', 'stressbalance.spcvx_base', 'Inf', 1, 'timeseries', 1)
             md = checkfield(md, 'fieldname', 'stressbalance.spcvy_base', 'Inf', 1, 'timeseries', 1)
             md = checkfield(md, 'fieldname', 'stressbalance.spcvx_shear', 'Inf', 1, 'timeseries', 1)
@@ -214 +214 @@
             outputscopy = outputs[0:max(0, indices[0] - 1)] + self.defaultoutputs(md) + outputs[indices[0] + 1:]
             outputs = outputscopy
         WriteData(fid, prefix, 'data', outputs, 'name', 'md.stressbalance.requested_outputs', 'format', 'StringArray')
-        # MLHO
-        if md.flowequation.isMLHO:
+        # MOLHO
+        if md.flowequation.isMOLHO:
             WriteData(fid, prefix, 'object', self, 'class', 'stressbalance', 'fieldname', 'spcvx_base', 'format', 'DoubleMat', 'mattype', 1, 'scale', 1. / yts, 'timeserieslength', md.mesh.numberofvertices + 1, 'yts', yts)
             WriteData(fid, prefix, 'object', self, 'class', 'stressbalance', 'fieldname', 'spcvy_base', 'format', 'DoubleMat', 'mattype', 1, 'scale', 1. / yts, 'timeserieslength', md.mesh.numberofvertices + 1, 'yts', yts)
             WriteData(fid, prefix, 'object', self, 'class', 'stressbalance', 'fieldname', 'spcvx_shear', 'format', 'DoubleMat', 'mattype', 1, 'scale', 1. / yts, 'timeserieslength', md.mesh.numberofvertices + 1, 'yts', yts)

../trunk-jpl/src/c/classes/Materials/Matestar.h

@@ -75 +75 @@
 
                 void       ViscosityFS(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input);
                 void       ViscosityHO(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input);
-                void       ViscosityMLHO(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input);
-                void       ViscosityMLHOAdjoint(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input);
+                void       ViscosityMOLHO(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input);
+                void       ViscosityMOLHOAdjoint(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input);
                 void       ViscosityL1L2(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* surf);
                 void       ViscositySSA(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input);
                 void       ViscosityBFS(IssmDouble* pmudB,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input,IssmDouble eps_eff);

../trunk-jpl/src/c/classes/Materials/Matice.cpp

@@ -733 +733 @@
         /*Assign output pointer*/
         *pviscosity=viscosity;
 }/*}}}*/
-void  Matice::ViscosityMLHO(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vxbase_input,Input* vybase_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input){/*{{{*/
+void  Matice::ViscosityMOLHO(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vxbase_input,Input* vybase_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input){/*{{{*/
 
         /*Intermediaries*/
         IssmDouble epsilon[5];  /* epsilon=[exx,eyy,exy,exz,eyz]; */
@@ -754 +754 @@
                 zeta=0.5*(gauss_seg->coord1+1); 
 
                 /* eps_eff^2 = exx^2 + eyy^2 + exy^2 + exz^2 + eyz^2 + exx*eyy (for a given zeta)*/
-                element->StrainRateMLHO(&epsilon[0],xyz_list,gauss,
+                element->StrainRateMOLHO(&epsilon[0],xyz_list,gauss,
                                                 vxbase_input,vybase_input,vxshear_input,vyshear_input,thickness_input,n_input,zeta);
                 epsilon_eff=sqrt(epsilon[0]*epsilon[0] + epsilon[1]*epsilon[1] + epsilon[2]*epsilon[2] 
                                                   +  epsilon[3]*epsilon[3] + epsilon[4]*epsilon[4] + epsilon[0]*epsilon[1]);
@@ -794 +794 @@
         /*Clean up*/
         delete gauss_seg;
 }/*}}}*/
-void  Matice::ViscosityMLHOAdjoint(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vxbase_input,Input* vybase_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input){/*{{{*/
+void  Matice::ViscosityMOLHOAdjoint(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vxbase_input,Input* vybase_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input){/*{{{*/
 
         /* To compute the additional 5 terms in the viscosity appear in the adjoint equation*/
         /*Intermediaries*/
@@ -816 +816 @@
                 zeta=0.5*(gauss_seg->coord1+1); 
 
                 /* eps_eff^2 = exx^2 + eyy^2 + exy^2 + exz^2 + eyz^2 + exx*eyy (for a given zeta)*/
-                element->StrainRateMLHO(&epsilon[0],xyz_list,gauss,
+                element->StrainRateMOLHO(&epsilon[0],xyz_list,gauss,
                                                 vxbase_input,vybase_input,vxshear_input,vyshear_input,thickness_input,n_input,zeta);
                 epsilon_eff=sqrt(epsilon[0]*epsilon[0] + epsilon[1]*epsilon[1] + epsilon[2]*epsilon[2] 
                                                   +  epsilon[3]*epsilon[3] + epsilon[4]*epsilon[4] + epsilon[0]*epsilon[1]);

../trunk-jpl/src/c/classes/Materials/Matice.h

@@ -77 +77 @@
 
                 void       ViscosityFS(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input);
                 void       ViscosityHO(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input);
-                void       ViscosityMLHO(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input);
-                void       ViscosityMLHOAdjoint(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input);
+                void       ViscosityMOLHO(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input);
+                void       ViscosityMOLHOAdjoint(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input);
                 void       ViscosityL1L2(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* surf);
                 void       ViscositySSA(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input);
                 void       ViscosityBFS(IssmDouble* pmudB,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input,IssmDouble eps_eff){_error_("not supported");};

../trunk-jpl/src/c/classes/Materials/Matlitho.h

@@ -71 +71 @@
 
                 void       ViscosityFS(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input){_error_("not supported");};
                 void       ViscosityHO(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input){_error_("not supported");};
-                void       ViscosityMLHO(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input){_error_("not supported");};
-                void       ViscosityMLHOAdjoint(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input){_error_("not supported");};
+                void       ViscosityMOLHO(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input){_error_("not supported");};
+                void       ViscosityMOLHOAdjoint(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vxshear_input,Input* vyshear_input,Input* thickness_input,Input* n_input){_error_("not supported");};
                 void       ViscosityL1L2(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* surf){_error_("not supported");};
                 void       ViscositySSA(IssmDouble* pviscosity,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input){_error_("not supported");};
                 void       ViscosityBFS(IssmDouble* pmudB,int dim,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input,IssmDouble epseff){_error_("not supported");};

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

@@ -117 +117 @@
                         if(in_approximation==L1L2ApproximationEnum && !reCast<int>(iomodel->Data("md.mesh.vertexonbase")[io_index])){
                                 this->HardDeactivate();
                         }
-                        if(in_approximation==MLHOApproximationEnum && !reCast<int>(iomodel->Data("md.mesh.vertexonbase")[io_index])){
+                        if(in_approximation==MOLHOApproximationEnum && !reCast<int>(iomodel->Data("md.mesh.vertexonbase")[io_index])){
                                 this->HardDeactivate();
                         }
                         if(in_approximation==SSAHOApproximationEnum && reCast<int>(iomodel->Data("md.flowequation.borderSSA")[io_index])){

../trunk-jpl/src/c/cores/stressbalance_core.cpp

@@ -18 +18 @@
         /*parameters: */
         bool       dakota_analysis,control_analysis;
         int        domaintype;
-        bool       isSIA,isSSA,isL1L2,isMLHO,isHO,isFS,isNitsche;
+        bool       isSIA,isSSA,isL1L2,isMOLHO,isHO,isFS,isNitsche;
         bool       save_results;
         int        solution_type;
         int        numoutputs        = 0;
@@ -30 +30 @@
         femmodel->parameters->FindParam(&isSIA,FlowequationIsSIAEnum);
         femmodel->parameters->FindParam(&isSSA,FlowequationIsSSAEnum);
         femmodel->parameters->FindParam(&isL1L2,FlowequationIsL1L2Enum);
-        femmodel->parameters->FindParam(&isMLHO,FlowequationIsMLHOEnum);
+        femmodel->parameters->FindParam(&isMOLHO,FlowequationIsMOLHOEnum);
         femmodel->parameters->FindParam(&isHO,FlowequationIsHOEnum);
         femmodel->parameters->FindParam(&isFS,FlowequationIsFSEnum);
         femmodel->parameters->FindParam(&isNitsche,FlowequationIsNitscheEnum);
@@ -74 +74 @@
         }
 
         /*Compute stressbalance for SSA L1L2 HO and FS*/
-        if(isSSA || isL1L2 || isMLHO || isHO || isFS){
+        if(isSSA || isL1L2 || isMOLHO || isHO || isFS){
                 analysis = new StressbalanceAnalysis();
                 analysis->Core(femmodel);
                 delete analysis;
@@ -81 +81 @@
         }
 
         /*Compute vertical velocities*/
-        if (domaintype==Domain3DEnum && (isSIA || isSSA || isL1L2 || isMLHO || isHO)){
+        if (domaintype==Domain3DEnum && (isSIA || isSSA || isL1L2 || isMOLHO || isHO)){
 
                 /*We need basal melt rates for vertical velocity*/
                 bmb_core(femmodel);

../trunk-jpl/src/c/modules/SurfaceAbsVelMisfitx/SurfaceAbsVelMisfitx.cpp

@@ -65 +65 @@
         Input* vy_input     = NULL;
         Input* vyobs_input  = NULL;
 
-        /*Read SurfaceEnum from 2D models:SSA, L1L2, MLHO*/
+        /*Read SurfaceEnum from 2D models:SSA, L1L2, MOLHO*/
         if (domaintype == Domain2DhorizontalEnum) {
                 vx_input = topelement->GetInput(VxSurfaceEnum);                                 _assert_(vx_input);
                 if(numcomponents==2){

../trunk-jpl/src/c/modules/SurfaceAverageVelMisfitx/SurfaceAverageVelMisfitx.cpp

@@ -73 +73 @@
    Input* vy_input     = NULL;
    Input* vyobs_input  = NULL;
 
-   /*Read SurfaceEnum from 2D models:SSA, L1L2, MLHO*/
+   /*Read SurfaceEnum from 2D models:SSA, L1L2, MOLHO*/
    if (domaintype == Domain2DhorizontalEnum) {
       vx_input = topelement->GetInput(VxSurfaceEnum);             _assert_(vx_input);
       if(numcomponents==2){

../trunk-jpl/src/c/modules/SurfaceLogVelMisfitx/SurfaceLogVelMisfitx.cpp

@@ -67 +67 @@
    Input* vy_input     = NULL;
    Input* vyobs_input  = NULL;
 
-   /*Read SurfaceEnum from 2D models:SSA, L1L2, MLHO*/
+   /*Read SurfaceEnum from 2D models:SSA, L1L2, MOLHO*/
    if (domaintype == Domain2DhorizontalEnum) {
       vx_input = topelement->GetInput(VxSurfaceEnum);             _assert_(vx_input);
       if(numcomponents==2){

../trunk-jpl/src/c/modules/SurfaceLogVxVyMisfitx/SurfaceLogVxVyMisfitx.cpp

@@ -67 +67 @@
    Input* vy_input     = NULL;
    Input* vyobs_input  = NULL;
 
-   /*Read SurfaceEnum from 2D models:SSA, L1L2, MLHO*/
+   /*Read SurfaceEnum from 2D models:SSA, L1L2, MOLHO*/
    if (domaintype == Domain2DhorizontalEnum) {
       vx_input = topelement->GetInput(VxSurfaceEnum);             _assert_(vx_input);
       if(numcomponents==2){

../trunk-jpl/src/c/modules/SurfaceRelVelMisfitx/SurfaceRelVelMisfitx.cpp

@@ -67 +67 @@
    Input* vy_input     = NULL;
    Input* vyobs_input  = NULL;
 
-   /*Read SurfaceEnum from 2D models:SSA, L1L2, MLHO*/
+   /*Read SurfaceEnum from 2D models:SSA, L1L2, MOLHO*/
    if (domaintype == Domain2DhorizontalEnum) {
       vx_input = topelement->GetInput(VxSurfaceEnum);             _assert_(vx_input);
       if(numcomponents==2){

../trunk-jpl/src/c/shared/Enum/EnumDefinitions.h

@@ -169 +169 @@
         FlowequationIsFSEnum,
         FlowequationIsHOEnum,
         FlowequationIsL1L2Enum,
-        FlowequationIsMLHOEnum,
+        FlowequationIsMOLHOEnum,
         FlowequationIsSIAEnum,
         FlowequationIsSSAEnum,
         FlowequationIsNitscheEnum,
@@ -1383 +1383 @@
         InversionVzObsEnum,
         JEnum,
         L1L2ApproximationEnum,
-        MLHOApproximationEnum,
+        MOLHOApproximationEnum,
         L2ProjectionBaseAnalysisEnum,
         L2ProjectionEPLAnalysisEnum,
         LACrouzeixRaviartEnum,

../trunk-jpl/src/c/shared/Enum/Enumjl.vim

@@ -105 +105 @@
 syn keyword juliaConstC CalvingLawEnum
 syn keyword juliaConstC CalvingMinthicknessEnum
 syn keyword juliaConstC CalvingTestSpeedfactorEnum
+syn keyword juliaConstC CalvingTestIndependentRateEnum
 syn keyword juliaConstC CalvingUseParamEnum
 syn keyword juliaConstC CalvingScaleThetaEnum
 syn keyword juliaConstC CalvingAmpAlphaEnum
@@ -167 +168 @@
 syn keyword juliaConstC FlowequationIsFSEnum
 syn keyword juliaConstC FlowequationIsHOEnum
 syn keyword juliaConstC FlowequationIsL1L2Enum
-syn keyword juliaConstC FlowequationIsMLHOEnum
+syn keyword juliaConstC FlowequationIsMOLHOEnum
 syn keyword juliaConstC FlowequationIsSIAEnum
 syn keyword juliaConstC FlowequationIsSSAEnum
 syn keyword juliaConstC FlowequationIsNitscheEnum
@@ -650 +651 @@
 syn keyword juliaConstC BottomPressureOldEnum
 syn keyword juliaConstC CalvingCalvingrateEnum
 syn keyword juliaConstC CalvingHabFractionEnum
+syn keyword juliaConstC CalvingAblationrateEnum
 syn keyword juliaConstC CalvingMeltingrateEnum
 syn keyword juliaConstC CalvingStressThresholdFloatingiceEnum
 syn keyword juliaConstC CalvingStressThresholdGroundediceEnum
@@ -1375 +1377 @@
 syn keyword juliaConstC InversionVzObsEnum
 syn keyword juliaConstC JEnum
 syn keyword juliaConstC L1L2ApproximationEnum
-syn keyword juliaConstC MLHOApproximationEnum
+syn keyword juliaConstC MOLHOApproximationEnum
 syn keyword juliaConstC L2ProjectionBaseAnalysisEnum
 syn keyword juliaConstC L2ProjectionEPLAnalysisEnum
 syn keyword juliaConstC LACrouzeixRaviartEnum

../trunk-jpl/src/c/shared/io/Marshalling/IoCodeConversions.cpp

@@ -344 +344 @@
                 case 1: return SIAApproximationEnum;
                 case 2: return SSAApproximationEnum;
                 case 3: return L1L2ApproximationEnum;
-                case 4: return MLHOApproximationEnum;
+                case 4: return MOLHOApproximationEnum;
                 case 5: return HOApproximationEnum;
                 case 6: return FSApproximationEnum;
                 case 7: return SSAHOApproximationEnum;
@@ -359 +359 @@
                 case 1: return SIAApproximationEnum;
                 case 2: return SSAApproximationEnum;
                 case 3: return L1L2ApproximationEnum;
-                case 4: return MLHOApproximationEnum;
+                case 4: return MOLHOApproximationEnum;
                 case 5: return HOApproximationEnum;
                 case 6: return FSApproximationEnum;
                 case 7: return SSAHOApproximationEnum;

../trunk-jpl/src/m/boundaryconditions/SetMOLHOBC.py

@@ +1 @@
+import numpy as np
+
+def SetMOLHOBC(md):
+    """
+    SETMOLHOBC - Create the boundary conditions for stressbalance for MOLHO: VxBase, VyBase, VxShear, VyShear
+
+       Usage:
+          md = SetIceShelfBC(md, varargin)
+
+       Example:
+          md = SetIceShelfBC(md)
+
+    """
+
+    #node on Dirichlet (boundary and ~icefront)
+    md.stressbalance.spcvx_base = md.stressbalance.spcvx
+    md.stressbalance.spcvy_base = md.stressbalance.spcvy
+    md.stressbalance.spcvx_shear = np.nan * md.stressbalance.spcvx
+    md.stressbalance.spcvy_shear = np.nan * md.stressbalance.spcvy
+
+    return md

../trunk-jpl/src/m/classes/flowequation.py

@@ -18 +18 @@
         self.isSIA = 0
         self.isSSA = 0
         self.isL1L2 = 0
-        self.isMLHO = 0
+        self.isMOLHO = 0
         self.isHO = 0
         self.isFS = 0
         self.isNitscheBC = 0
@@ -45 +45 @@
         s += '{}\n'.format(fielddisplay(self, 'isSIA', "is the Shallow Ice Approximation (SIA) used?"))
         s += '{}\n'.format(fielddisplay(self, 'isSSA', "is the Shelfy-Stream Approximation (SSA) used?"))
         s += '{}\n'.format(fielddisplay(self, 'isL1L2', "are L1L2 equations used?"))
-        s += '{}\n'.format(fielddisplay(self, 'isMLHO', "are Mono-layer Higher-Order equations used?"))
+        s += '{}\n'.format(fielddisplay(self, 'isMOLHO', "are MOno-layer Higher-Order (MOLHO) equations used?"))
         s += '{}\n'.format(fielddisplay(self, 'isHO', "is the Higher-Order (HO) approximation used?"))
         s += '{}\n'.format(fielddisplay(self, 'isFS', "are the Full-FS (FS) equations used?"))
         s += '{}\n'.format(fielddisplay(self, 'isNitscheBC', "is weakly imposed condition used?"))
@@ -89 +89 @@
         md = checkfield(md, 'fieldname', 'flowequation.isSIA', 'numel', [1], 'values', [0, 1])
         md = checkfield(md, 'fieldname', 'flowequation.isSSA', 'numel', [1], 'values', [0, 1])
         md = checkfield(md, 'fieldname', 'flowequation.isL1L2', 'numel', [1], 'values', [0, 1])
-        md = checkfield(md, 'fieldname', 'flowequation.isMLHO', 'numel', [1], 'values', [0, 1])
+        md = checkfield(md, 'fieldname', 'flowequation.isMOLHO', 'numel', [1], 'values', [0, 1])
         md = checkfield(md, 'fieldname', 'flowequation.isHO', 'numel', [1], 'values', [0, 1])
         md = checkfield(md, 'fieldname', 'flowequation.isFS', 'numel', [1], 'values', [0, 1])
         md = checkfield(md, 'fieldname', 'flowequation.isNitscheBC', 'numel', [1], 'values', [0, 1])
@@ -120 +120 @@
         else:
             raise RuntimeError('Case not supported yet')
 
-        if not (self.isSIA or self.isSSA or self.isL1L2 or self.isMLHO or self.isHO or self.isFS):
+        if not (self.isSIA or self.isSSA or self.isL1L2 or self.isMOLHO or self.isHO or self.isFS):
             md.checkmessage("no element types set for this model")
         if 'StressbalanceSIAAnalysis' in analyses:
             if any(self.element_equation == 1):
@@ -133 +133 @@
         WriteData(fid, prefix, 'object', self, 'fieldname', 'isSIA', 'format', 'Boolean')
         WriteData(fid, prefix, 'object', self, 'fieldname', 'isSSA', 'format', 'Boolean')
         WriteData(fid, prefix, 'object', self, 'fieldname', 'isL1L2', 'format', 'Boolean')
-        WriteData(fid, prefix, 'object', self, 'fieldname', 'isMLHO', 'format', 'Boolean')
+        WriteData(fid, prefix, 'object', self, 'fieldname', 'isMOLHO', 'format', 'Boolean')
         WriteData(fid, prefix, 'object', self, 'fieldname', 'isHO', 'format', 'Boolean')
         WriteData(fid, prefix, 'object', self, 'fieldname', 'isFS', 'format', 'Boolean')
         WriteData(fid, prefix, 'object', self, 'fieldname', 'isNitscheBC', 'format', 'Boolean')

../trunk-jpl/src/m/classes/stressbalance.m

@@ -35 +35 @@
                         self.referential=project3d(md,'vector',self.referential,'type','node');
                         self.loadingforce=project3d(md,'vector',self.loadingforce,'type','node');
 
-                        % for MLHO
-                        if md.flowequation.isMLHO
+                        % for MOLHO
+                        if md.flowequation.isMOLHO
                                 self.spcvx_base=project3d(md,'vector',self.spcvx_base,'type','node');
                                 self.spcvy_base=project3d(md,'vector',self.spcvy_base,'type','node');
                                 self.spcvx_shear=project3d(md,'vector',self.spcvx_shear,'type','poly','degree',4);
@@ -133 +133 @@
                                 end
                                 md = checkfield(md,'fieldname','stressbalance.FSreconditioning','>',0);
                         end
-                        % CHECK THIS ONLY WORKS FOR MLHO
-                        if md.flowequation.isMLHO
+                        % CHECK THIS ONLY WORKS FOR MOLHO
+                        if md.flowequation.isMOLHO
                                 md = checkfield(md,'fieldname','stressbalance.spcvx_base','Inf',1,'timeseries',1);
                                 md = checkfield(md,'fieldname','stressbalance.spcvy_base','Inf',1,'timeseries',1);
                                 md = checkfield(md,'fieldname','stressbalance.spcvx_shear','Inf',1,'timeseries',1);
@@ -168 +168 @@
                         fielddisplay(self,'spcvy','y-axis velocity constraint (NaN means no constraint) [m/yr]');
                         fielddisplay(self,'spcvz','z-axis velocity constraint (NaN means no constraint) [m/yr]');
 
-                        disp(sprintf('\n      %s','MLHO boundary conditions:'));
+                        disp(sprintf('\n      %s','MOLHO boundary conditions:'));
                         fielddisplay(self,'spcvx_base','x-axis basal velocity constraint (NaN means no constraint) [m/yr]');
                         fielddisplay(self,'spcvy_base','y-axis basal velocity constraint (NaN means no constraint) [m/yr]');
                         fielddisplay(self,'spcvx_shear','x-axis shear velocity constraint (NaN means no constraint) [m/yr]');
@@ -225 +225 @@
                                 outputs      = [outputs defaultoutputs(self,md)]; %add defaults
                         end
                         WriteData(fid,prefix,'data',outputs,'name','md.stressbalance.requested_outputs','format','StringArray');
-                        % for MLHO
-                        if (md.flowequation.isMLHO)
+                        % for MOLHO
+                        if (md.flowequation.isMOLHO)
                                 WriteData(fid,prefix,'object',self,'class','stressbalance','fieldname','spcvx_base','format','DoubleMat','mattype',1,'scale',1./yts,'timeserieslength',md.mesh.numberofvertices+1,'yts',md.constants.yts);
                                 WriteData(fid,prefix,'object',self,'class','stressbalance','fieldname','spcvy_base','format','DoubleMat','mattype',1,'scale',1./yts,'timeserieslength',md.mesh.numberofvertices+1,'yts',md.constants.yts);
                                 WriteData(fid,prefix,'object',self,'class','stressbalance','fieldname','spcvx_shear','format','DoubleMat','mattype',1,'scale',1./yts,'timeserieslength',md.mesh.numberofvertices+1,'yts',md.constants.yts);

../trunk-jpl/src/m/parameterization/setflowequation.m

@@ -2 +2 @@
 %SETFLOWEQUATION - associate a solution type to each element
 %
 %   This routine works like plotmodel: it works with an even number of inputs
-%   'SIA','SSA','L1L2','MLHO','HO','FS' and 'fill' are the possible options
+%   'SIA','SSA','L1L2','MOLHO','HO','FS' and 'fill' are the possible options
 %   that must be followed by the corresponding exp file or flags list
 %   It can either be a domain file (argus type, .exp extension), or an array of element flags. 
 %   If user wants every element outside the domain to be 
@@ -10 +10 @@
 %   an empty string '' will be considered as an empty domain
 %   a string 'all' will be considered as the entire domain
 %   You can specify the type of coupling, 'penalties' or 'tiling', to use with the input 'coupling'
-%   NB: L1L2 and MLHO cannot currently be coupled to any other ice flow model
+%   NB: L1L2 and MOLHO cannot currently be coupled to any other ice flow model
 %
 %   Usage:
 %      md=setflowequation(md,varargin)
@@ -38 +38 @@
 SSAflag  = FlagElements(md,getfieldvalue(options,'SSA',''));
 HOflag   = FlagElements(md,getfieldvalue(options,'HO',''));
 L1L2flag = FlagElements(md,getfieldvalue(options,'L1L2',''));
-MLHOflag = FlagElements(md,getfieldvalue(options,'MLHO',''));
+MOLHOflag = FlagElements(md,getfieldvalue(options,'MOLHO',''));
 FSflag   = FlagElements(md,getfieldvalue(options,'FS',''));
 filltype = getfieldvalue(options,'fill','none');
 displayunused(options);
@@ -53 +53 @@
 end
 
 %check that each element has at least one flag
-if any(SIAflag+SSAflag+HOflag+L1L2flag+MLHOflag+FSflag==0),
+if any(SIAflag+SSAflag+HOflag+L1L2flag+MOLHOflag+FSflag==0),
         error('elements type not assigned, supported models are ''SIA'',''SSA'',''HO'' and ''FS''')
 end
 
 %check that each element has only one flag
-if any(SIAflag+SSAflag+HOflag+L1L2flag+MLHOflag+FSflag>1),
+if any(SIAflag+SSAflag+HOflag+L1L2flag+MOLHOflag+FSflag>1),
         disp('setflowequation.m: Warning: some elements have several types, higher order type is used for them')
         SIAflag(find(SIAflag & SSAflag))=0;
         SIAflag(find(SIAflag & HOflag))=0;
@@ -69 +69 @@
 if any(L1L2flag) & any(SIAflag | SSAflag | HOflag | FSflag)
         error('L1L2 cannot be coupled to any other model');
 end
-if any(MLHOflag) & any(SIAflag | SSAflag | HOflag | FSflag)
-        error('MLHO cannot be coupled to any other model');
+if any(MOLHOflag) & any(SIAflag | SSAflag | HOflag | FSflag)
+        error('MOLHO cannot be coupled to any other model');
 end
 
 %Check that no HO or FS for 2d mesh
@@ -90 +90 @@
 nodeonSSA=zeros(md.mesh.numberofvertices,1);  nodeonSSA(md.mesh.elements(find(SSAflag),:))=1;
 nodeonHO=zeros(md.mesh.numberofvertices,1);   nodeonHO(md.mesh.elements(find(HOflag),:))=1;
 nodeonL1L2=zeros(md.mesh.numberofvertices,1); nodeonL1L2(md.mesh.elements(find(L1L2flag),:))=1;
-nodeonMLHO=zeros(md.mesh.numberofvertices,1); nodeonMLHO(md.mesh.elements(find(MLHOflag),:))=1;
+nodeonMOLHO=zeros(md.mesh.numberofvertices,1); nodeonMOLHO(md.mesh.elements(find(MOLHOflag),:))=1;
 nodeonFS=zeros(md.mesh.numberofvertices,1);
 noneflag=zeros(md.mesh.numberofelements,1);
 
@@ -247 +247 @@
 md.flowequation.element_equation(find(SIAflag))=1;
 md.flowequation.element_equation(find(SSAflag))=2;
 md.flowequation.element_equation(find(L1L2flag))=3;
-md.flowequation.element_equation(find(MLHOflag))=4;
+md.flowequation.element_equation(find(MOLHOflag))=4;
 md.flowequation.element_equation(find(HOflag))=5;
 md.flowequation.element_equation(find(FSflag))=6;
 md.flowequation.element_equation(find(SSAHOflag))=7;
@@ -263 +263 @@
 md.flowequation.vertex_equation=zeros(md.mesh.numberofvertices,1);
 pos=find(nodeonSSA);  md.flowequation.vertex_equation(pos)=2;
 pos=find(nodeonL1L2); md.flowequation.vertex_equation(pos)=3;
-pos=find(nodeonMLHO); md.flowequation.vertex_equation(pos)=4;
+pos=find(nodeonMOLHO); md.flowequation.vertex_equation(pos)=4;
 pos=find(nodeonHO);   md.flowequation.vertex_equation(pos)=5;
 pos=find(nodeonFS);   md.flowequation.vertex_equation(pos)=6;
 %DO SIA LAST! Otherwise spcs might not be set up correctly (SIA should have priority)
@@ -286 +286 @@
 md.flowequation.isSIA  = double(any(md.flowequation.element_equation == 1));
 md.flowequation.isSSA  = double(any(md.flowequation.element_equation == 2));
 md.flowequation.isL1L2 = double(any(md.flowequation.element_equation == 3));
-md.flowequation.isMLHO = double(any(md.flowequation.element_equation == 4));
+md.flowequation.isMOLHO = double(any(md.flowequation.element_equation == 4));
 md.flowequation.isHO   = double(any(md.flowequation.element_equation == 5));
 md.flowequation.isFS   = double(any(md.flowequation.element_equation == 6));
 

../trunk-jpl/src/m/plot/plot_elementstype.m

@@ -36 +36 @@
                 patch( 'Faces', [C A D F],'Vertices', [x y z],'CData', i,'FaceColor','flat','EdgeColor',edgecolor);
         end
 end
-legend(p,'None','SIA','SSA','L1L2','MLHO','HO',...
+legend(p,'None','SIA','SSA','L1L2','MOLHO','HO',...
                 'SSAHO','FS','SSAFS','HOFS');
 
 %apply options

../trunk-jpl/test/NightlyRun/test128.py

@@ -1 @@
-#Test Name: SquareShelfConstrainedTranMLHO2d
+#Test Name: SquareShelfConstrainedTranMOLHO2d
 from model import *
-from SetMLHOBC import SetMLHOBC
+from SetMOLHOBC import SetMOLHOBC
 from socket import gethostname
 from triangle import *
 from setmask import *
@@ -12 +12 @@
 md = triangle(model(), '../Exp/Square.exp', 150000)
 md = setmask(md, 'all', '')
 md = parameterize(md, '../Par/SquareShelfConstrained.py')
-md = setflowequation(md, 'MLHO', 'all')
+md = setflowequation(md, 'MOLHO', 'all')
 md.cluster = generic('name', gethostname(), 'np', 3)
 md.transient.requested_outputs = ['IceVolume','VxSurface','VySurface','VxShear','VyShear','VxBase','VyBase']
 
-md = SetMLHOBC(md)
+md = SetMOLHOBC(md)
 md = solve(md, 'Transient')
 
 #Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test248.py

@@ -1 @@
-#Test Name: SquareShelfStressMLHO2d
+#Test Name: SquareShelfStressMOLHO2d
 from model import *
 from socket import gethostname
 from triangle import triangle
@@ -6 +6 @@
 from parameterize import parameterize
 from setflowequation import setflowequation
 from solve import solve
-from SetMLHOBC import SetMLHOBC
+from SetMOLHOBC import SetMOLHOBC
 
 md = triangle(model(), '../Exp/Square.exp', 150000)
 md = setmask(md, 'all', '')
 md = parameterize(md, '../Par/SquareShelf.py')
-md = setflowequation(md, 'MLHO', 'all')
+md = setflowequation(md, 'MOLHO', 'all')
 md.cluster = generic('name', gethostname(), 'np', 3)
 md.stressbalance.requested_outputs = ['default', 'VxSurface', 'VySurface', 'VxShear', 'VyShear', 'VxBase', 'VyBase']
-md = SetMLHOBC(md)
+md = SetMOLHOBC(md)
 md = solve(md, 'Stressbalance')
 
 #Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test254.py

@@ -1 @@
-#Test Name: SquareShelfConstrainedMLHO
+#Test Name: SquareShelfConstrainedMOLHO
 from model import *
 from socket import gethostname
 import numpy as np
@@ -8 +8 @@
 from setflowequation import setflowequation
 from paterson import paterson
 from solve import solve
-from SetMLHOBC import SetMLHOBC
+from SetMOLHOBC import SetMOLHOBC
 from generic import generic
 
 md = triangle(model(), '../Exp/Square.exp', 150000)
 md = setmask(md, 'all', '')
 md = parameterize(md, '../Par/SquareShelf.py')
-md = setflowequation(md, 'MLHO', 'all')
+md = setflowequation(md, 'MOLHO', 'all')
 md.cluster = generic('name', gethostname(), 'np', 3)
 
 # redo the parameter file for this special shelf.
@@ -70 +70 @@
 md.mask.ice_levelset = -1 + nodeonicefront
 
 md.stressbalance.requested_outputs = ['default', 'VySurface', 'VyShear', 'VyBase']
-md = SetMLHOBC(md)
+md = SetMOLHOBC(md)
 md = solve(md, 'Stressbalance')
 
 # create analytical solution: strain rate is constant = ((rho_ice * g * h) / 4B)^3 (Paterson, 4th Edition, page 292.

../trunk-jpl/src/m/parameterization/setflowequation.py

@@ -8 +8 @@
     """SETFLOWEQUATION - associate a solution type to each element
 
     This routine works like plotmodel: it works with an even number of inputs
-    'SIA', 'SSA', 'HO', 'L1L2', 'MLHO', 'FS' and 'fill' are the possible 
+    'SIA', 'SSA', 'HO', 'L1L2', 'MOLHO', 'FS' and 'fill' are the possible 
     options that must be followed by the corresponding exp file or flags list. 
     It can either be a domain file (argus type, .exp extension), or an array of 
     element flags.
@@ -44 +44 @@
     SSAflag = FlagElements(md, options.getfieldvalue('SSA', ''))
     HOflag = FlagElements(md, options.getfieldvalue('HO', ''))
     L1L2flag = FlagElements(md, options.getfieldvalue('L1L2', ''))
-    MLHOflag = FlagElements(md, options.getfieldvalue('MLHO', ''))
+    MOLHOflag = FlagElements(md, options.getfieldvalue('MOLHO', ''))
     FSflag = FlagElements(md, options.getfieldvalue('FS', ''))
     filltype = options.getfieldvalue('fill', 'none')
 
@@ -56 +56 @@
     elif 'HO' in filltype:
         HOflag = ~SIAflag & ~SSAflag & ~FSflag
     #check that each element has at least one flag
-    if not any(SIAflag + SSAflag + L1L2flag + MLHOflag + HOflag + FSflag):
+    if not any(SIAflag + SSAflag + L1L2flag + MOLHOflag + HOflag + FSflag):
         raise TypeError("elements type not assigned, supported models are 'SIA', 'SSA', 'HO' and 'FS'")
 
     #check that each element has only one flag
-    if any(SIAflag + SSAflag + L1L2flag + MLHOflag + HOflag + FSflag > 1):
+    if any(SIAflag + SSAflag + L1L2flag + MOLHOflag + HOflag + FSflag > 1):
         print('Warning: setflowequation.py: some elements have several types, higher order type is used for them')
         SIAflag[np.where(np.logical_and(SIAflag, SSAflag))] = False
         SIAflag[np.where(np.logical_and(SIAflag, HOflag))] = False
         SSAflag[np.where(np.logical_and(SSAflag, HOflag))] = False
 
-        #check that L1L2 and MLHO is not coupled to any other model for now
+        #check that L1L2 and MOLHO is not coupled to any other model for now
         if any(L1L2flag) and any(SIAflag + SSAflag + HOflag + FSflag):
             raise TypeError('L1L2 cannot be coupled to any other model')
-        if any(MLHOflag) and any(SIAflag + SSAflag + HOflag + FSflag):
-            raise TypeError('MLHO cannot be coupled to any other model')
+        if any(MOLHOflag) and any(SIAflag + SSAflag + HOflag + FSflag):
+            raise TypeError('MOLHO cannot be coupled to any other model')
 
         #Check that no HO or FS for 2d mesh
         if md.mesh.domaintype == '2Dhorizontal':
@@ -88 +88 @@
     nodeonSSA[md.mesh.elements[np.where(SSAflag), :] - 1] = True
     nodeonL1L2 = np.zeros(md.mesh.numberofvertices, bool)
     nodeonL1L2[md.mesh.elements[np.where(L1L2flag), :] - 1] = True
-    nodeonMLHO = np.zeros(md.mesh.numberofvertices, bool)
-    nodeonMLHO[md.mesh.elements[np.where(MLHOflag), :] - 1] = True
+    nodeonMOLHO = np.zeros(md.mesh.numberofvertices, bool)
+    nodeonMOLHO[md.mesh.elements[np.where(MOLHOflag), :] - 1] = True
     nodeonHO = np.zeros(md.mesh.numberofvertices, bool)
     nodeonHO[md.mesh.elements[np.where(HOflag), :] - 1] = True
     nodeonFS = np.zeros(md.mesh.numberofvertices, bool)
@@ -240 +240 @@
     md.flowequation.element_equation[np.where(SIAflag)] = 1
     md.flowequation.element_equation[np.where(SSAflag)] = 2
     md.flowequation.element_equation[np.where(L1L2flag)] = 3
-    md.flowequation.element_equation[np.where(MLHOflag)] = 4
+    md.flowequation.element_equation[np.where(MOLHOflag)] = 4
     md.flowequation.element_equation[np.where(HOflag)] = 5
     md.flowequation.element_equation[np.where(FSflag)] = 6
     md.flowequation.element_equation[np.where(SSAHOflag)] = 7
@@ -258 +258 @@
     md.flowequation.vertex_equation[pos] = 2
     pos = np.where(nodeonL1L2)
     md.flowequation.vertex_equation[pos] = 3
-    pos = np.where(nodeonMLHO)
+    pos = np.where(nodeonMOLHO)
     md.flowequation.vertex_equation[pos] = 4
     pos = np.where(nodeonHO)
     md.flowequation.vertex_equation[pos] = 5
@@ -282 +282 @@
     md.flowequation.isSIA = any(md.flowequation.element_equation == 1)
     md.flowequation.isSSA = any(md.flowequation.element_equation == 2)
     md.flowequation.isL1L2= any(md.flowequation.element_equation == 3)
-    md.flowequation.isMLHO= any(md.flowequation.element_equation == 4)
+    md.flowequation.isMOLHO= any(md.flowequation.element_equation == 4)
     md.flowequation.isHO = any(md.flowequation.element_equation == 5)
     md.flowequation.isFS = any(md.flowequation.element_equation == 6)
 

../trunk-jpl/test/NightlyRun/test128.m

@@ -1 @@
-%Test Name: SquareShelfConstrainedTranMLHO2d
+%Test Name: SquareShelfConstrainedTranMOLHO2d
 md=triangle(model(),'../Exp/Square.exp',150000.);
 md=setmask(md,'all','');
 md=parameterize(md,'../Par/SquareShelfConstrained.par');
-md=setflowequation(md,'MLHO','all');
+md=setflowequation(md,'MOLHO','all');
 md.cluster=generic('name',oshostname(),'np',3);
 md.transient.requested_outputs={'IceVolume','VxShear','VyShear','VxBase','VyBase','VxSurface','VySurface'};
 
-md=SetMLHOBC(md);
+md=SetMOLHOBC(md);
 md=solve(md,'Transient');
 
 %Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test248.m

@@ -1 @@
-%Test Name: SquareShelfStressMLHO2d
+%Test Name: SquareShelfStressMOLHO2d
 md=triangle(model(),'../Exp/Square.exp',150000.);
 md=setmask(md,'all','');
 md=parameterize(md,'../Par/SquareShelf.par');
-md=setflowequation(md,'MLHO','all');
+md=setflowequation(md,'MOLHO','all');
 md.cluster=generic('name',oshostname(),'np',3);
 md.stressbalance.requested_outputs={'default','VxSurface','VySurface','VxShear','VyShear','VxBase','VyBase'};
-md=SetMLHOBC(md);
+md=SetMOLHOBC(md);
 md=solve(md,'Stressbalance');
 
 %Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test254.m

@@ -1 @@
-%Test Name: SquareShelfConstrainedMLHO
+%Test Name: SquareShelfConstrainedMOLHO
 md=triangle(model(),'../Exp/Square.exp',150000.);
 md=setmask(md,'all','');
 md=parameterize(md,'../Par/SquareShelf.par');
-md=setflowequation(md,'MLHO','all');
+md=setflowequation(md,'MOLHO','all');
 md.cluster=generic('name',oshostname(),'np',3);
 
 %redo the parameter file for this special shelf. 
@@ -57 +57 @@
 md.mask.ice_levelset=-1+nodeonicefront;
 
 md.stressbalance.requested_outputs={'default','VySurface','VyShear','VyBase'};
-md=SetMLHOBC(md);
+md=SetMOLHOBC(md);
 md=solve(md,'Stressbalance');
 
 %create analytical solution: strain rate is constant = ((rho_ice*g*h)/4B)^3 (Paterson, 4th Edition, page 292.

../trunk-jpl/test/NightlyRun/test256.m

@@ -1 @@
-%Test Name: SquareShelfStressMLHO2dTransientIncrNonHydro
+%Test Name: SquareShelfStressMOLHO2dTransientIncrNonHydro
 md=triangle(model(),'../Exp/Square.exp',150000.);
 md=setmask(md,'all','');
 md=parameterize(md,'../Par/SquareShelf.par');
-md=setflowequation(md,'MLHO','all');
+md=setflowequation(md,'MOLHO','all');
 md.geometry.base=md.geometry.base+50.; md.geometry.surface=md.geometry.surface+50.;
 md.cluster=generic('name',oshostname(),'np',1);
 md.masstransport.hydrostatic_adjustment='Incremental';
-md=SetMLHOBC(md);
+md=SetMOLHOBC(md);
 md=solve(md,'Transient');
 
 %Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test127.py

@@ -1 @@
-#Test Name: SquareShelfConstrainedStressMLHO2d
+#Test Name: SquareShelfConstrainedStressMOLHO2d
 from model import *
 from socket import gethostname
 from triangle import triangle
@@ -8 +8 @@
 from solve import solve
 from massfluxatgate import massfluxatgate
 from generic import generic
-from SetMLHOBC import SetMLHOBC
+from SetMOLHOBC import SetMOLHOBC
 
 md = triangle(model(), '../Exp/Square.exp', 50000)
 md = setmask(md, 'all', '')
 md = parameterize(md, '../Par/SquareShelfConstrained.py')
-md = setflowequation(md, 'MLHO', 'all')
+md = setflowequation(md, 'MOLHO', 'all')
 md.cluster = generic('name', gethostname(), 'np', 2)
 #outputs
-#FIXME compute the stress components for MLHO
+#FIXME compute the stress components for MOLHO
 md.stressbalance.requested_outputs = ['default', 'VxSurface', 'VySurface', 'VxShear', 'VyShear', 'VxBase', 'VyBase', 'MassFlux1', 'MassFlux2', 'MassFlux3', 'MassFlux4', 'MassFlux5', 'MassFlux6']
 #md.stressbalance.requested_outputs = ['default', 'DeviatoricStressxx', 'DeviatoricStressyy', 'DeviatoricStressxy', 'MassFlux1', 'MassFlux2', 'MassFlux3', 'MassFlux4', 'MassFlux5', 'MassFlux6']
 md.outputdefinition.definitions = [massfluxatgate('name', 'MassFlux1', 'profilename', '../Exp/MassFlux1.exp', 'definitionstring', 'Outputdefinition1'),
@@ -25 +25 @@
                                    massfluxatgate('name', 'MassFlux4', 'profilename', '../Exp/MassFlux4.exp', 'definitionstring', 'Outputdefinition4'),
                                    massfluxatgate('name', 'MassFlux5', 'profilename', '../Exp/MassFlux5.exp', 'definitionstring', 'Outputdefinition5'),
                                    massfluxatgate('name', 'MassFlux6', 'profilename', '../Exp/MassFlux6.exp', 'definitionstring', 'Outputdefinition6')]
-md = SetMLHOBC(md)
+md = SetMOLHOBC(md)
 md = solve(md, 'Stressbalance')
 
 #Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test129.py

@@ -1 @@
-#Test Name: SquareShelfConstrainedRestartTranMLHO2d
+#Test Name: SquareShelfConstrainedRestartTranMOLHO2d
 from model import *
 from socket import gethostname
-from SetMLHOBC import SetMLHOBC
+from SetMOLHOBC import SetMOLHOBC
 from triangle import *
 from setmask import *
 from parameterize import *
@@ -13 +13 @@
 md = triangle(model(), '../Exp/Square.exp', 150000.)
 md = setmask(md, 'all', '')
 md = parameterize(md, '../Par/SquareShelfConstrained.py')
-md = setflowequation(md, 'MLHO', 'all')
+md = setflowequation(md, 'MOLHO', 'all')
 md.cluster = generic('name', gethostname(), 'np', 1)
 md.transient.requested_outputs = ['IceVolume', 'TotalSmb', 'VxShear','VyShear','VxBase','VyBase','VxSurface','VySurface']
 
@@ -23 +23 @@
 # time steps and resolution
 md.timestepping.final_time = 19
 md.settings.output_frequency = 2
-md = SetMLHOBC(md)
+md = SetMOLHOBC(md)
 
 md = solve(md, 'Transient')
 md2 = copy.deepcopy(md)

../trunk-jpl/test/NightlyRun/test249.py

@@ -1 @@
-#Test Name: SquareShelfTranMLHO2d
+#Test Name: SquareShelfTranMOLHO2d
 
 from model import *
 from socket import gethostname
@@ -8 +8 @@
 from setflowequation import setflowequation
 from solve import solve
 from generic import generic
-from SetMLHOBC import SetMLHOBC
+from SetMOLHOBC import SetMOLHOBC
 
 md = triangle(model(), '../Exp/Square.exp', 150000)
 md = setmask(md, 'all', '')
 md = parameterize(md, '../Par/SquareShelf.py')
 md.basalforcings.floatingice_melting_rate[:] = 1.
-md = setflowequation(md, 'MLHO', 'all')
+md = setflowequation(md, 'MOLHO', 'all')
 md.cluster = generic('name', gethostname(), 'np', 3)
 md.transient.requested_outputs = ['default', 'FloatingArea', 'GroundedArea', 'TotalFloatingBmb', 'TotalGroundedBmb']
-md = SetMLHOBC(md)
+md = SetMOLHOBC(md)
 md = solve(md, 'Transient')
 
 

../trunk-jpl/test/NightlyRun/test255.py

@@ -1 @@
-#Test Name: SquareShelfStressMLHO2dTransientIncrHydro
+#Test Name: SquareShelfStressMOLHO2dTransientIncrHydro
 from model import *
 from socket import gethostname
 from triangle import triangle
@@ -6 +6 @@
 from parameterize import parameterize
 from setflowequation import setflowequation
 from solve import solve
-from SetMLHOBC import SetMLHOBC
+from SetMOLHOBC import SetMOLHOBC
 
 md = triangle(model(), '../Exp/Square.exp', 150000.)
 md = setmask(md, 'all', '')
 md = parameterize(md, '../Par/SquareShelf.py')
-md = setflowequation(md, 'MLHO', 'all')
+md = setflowequation(md, 'MOLHO', 'all')
 md.cluster = generic('name', gethostname(), 'np', 3)
 md.masstransport.hydrostatic_adjustment = 'Incremental'
-md = SetMLHOBC(md)
+md = SetMOLHOBC(md)
 md = solve(md, 'Transient')
 
 #Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test256.py

@@ -1 @@
-#Test Name: SquareShelfStressMLHO2dTransientIncrNonHydro
+#Test Name: SquareShelfStressMOLHO2dTransientIncrNonHydro
 from model import *
 from socket import gethostname
 from triangle import triangle
@@ -6 +6 @@
 from parameterize import parameterize
 from setflowequation import setflowequation
 from solve import solve
-from SetMLHOBC import SetMLHOBC
+from SetMOLHOBC import SetMOLHOBC
 
 
 md = triangle(model(), '../Exp/Square.exp', 150000.)
 md = setmask(md, 'all', '')
 md = parameterize(md, '../Par/SquareShelf.py')
-md = setflowequation(md, 'MLHO', 'all')
+md = setflowequation(md, 'MOLHO', 'all')
 md.geometry.base = md.geometry.base + 50.
 md.geometry.surface = md.geometry.surface + 50.
 md.cluster = generic('name', gethostname(), 'np', 1)
 md.masstransport.hydrostatic_adjustment = 'Incremental'
-md = SetMLHOBC(md)
+md = SetMOLHOBC(md)
 md = solve(md, 'Transient')
 
 #Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test332.py

@@ -1 @@
-#Test Name: SquareSheetConstrainedTranMLHO2d
+#Test Name: SquareSheetConstrainedTranMOLHO2d
 from model import *
 from socket import gethostname
 from triangle import triangle
@@ -6 +6 @@
 from parameterize import parameterize
 from setflowequation import setflowequation
 from solve import solve
-from SetMLHOBC import SetMLHOBC
+from SetMOLHOBC import SetMOLHOBC
 
 md = triangle(model(), '../Exp/Square.exp', 150000.)
 md = setmask(md, '', '')
 md = parameterize(md, '../Par/SquareSheetConstrained.py')
-md = setflowequation(md, 'MLHO', 'all')
+md = setflowequation(md, 'MOLHO', 'all')
 md.cluster = generic('name', gethostname(), 'np', 3)
-md = SetMLHOBC(md)
+md = SetMOLHOBC(md)
 md = solve(md, 'Transient')
 
 #Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test334.py

@@ -1 @@
-#Test Name: SquareShelfCMBMLHO
+#Test Name: SquareShelfCMBMOLHO
 from model import *
 from socket import gethostname
 import numpy as np
@@ -8 +8 @@
 from setflowequation import setflowequation
 from solve import solve
 from generic import generic
-from SetMLHOBC import SetMLHOBC
+from SetMOLHOBC import SetMOLHOBC
 
 md = triangle(model(), '../Exp/Square.exp', 200000)
 md = setmask(md, 'all', '')
 md = parameterize(md, '../Par/SquareShelf.py')
-md = setflowequation(md, 'MLHO', 'all')
+md = setflowequation(md, 'MOLHO', 'all')
 
 # control parameters
 
@@ -32 +32 @@
 
 
 md.cluster = generic('name', gethostname(), 'np', 3)
-md = SetMLHOBC(md)
+md = SetMOLHOBC(md)
 md = solve(md, 'Stressbalance')
 
 

../trunk-jpl/test/NightlyRun/test446.py

@@ -1 @@
-#Test Name: SquareSheetShelfStressMLHO2d
+#Test Name: SquareSheetShelfStressMOLHO2d
 from model import *
 from socket import gethostname
 from triangle import triangle
@@ -6 +6 @@
 from parameterize import parameterize
 from setflowequation import setflowequation
 from solve import solve
-from SetMLHOBC import SetMLHOBC
+from SetMOLHOBC import SetMOLHOBC
 
 md = triangle(model(), '../Exp/Square.exp', 150000.)
 md = setmask(md, '../Exp/SquareShelf.exp', '')
 md = parameterize(md, '../Par/SquareSheetShelf.py')
-md = setflowequation(md, 'MLHO', 'all')
+md = setflowequation(md, 'MOLHO', 'all')
 md.cluster = generic('name', gethostname(), 'np', 3)
 md.stressbalance.requested_outputs = ['default', 'VxSurface', 'VySurface', 'VxShear', 'VyShear', 'VxBase', 'VyBase']
-md = SetMLHOBC(md);
+md = SetMOLHOBC(md);
 md = solve(md, 'Stressbalance')
 
 #Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test448.py

@@ -1 @@
-#Test Name: RoundSheetShelfGLMigrationMLHO2d
+#Test Name: RoundSheetShelfGLMigrationMOLHO2d
 import numpy as np
 from model import *
 from socket import gethostname
@@ -7 +7 @@
 from parameterize import parameterize
 from setflowequation import setflowequation
 from solve import solve
-from SetMLHOBC import SetMLHOBC
+from SetMOLHOBC import SetMOLHOBC
 
 
 radius = 1.e6
@@ -27 +27 @@
 flags[pos] = 1
 md = setmask(md, flags, '')
 md = parameterize(md, '../Par/RoundSheetShelf.py')
-md = setflowequation(md, 'MLHO', 'all')
+md = setflowequation(md, 'MOLHO', 'all')
 md.cluster = generic('name', gethostname(), 'np', 3)
 
 md.transient.isthermal = False
@@ -38 +38 @@
 
 #test different grounding line dynamics.
 md.groundingline.migration = 'AggressiveMigration'
-md = SetMLHOBC(md)
+md = SetMOLHOBC(md)
 md = solve(md, 'Transient')
 element_on_iceshelf_agressive = md.results.TransientSolution[0].MaskOceanLevelset
 vel_agressive = md.results.TransientSolution[0].Vel

../trunk-jpl/test/NightlyRun/test518.py

@@ -1 +1 @@
 
-#Test Name: PigStressMLHO2d
+#Test Name: PigStressMOLHO2d
 from model import *
 from socket import gethostname
 from triangle import triangle
@@ -7 +7 @@
 from parameterize import parameterize
 from setflowequation import setflowequation
 from solve import solve
-from SetMLHOBC import SetMLHOBC
+from SetMOLHOBC import SetMOLHOBC
 
 
 md = triangle(model(), '../Exp/Pig.exp', 20000.)
 md = setmask(md, '../Exp/PigShelves.exp', '../Exp/PigIslands.exp')
 md = parameterize(md, '../Par/Pig.py')
-md = setflowequation(md, 'MLHO', 'all')
+md = setflowequation(md, 'MOLHO', 'all')
 md.cluster = generic('name', gethostname(), 'np', 3)
 md.stressbalance.requested_outputs = ['default', 'VxSurface', 'VySurface', 'VxShear', 'VyShear', 'VxBase', 'VyBase']
-md = SetMLHOBC(md)
+md = SetMOLHOBC(md)
 md = solve(md, 'Stressbalance')
 
 #Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test332.m

@@ -1 @@
-%Test Name: SquareSheetConstrainedTranMLHO2d
+%Test Name: SquareSheetConstrainedTranMOLHO2d
 md=triangle(model(),'../Exp/Square.exp',150000.);
 md=setmask(md,'','');
 md=parameterize(md,'../Par/SquareSheetConstrained.par');
-md=setflowequation(md,'MLHO','all');
+md=setflowequation(md,'MOLHO','all');
 md.cluster=generic('name',oshostname(),'np',3);
-md=SetMLHOBC(md);
+md=SetMOLHOBC(md);
 md=solve(md,'Transient');
 
 %Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test334.m

@@ -1 @@
-%Test Name: SquareShelfCMBMLHO
+%Test Name: SquareShelfCMBMOLHO
 md=triangle(model(),'../Exp/Square.exp',200000.);
 md=setmask(md,'all','');
 md=parameterize(md,'../Par/SquareShelf.par');
-md=setflowequation(md,'MLHO','all');
+md=setflowequation(md,'MOLHO','all');
 
 %control parameters
 md.inversion.iscontrol=1;
@@ -19 +19 @@
 md.verbose.control=true;
 
 md.cluster=generic('name',oshostname(),'np',3);
-md=SetMLHOBC(md);
+md=SetMOLHOBC(md);
 md=solve(md,'Stressbalance');
 
 %Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test446.m

@@ -1 @@
-%Test Name: SquareSheetShelfStressMLHO2d
+%Test Name: SquareSheetShelfStressMOLHO2d
 md=triangle(model(),'../Exp/Square.exp',150000.);
 md=setmask(md,'../Exp/SquareShelf.exp','');
 md=parameterize(md,'../Par/SquareSheetShelf.par');
-md=setflowequation(md,'MLHO','all');
+md=setflowequation(md,'MOLHO','all');
 md.cluster=generic('name',oshostname(),'np',3);
 md.stressbalance.requested_outputs={'default','VxSurface','VySurface','VxShear','VyShear','VxBase','VyBase'};
-md=SetMLHOBC(md);
+md=SetMOLHOBC(md);
 md=solve(md,'Stressbalance');
 
 %Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test448.m

@@ -1 @@
-%Test Name: RoundSheetShelfGLMigrationMLHO2d
+%Test Name: RoundSheetShelfGLMigrationMOLHO2d
 radius=1.e6;
 shelfextent=2.e5;
 
@@ -15 +15 @@
 flags(pos)=1;
 md=setmask(md,flags,''); 
 md=parameterize(md,'../Par/RoundSheetShelf.par');
-md=setflowequation(md,'MLHO','all');
+md=setflowequation(md,'MOLHO','all');
 md.cluster=generic('name',oshostname(),'np',3);
 
 md.transient.isthermal=0;
@@ -24 +24 @@
 md.transient.isstressbalance=1;
 md.transient.isgroundingline=1;
 
-md=SetMLHOBC(md);
+md=SetMOLHOBC(md);
 %test different grounding line dynamics.
 md.groundingline.migration='AggressiveMigration';
 md=solve(md,'Transient');

../trunk-jpl/test/NightlyRun/test518.m

@@ -1 @@
-%Test Name: PigStressMLHO2d
+%Test Name: PigStressMOLHO2d
 md=triangle(model(),'../Exp/Pig.exp',20000.);
 md=setmask(md,'../Exp/PigShelves.exp','../Exp/PigIslands.exp');
 md=parameterize(md,'../Par/Pig.par');
-md=setflowequation(md,'MLHO','all');
+md=setflowequation(md,'MOLHO','all');
 md.cluster=generic('name',oshostname(),'np',3);
 md.stressbalance.requested_outputs={'default','VxSurface','VySurface','VxShear','VyShear','VxBase','VyBase'};
-md=SetMLHOBC(md);
+md=SetMOLHOBC(md);
 md=solve(md,'Stressbalance');
 
 %Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test810.m

@@ -1 @@
-%Test Name: ValleyGlacierLevelsetMLHO2d
+%Test Name: ValleyGlacierLevelsetMOLHO2d
 md=triangle(model(),'../Exp/Square.exp',50000);
 md=setmask(md,'','');
 md=parameterize(md,'../Par/ValleyGlacierShelf.par');
-md=setflowequation(md,'MLHO','all');
+md=setflowequation(md,'MOLHO','all');
 md.levelset.stabilization=2;
 md.cluster=generic('name',oshostname(),'np',3);
 
@@ -14 +14 @@
 md.transient.isthermal=0;
 md.transient.isgroundingline=1;
 
-md=SetMLHOBC(md);
+md=SetMOLHOBC(md);
 md=solve(md,'Transient');
 
 %Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test330.py

@@ -1 @@
-#Test Name: SquareSheetConstrainedStressMLHO2d
+#Test Name: SquareSheetConstrainedStressMOLHO2d
 from model import *
 from socket import gethostname
 from triangle import triangle
@@ -6 +6 @@
 from parameterize import parameterize
 from setflowequation import setflowequation
 from solve import solve
-from SetMLHOBC import SetMLHOBC
+from SetMOLHOBC import SetMOLHOBC
 
 
 md = triangle(model(), '../Exp/Square.exp', 150000.)
 md = setmask(md, '', '')
 md = parameterize(md, '../Par/SquareSheetConstrained.py')
-md = setflowequation(md, 'MLHO', 'all')
+md = setflowequation(md, 'MOLHO', 'all')
 md.cluster = generic('name', gethostname(), 'np', 3)
 md.stressbalance.requested_outputs = ['default', 'VxSurface', 'VySurface', 'VxShear', 'VyShear', 'VxBase', 'VyBase']
-md = SetMLHOBC(md)
+md = SetMOLHOBC(md)
 md = solve(md, 'Stressbalance')
 
 #Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test333.py

@@ -1 @@
-#Test Name: SquareSheetConstrainedCMDragMLHO
+#Test Name: SquareSheetConstrainedCMDragMOLHO
 import numpy as np
 from model import *
 from socket import gethostname
@@ -7 +7 @@
 from parameterize import parameterize
 from setflowequation import setflowequation
 from solve import solve
-from SetMLHOBC import SetMLHOBC
+from SetMOLHOBC import SetMOLHOBC
 
 
 md = triangle(model(), '../Exp/Square.exp', 200000.)
 md = setmask(md, '', '')
 md = parameterize(md, '../Par/SquareSheetConstrained.py')
-md = setflowequation(md, 'MLHO', 'all')
+md = setflowequation(md, 'MOLHO', 'all')
 
 #control parameters
 
@@ -32 +32 @@
 md.inversion.vy_obs = md.initialization.vy
 
 md.cluster = generic('name', gethostname(), 'np', 3)
-md = SetMLHOBC(md)
+md = SetMOLHOBC(md)
 md = solve(md, 'Stressbalance')
 
 #Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test335.py

@@ -1 @@
-#Test Name: SquareSheetConstrainedStressMLHO2d
+#Test Name: SquareSheetConstrainedStressMOLHO2d
 from model import *
 from socket import gethostname
 from triangle import triangle
@@ -6 +6 @@
 from parameterize import parameterize
 from setflowequation import setflowequation
 from solve import solve
-from SetMLHOBC import SetMLHOBC
+from SetMOLHOBC import SetMOLHOBC
 
 
 md = triangle(model(), '../Exp/Square.exp', 200000.)
 md = setmask(md, '', '')
 md = parameterize(md, '../Par/SquareSheetConstrained.py')
-md = setflowequation(md, 'MLHO', 'all')
-md = SetMLHOBC(md)
+md = setflowequation(md, 'MOLHO', 'all')
+md = SetMOLHOBC(md)
 md.extrude(5, 1.)
 md.cluster = generic('name', gethostname(), 'np', 3)
 md.stressbalance.requested_outputs = ['default', 'VxSurface', 'VySurface', 'VxShear', 'VyShear', 'VxBase', 'VyBase']

../trunk-jpl/test/NightlyRun/test447.py

@@ -1 @@
-#Test Name: SquareSheetShelfTranMLHO2d
+#Test Name: SquareSheetShelfTranMOLHO2d
 import numpy as np
 from model import *
 from socket import gethostname
@@ -7 +7 @@
 from parameterize import parameterize
 from setflowequation import setflowequation
 from solve import solve
-from SetMLHOBC import SetMLHOBC
+from SetMOLHOBC import SetMOLHOBC
 from frictioncoulomb import frictioncoulomb
 from generic import generic
 
@@ -24 +24 @@
 md.friction.coefficientcoulomb = 0.02 * np.ones(md.mesh.numberofvertices)
 md.transient.isthermal = False
 md.transient.isgroundingline = True
-md = setflowequation(md, 'MLHO', 'all')
+md = setflowequation(md, 'MOLHO', 'all')
 md.cluster = generic('name', gethostname(), 'np', 3)
 md.transient.requested_outputs = ['default', 'GroundedArea', 'FloatingArea', 'TotalFloatingBmb', 'TotalGroundedBmb', 'TotalSmb']
-md = SetMLHOBC(md)
+md = SetMOLHOBC(md)
 md = solve(md, 'Transient')
 
 #Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test449.py

@@ -1 @@
-#Test Name: MISMIP3DMLHO
+#Test Name: MISMIP3DMOLHO
 import numpy as np
 from model import *
 from socket import gethostname
@@ -7 +7 @@
 from parameterize import parameterize
 from setflowequation import setflowequation
 from solve import solve
-from SetMLHOBC import SetMLHOBC
+from SetMOLHOBC import SetMOLHOBC
 from generic import generic
 
 
@@ -23 +23 @@
 pos = np.where(md.mask.ocean_levelset >= 0.)
 md.geometry.base[pos] = md.geometry.bed[pos]
 md.geometry.surface = md.geometry.base + md.geometry.thickness
-md = setflowequation(md, 'MLHO', 'all')
+md = setflowequation(md, 'MOLHO', 'all')
 
 #Boundary conditions:
 md.mask.ice_levelset = -np.ones((md.mesh.numberofvertices, ))
@@ -58 +58 @@
 md.timestepping.time_step = 10
 
 md.cluster = generic('name', gethostname(), 'np', 3)
-md = SetMLHOBC(md)
+md = SetMOLHOBC(md)
 md = solve(md, 'Transient')
 #print md.results.TransientSolution[0].BasalforcingsFloatingiceMeltingRate
 #print md.results.TransientSolution[1].BasalforcingsFloatingiceMeltingRate

../trunk-jpl/test/NightlyRun/test519.py

@@ -1 @@
-#Test Name: PigTranMLHO2d
+#Test Name: PigTranMOLHO2d
 import numpy as np
 from model import *
 from socket import gethostname
@@ -7 +7 @@
 from parameterize import parameterize
 from setflowequation import setflowequation
 from solve import solve
-from SetMLHOBC import SetMLHOBC
+from SetMOLHOBC import SetMOLHOBC
 from generic import generic
 
 md = triangle(model(), '../Exp/Pig.exp', 20000.)
 md = setmask(md, '../Exp/PigShelves.exp', '../Exp/PigIslands.exp')
 md = parameterize(md, '../Par/Pig.py')
-md = setflowequation(md, 'MLHO', 'all')
+md = setflowequation(md, 'MOLHO', 'all')
 md.mesh.scale_factor = 0.9 * np.ones((md.mesh.numberofvertices))
 md.transient.requested_outputs = ['default', 'IceVolume', 'IceVolumeScaled', 'GroundedArea', 'GroundedAreaScaled', 'FloatingArea', 'FloatingAreaScaled', 'TotalSmb', 'TotalSmbScaled', 'TotalFloatingBmb', 'TotalFloatingBmbScaled']
 md.cluster = generic('name', gethostname(), 'np', 3)
-md = SetMLHOBC(md)
+md = SetMOLHOBC(md)
 md = solve(md, 'Transient')
 
 # Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test810.py

@@ -1 @@
-#Test Name: ValleyGlacierLevelsetMLHO2d
+#Test Name: ValleyGlacierLevelsetMOLHO2d
 from model import *
 from socket import gethostname
 from triangle import triangle
@@ -6 +6 @@
 from parameterize import parameterize
 from setflowequation import setflowequation
 from solve import solve
-from SetMLHOBC import SetMLHOBC
+from SetMOLHOBC import SetMOLHOBC
 
 md = triangle(model(), '../Exp/Square.exp', 50000)
 md = setmask(md, '', '')
 md = parameterize(md, '../Par/ValleyGlacierShelf.py')
 md.levelset.stabilization = 2
-md = setflowequation(md, 'MLHO', 'all')
+md = setflowequation(md, 'MOLHO', 'all')
 md.cluster = generic('name', gethostname(), 'np', 3)
 
 #Transient
@@ -22 +22 @@
 md.transient.issmb = True
 md.transient.isthermal = False
 md.transient.isgroundingline = True
-md = SetMLHOBC(md)
+md = SetMOLHOBC(md)
 
 md = solve(md, 'Transient')
 

../trunk-jpl/test/NightlyRun/test812.py

@@ -1 @@
-#Test Name: SquareShelfLevelsetCalvingMLHO2dLevermann
+#Test Name: SquareShelfLevelsetCalvingMOLHO2dLevermann
 import numpy as np
 from socket import gethostname
 from model import *
@@ -7 +7 @@
 from setmask import setmask
 from solve import solve
 from triangle import triangle
-from SetMLHOBC import SetMLHOBC
+from SetMOLHOBC import SetMOLHOBC
 
 md = triangle(model(), '../Exp/Square.exp', 50000.)
 md = setmask(md, 'all', '')
 md = parameterize(md, '../Par/SquareShelf.py')
-md = setflowequation(md, 'MLHO', 'all')
+md = setflowequation(md, 'MOLHO', 'all')
 md.cluster = generic('name', gethostname(), 'np', 3)
 
 x = md.mesh.x
@@ -44 +44 @@
 
 md.transient.requested_outputs = ['default', 'StrainRateparallel', 'StrainRateperpendicular', 'Calvingratex', 'Calvingratey', 'CalvingCalvingrate']
 
-md = SetMLHOBC(md);
+md = SetMOLHOBC(md);
 md = solve(md, 'Transient')
 
 #Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test812.m

@@ -1 @@
-%Test Name: SquareShelfLevelsetCalvingMLHO2dLevermann
+%Test Name: SquareShelfLevelsetCalvingMOLHO2dLevermann
 md=triangle(model(),'../Exp/Square.exp',50000.);
 md=setmask(md,'all','');
 md=parameterize(md,'../Par/SquareShelf.par');
-md=setflowequation(md,'MLHO','all');
+md=setflowequation(md,'MOLHO','all');
 md.cluster=generic('name',oshostname(),'np',3);
 
 %Do not kill ice bergs as all is floating
@@ -34 +34 @@
 
 md.transient.requested_outputs={'default','StrainRateparallel','StrainRateperpendicular','Calvingratex','Calvingratey','CalvingCalvingrate'};
 
-md=SetMLHOBC(md);
+md=SetMOLHOBC(md);
 md=solve(md,'Transient');
 
 %Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test811.py

@@ -1 @@
-#Test Name: ValleyGlacierLevelsetCalvingMLHO2d
+#Test Name: ValleyGlacierLevelsetCalvingMOLHO2d
 import numpy as np
 from model import *
 from socket import gethostname
@@ -7 +7 @@
 from parameterize import parameterize
 from setflowequation import setflowequation
 from solve import solve
-from SetMLHOBC import SetMLHOBC
+from SetMOLHOBC import SetMOLHOBC
 from generic import generic
 
 md = triangle(model(), '../Exp/Square.exp', 50000)
 md = setmask(md, '', '')
 md = parameterize(md, '../Par/ValleyGlacierShelf.py')
-md = setflowequation(md, 'MLHO', 'all')
+md = setflowequation(md, 'MOLHO', 'all')
 md.cluster = generic('name', gethostname(), 'np', 3)
 
 #Transient
@@ -28 +28 @@
 md.frontalforcings.meltingrate = np.zeros((md.mesh.numberofvertices))
 md.levelset.migration_max = 1e10
 
-md = SetMLHOBC(md)
+md = SetMOLHOBC(md)
 md = solve(md, 'Transient')
 
 #Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test127.m

@@ -1 @@
-%Test Name: SquareShelfConstrainedStressMLHO2d
+%Test Name: SquareShelfConstrainedStressMOLHO2d
 md=triangle(model(),'../Exp/Square.exp',50000.);
 md=setmask(md,'all','');
 md=parameterize(md,'../Par/SquareShelfConstrained.par');
-md=setflowequation(md,'MLHO','all');
+md=setflowequation(md,'MOLHO','all');
 md.cluster=generic('name',oshostname(),'np',2);
 
 %output
-%FIXME compute the stress components for MLHO
+%FIXME compute the stress components for MOLHO
 md.stressbalance.requested_outputs={'default','VxSurface','VySurface','VxShear','VyShear','VxBase','VyBase','MassFlux1','MassFlux2','MassFlux3','MassFlux4','MassFlux5','MassFlux6'};
 %md.stressbalance.requested_outputs={'default','DeviatoricStressxx','DeviatoricStressyy','DeviatoricStressxy','MassFlux1','MassFlux2','MassFlux3','MassFlux4','MassFlux5','MassFlux6'};
 md.outputdefinition.definitions={...
@@ -17 +17 @@
         massfluxatgate('name','MassFlux5','profilename',['../Exp/MassFlux5.exp'],'definitionstring','Outputdefinition5'),...
         massfluxatgate('name','MassFlux6','profilename',['../Exp/MassFlux6.exp'],'definitionstring','Outputdefinition6')...
         };
-md=SetMLHOBC(md);
+md=SetMOLHOBC(md);
 md=solve(md,'Stressbalance');
 
 %Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test129.m

@@ -1 @@
-%Test Name: SquareShelfConstrainedRestartTranMLHO2d
+%Test Name: SquareShelfConstrainedRestartTranMOLHO2d
 md=triangle(model(),'../Exp/Square.exp',150000.);
 md=setmask(md,'all','');
 md=parameterize(md,'../Par/SquareShelfConstrained.par');
-md=setflowequation(md,'MLHO','all');
+md=setflowequation(md,'MOLHO','all');
 md.cluster=generic('name',oshostname(),'np',1);
 md.transient.requested_outputs={'IceVolume','TotalSmb','VxShear','VyShear','VxBase','VyBase','VxSurface','VySurface'};
 
@@ -13 +13 @@
 md.timestepping.final_time=19;
 md.settings.output_frequency=2;
 
-md=SetMLHOBC(md);
+md=SetMOLHOBC(md);
 md=solve(md,'Transient');
 md2=solve(md,'Transient','restart',1);
 

../trunk-jpl/test/NightlyRun/test249.m

@@ -1 @@
-%Test Name: SquareShelfTranMLHO2d
+%Test Name: SquareShelfTranMOLHO2d
 md=triangle(model(),'../Exp/Square.exp',150000.);
 md=setmask(md,'all','');
 md=parameterize(md,'../Par/SquareShelf.par');
-md=setflowequation(md,'MLHO','all');
+md=setflowequation(md,'MOLHO','all');
 md.cluster=generic('name',oshostname(),'np',3);
 md.transient.requested_outputs={'default','FloatingArea','GroundedArea','TotalGroundedBmb','TotalFloatingBmb'};
 md.basalforcings.floatingice_melting_rate(:)=1;
-md=SetMLHOBC(md);
+md=SetMOLHOBC(md);
 md=solve(md,'Transient');
 
 %Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test255.m

@@ -1 @@
-%Test Name: SquareShelfStressMLHO2dTransientIncrHydro
+%Test Name: SquareShelfStressMOLHO2dTransientIncrHydro
 md=triangle(model(),'../Exp/Square.exp',150000.);
 md=setmask(md,'all','');
 md=parameterize(md,'../Par/SquareShelf.par');
-md=setflowequation(md,'MLHO','all');
+md=setflowequation(md,'MOLHO','all');
 md.cluster=generic('name',oshostname(),'np',3);
 md.masstransport.hydrostatic_adjustment='Incremental';
-md=SetMLHOBC(md);
+md=SetMOLHOBC(md);
 md=solve(md,'Transient');
 
 %Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test330.m

@@ -1 @@
-%Test Name: SquareSheetConstrainedStressMLHO2d
+%Test Name: SquareSheetConstrainedStressMOLHO2d
 md=triangle(model(),'../Exp/Square.exp',150000.);
 md=setmask(md,'','');
 md=parameterize(md,'../Par/SquareSheetConstrained.par');
-md=setflowequation(md,'MLHO','all');
+md=setflowequation(md,'MOLHO','all');
 md.cluster=generic('name',oshostname(),'np',3);
 md.stressbalance.requested_outputs={'default','VxSurface','VySurface','VxShear','VyShear','VxBase','VyBase'};
-md=SetMLHOBC(md);
+md=SetMOLHOBC(md);
 md=solve(md,'Stressbalance');
 
 %Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test333.m

@@ -1 @@
-%Test Name: SquareSheetConstrainedCMDragMLHO
+%Test Name: SquareSheetConstrainedCMDragMOLHO
 md=triangle(model(),'../Exp/Square.exp',200000.);
 md=setmask(md,'','');
 md=parameterize(md,'../Par/SquareSheetConstrained.par');
-md=setflowequation(md,'MLHO','all');
+md=setflowequation(md,'MOLHO','all');
 
 %control parameters
 md.inversion.iscontrol=1;
@@ -18 +18 @@
 md.inversion.vx_obs=md.initialization.vx; md.inversion.vy_obs=md.initialization.vy;
 
 md.cluster=generic('name',oshostname(),'np',3);
-md=SetMLHOBC(md);
+md=SetMOLHOBC(md);
 md=solve(md,'Stressbalance');
 
 %Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test335.m

@@ -1 @@
-%Test Name: SquareSheetConstrainedStressMLHO3d
+%Test Name: SquareSheetConstrainedStressMOLHO3d
 md=triangle(model(),'../Exp/Square.exp',200000.);
 md=setmask(md,'','');
 md=parameterize(md,'../Par/SquareSheetConstrained.par');
-md=setflowequation(md,'MLHO','all');
+md=setflowequation(md,'MOLHO','all');
 md = extrude(md, 5, 1);
 md.cluster=generic('name',oshostname(),'np',3);
 md.stressbalance.requested_outputs={'default','VxSurface','VySurface','VxShear','VyShear','VxBase','VyBase'};
-md=SetMLHOBC(md);
+md=SetMOLHOBC(md);
 md=solve(md,'Stressbalance');
 
 %Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test447.m

@@ -1 @@
-%Test Name: SquareSheetShelfTranMLHO2d
+%Test Name: SquareSheetShelfTranMOLHO2d
 md=triangle(model(),'../Exp/Square.exp',150000.);
 md=setmask(md,'../Exp/SquareShelf.exp','');
 md=parameterize(md,'../Par/SquareSheetShelf.par');
@@ -12 +12 @@
 md.friction.coefficientcoulomb=0.02*ones(md.mesh.numberofvertices,1);
 md.transient.isthermal=0;
 md.transient.isgroundingline=1;
-md=setflowequation(md,'MLHO','all');
+md=setflowequation(md,'MOLHO','all');
 md.cluster=generic('name',oshostname(),'np',3);
 md.transient.requested_outputs={'default','GroundedArea','FloatingArea','TotalFloatingBmb','TotalGroundedBmb','TotalSmb'};
-md=SetMLHOBC(md);
+md=SetMOLHOBC(md);
 md=solve(md,'Transient');
 
 %Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test449.m

@@ -1 @@
-%Test Name: MISMIP3DMLHO
+%Test Name: MISMIP3DMOLHO
 md=triangle(model(),'../Exp/Square.exp',100000.);
 md=setmask(md,'../Exp/SquareShelf.exp','');
 md=parameterize(md,'../Par/SquareSheetShelf.par');
@@ -11 +11 @@
 pos=find(md.mask.ocean_levelset>=0);
 md.geometry.base(pos)=md.geometry.bed(pos);
 md.geometry.surface=md.geometry.base+md.geometry.thickness;
-md=setflowequation(md,'MLHO','all');
+md=setflowequation(md,'MOLHO','all');
 
 %Boundary conditions:
 md.mask.ice_levelset=-ones(md.mesh.numberofvertices,1);
@@ -44 +44 @@
 md.timestepping.time_step=10;
 
 md.cluster=generic('name',oshostname(),'np',3);
-md=SetMLHOBC(md);
+md=SetMOLHOBC(md);
 md=solve(md,'Transient');
 
 %Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test519.m

@@ -1 @@
-%Test Name: PigTranMLHO2d
+%Test Name: PigTranMOLHO2d
 md=triangle(model(),'../Exp/Pig.exp',20000.);
 md=setmask(md,'../Exp/PigShelves.exp','../Exp/PigIslands.exp');
 md=parameterize(md,'../Par/Pig.par');
-md=setflowequation(md,'MLHO','all');
+md=setflowequation(md,'MOLHO','all');
 md.cluster=generic('name',oshostname(),'np',3);
 md.geometry.bed=md.geometry.base;
 [md.mesh.lat,md.mesh.long] = xy2ll(md.mesh.x,md.mesh.y,-1);
 md.mesh.scale_factor=0.9*ones(md.mesh.numberofvertices,1);
 md.transient.requested_outputs={'default','IceVolume','IceVolumeScaled','GroundedArea','GroundedAreaScaled','FloatingArea','FloatingAreaScaled','TotalSmb','TotalSmbScaled','TotalFloatingBmb','TotalFloatingBmbScaled'};
-md=SetMLHOBC(md);
+md=SetMOLHOBC(md);
 md=solve(md,'Transient');
 
 %Fields and tolerances to track changes

../trunk-jpl/test/NightlyRun/test811.m

@@ -1 @@
-%Test Name: ValleyGlacierLevelsetCalvingMLHO2d
+%Test Name: ValleyGlacierLevelsetCalvingMOLHO2d
 md=triangle(model(),'../Exp/Square.exp',50000);
 md=setmask(md,'','');
 md=parameterize(md,'../Par/ValleyGlacierShelf.par');
-md=setflowequation(md,'MLHO','all');
+md=setflowequation(md,'MOLHO','all');
 md.cluster=generic('name',oshostname(),'np',3);
 
 %Transient
@@ -17 +17 @@
 md.frontalforcings.meltingrate=zeros(md.mesh.numberofvertices,1);
 md.levelset.migration_max = 1e10;
 
-md=SetMLHOBC(md);
+md=SetMOLHOBC(md);
 md=solve(md,'Transient');
 
 %Fields and tolerances to track changes