Index: /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp	(revision 15419)
+++ /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp	(revision 15420)
@@ -6822,8 +6822,26 @@
 ElementMatrix* Penta::CreateKMatrixDiagnosticStokes(void){
 
-	/*compute all stiffness matrices for this element*/
-	ElementMatrix* Ke1=CreateKMatrixDiagnosticStokesViscous();
-	ElementMatrix* Ke2=CreateKMatrixDiagnosticStokesFriction();
-	ElementMatrix* Ke =new ElementMatrix(Ke1,Ke2);
+	int fe_stokes;
+	ElementMatrix* Ke1;
+	ElementMatrix* Ke2;
+	ElementMatrix* Ke;
+	parameters->FindParam(&fe_stokes,FlowequationFeStokesEnum);
+
+	switch(fe_stokes){
+		case 0:
+			/*compute all stiffness matrices for this element*/
+			Ke1=CreateKMatrixDiagnosticStokesViscous();
+			Ke2=CreateKMatrixDiagnosticStokesFriction();
+			Ke =new ElementMatrix(Ke1,Ke2);
+
+		case 1:
+			/*compute all stiffness matrices for this element*/
+			Ke1=CreateKMatrixDiagnosticStokesViscous();
+			Ke2=CreateKMatrixDiagnosticStokesFriction();
+			Ke =new ElementMatrix(Ke1,Ke2);
+
+		default:
+			_error_("Finite element" << fe_stokes << " not supported yet");
+	}
 
 	/*clean-up and return*/
@@ -6831,4 +6849,5 @@
 	delete Ke2;
 	return Ke;
+
 }
 /*}}}*/
Index: /issm/trunk-jpl/src/c/modules/ModelProcessorx/CreateParameters.cpp
===================================================================
--- /issm/trunk-jpl/src/c/modules/ModelProcessorx/CreateParameters.cpp	(revision 15419)
+++ /issm/trunk-jpl/src/c/modules/ModelProcessorx/CreateParameters.cpp	(revision 15420)
@@ -49,4 +49,5 @@
 	parameters->AddObject(iomodel->CopyConstantObject(FlowequationIsl1l2Enum));
 	parameters->AddObject(iomodel->CopyConstantObject(FlowequationIsstokesEnum));
+	parameters->AddObject(iomodel->CopyConstantObject(FlowequationFeStokesEnum));
 	parameters->AddObject(iomodel->CopyConstantObject(SettingsOutputFrequencyEnum));
 	parameters->AddObject(iomodel->CopyConstantObject(DiagnosticRestolEnum));
Index: /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h
===================================================================
--- /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h	(revision 15419)
+++ /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h	(revision 15420)
@@ -71,4 +71,5 @@
 	FlowequationIsmacayealpattynEnum,
 	FlowequationIsstokesEnum,
+	FlowequationFeStokesEnum,
 	FlowequationVertexEquationEnum,
 	FrictionCoefficientEnum,
Index: /issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp
===================================================================
--- /issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp	(revision 15419)
+++ /issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp	(revision 15420)
@@ -79,4 +79,5 @@
 		case FlowequationIsmacayealpattynEnum : return "FlowequationIsmacayealpattyn";
 		case FlowequationIsstokesEnum : return "FlowequationIsstokes";
+		case FlowequationFeStokesEnum : return "FlowequationFeStokes";
 		case FlowequationVertexEquationEnum : return "FlowequationVertexEquation";
 		case FrictionCoefficientEnum : return "FrictionCoefficient";
Index: /issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp
===================================================================
--- /issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp	(revision 15419)
+++ /issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp	(revision 15420)
@@ -79,4 +79,5 @@
 	      else if (strcmp(name,"FlowequationIsmacayealpattyn")==0) return FlowequationIsmacayealpattynEnum;
 	      else if (strcmp(name,"FlowequationIsstokes")==0) return FlowequationIsstokesEnum;
+	      else if (strcmp(name,"FlowequationFeStokes")==0) return FlowequationFeStokesEnum;
 	      else if (strcmp(name,"FlowequationVertexEquation")==0) return FlowequationVertexEquationEnum;
 	      else if (strcmp(name,"FrictionCoefficient")==0) return FrictionCoefficientEnum;
@@ -136,9 +137,9 @@
 	      else if (strcmp(name,"InversionMaxParameters")==0) return InversionMaxParametersEnum;
 	      else if (strcmp(name,"InversionMaxiterPerStep")==0) return InversionMaxiterPerStepEnum;
-	      else if (strcmp(name,"InversionMinParameters")==0) return InversionMinParametersEnum;
          else stage=2;
    }
    if(stage==2){
-	      if (strcmp(name,"InversionNsteps")==0) return InversionNstepsEnum;
+	      if (strcmp(name,"InversionMinParameters")==0) return InversionMinParametersEnum;
+	      else if (strcmp(name,"InversionNsteps")==0) return InversionNstepsEnum;
 	      else if (strcmp(name,"InversionNumControlParameters")==0) return InversionNumControlParametersEnum;
 	      else if (strcmp(name,"InversionNumCostFunctions")==0) return InversionNumCostFunctionsEnum;
@@ -259,9 +260,9 @@
 	      else if (strcmp(name,"ThermalIsenthalpy")==0) return ThermalIsenthalpyEnum;
 	      else if (strcmp(name,"GiaMantleViscosity")==0) return GiaMantleViscosityEnum;
-	      else if (strcmp(name,"GiaLithosphereThickness")==0) return GiaLithosphereThicknessEnum;
          else stage=3;
    }
    if(stage==3){
-	      if (strcmp(name,"Thickness")==0) return ThicknessEnum;
+	      if (strcmp(name,"GiaLithosphereThickness")==0) return GiaLithosphereThicknessEnum;
+	      else if (strcmp(name,"Thickness")==0) return ThicknessEnum;
 	      else if (strcmp(name,"TimesteppingStartTime")==0) return TimesteppingStartTimeEnum;
 	      else if (strcmp(name,"TimesteppingFinalTime")==0) return TimesteppingFinalTimeEnum;
@@ -382,9 +383,9 @@
 	      else if (strcmp(name,"VectorParam")==0) return VectorParamEnum;
 	      else if (strcmp(name,"Riftfront")==0) return RiftfrontEnum;
-	      else if (strcmp(name,"RiftfrontType")==0) return RiftfrontTypeEnum;
          else stage=4;
    }
    if(stage==4){
-	      if (strcmp(name,"Segment")==0) return SegmentEnum;
+	      if (strcmp(name,"RiftfrontType")==0) return RiftfrontTypeEnum;
+	      else if (strcmp(name,"Segment")==0) return SegmentEnum;
 	      else if (strcmp(name,"SegmentRiftfront")==0) return SegmentRiftfrontEnum;
 	      else if (strcmp(name,"SpcDynamic")==0) return SpcDynamicEnum;
@@ -505,9 +506,9 @@
 	      else if (strcmp(name,"BoolElementResult")==0) return BoolElementResultEnum;
 	      else if (strcmp(name,"BoolExternalResult")==0) return BoolExternalResultEnum;
-	      else if (strcmp(name,"DoubleElementResult")==0) return DoubleElementResultEnum;
          else stage=5;
    }
    if(stage==5){
-	      if (strcmp(name,"DoubleExternalResult")==0) return DoubleExternalResultEnum;
+	      if (strcmp(name,"DoubleElementResult")==0) return DoubleElementResultEnum;
+	      else if (strcmp(name,"DoubleExternalResult")==0) return DoubleExternalResultEnum;
 	      else if (strcmp(name,"DoubleMatExternalResult")==0) return DoubleMatExternalResultEnum;
 	      else if (strcmp(name,"IntExternalResult")==0) return IntExternalResultEnum;
Index: /issm/trunk-jpl/src/m/classes/flowequation.m
===================================================================
--- /issm/trunk-jpl/src/m/classes/flowequation.m	(revision 15419)
+++ /issm/trunk-jpl/src/m/classes/flowequation.m	(revision 15420)
@@ -10,4 +10,5 @@
 		isl1l2               = 0;
 		isstokes             = 0;
+		fe_stokes            = 0;
 		vertex_equation      = NaN;
 		element_equation     = NaN;
@@ -78,4 +79,5 @@
 			WriteData(fid,'object',obj,'fieldname','isl1l2','format','Boolean');
 			WriteData(fid,'object',obj,'fieldname','isstokes','format','Boolean');
+			WriteData(fid,'object',obj,'fieldname','fe_stokes','format','Integer');
 			WriteData(fid,'object',obj,'fieldname','bordermacayeal','format','DoubleMat','mattype',1);
 			WriteData(fid,'object',obj,'fieldname','borderpattyn','format','DoubleMat','mattype',1);
Index: /issm/trunk-jpl/src/m/classes/flowequation.py
===================================================================
--- /issm/trunk-jpl/src/m/classes/flowequation.py	(revision 15419)
+++ /issm/trunk-jpl/src/m/classes/flowequation.py	(revision 15420)
@@ -20,4 +20,5 @@
 		self.isl1l2               = 0
 		self.isstokes             = 0
+		self.fe_stokes            = 0
 		self.vertex_equation      = float('NaN')
 		self.element_equation     = float('NaN')
@@ -78,4 +79,5 @@
 		WriteData(fid,'object',self,'fieldname','isl1l2','format','Boolean')
 		WriteData(fid,'object',self,'fieldname','isstokes','format','Boolean')
+		WriteData(fid,'object',self,'fieldname','fe_stokes','format','Integer')
 		WriteData(fid,'object',self,'fieldname','bordermacayeal','format','DoubleMat','mattype',1)
 		WriteData(fid,'object',self,'fieldname','borderpattyn','format','DoubleMat','mattype',1)
Index: /issm/trunk-jpl/src/m/enum/EnumDefinitions.py
===================================================================
--- /issm/trunk-jpl/src/m/enum/EnumDefinitions.py	(revision 15419)
+++ /issm/trunk-jpl/src/m/enum/EnumDefinitions.py	(revision 15420)
@@ -877,4 +877,18 @@
 	return StringToEnum('FlowequationIsstokes')[0]
 
+def FlowequationFeStokesEnum():
+	"""
+	FLOWEQUATIONFESTOKESENUM - Enum of FlowequationFeStokes
+
+	WARNING: DO NOT MODIFY THIS FILE
+				this file has been automatically generated by src/c/shared/Enum/Synchronize.sh
+				Please read src/c/shared/Enum/README for more information
+
+	   Usage:
+	      macro=FlowequationFeStokesEnum()
+	"""
+
+	return StringToEnum('FlowequationFeStokes')[0]
+
 def FlowequationVertexEquationEnum():
 	"""
@@ -7875,4 +7889,4 @@
 	"""
 
-	return 561
-
+	return 562
+
Index: /issm/trunk-jpl/src/m/enum/FlowequationFeStokesEnum.m
===================================================================
--- /issm/trunk-jpl/src/m/enum/FlowequationFeStokesEnum.m	(revision 15420)
+++ /issm/trunk-jpl/src/m/enum/FlowequationFeStokesEnum.m	(revision 15420)
@@ -0,0 +1,11 @@
+function macro=FlowequationFeStokesEnum()
+%FLOWEQUATIONFESTOKESENUM - Enum of FlowequationFeStokes
+%
+%   WARNING: DO NOT MODIFY THIS FILE
+%            this file has been automatically generated by src/c/shared/Enum/Synchronize.sh
+%            Please read src/c/shared/Enum/README for more information
+%
+%   Usage:
+%      macro=FlowequationFeStokesEnum()
+
+macro=StringToEnum('FlowequationFeStokes');
Index: /issm/trunk-jpl/src/m/enum/MaximumNumberOfEnums.m
===================================================================
--- /issm/trunk-jpl/src/m/enum/MaximumNumberOfEnums.m	(revision 15419)
+++ /issm/trunk-jpl/src/m/enum/MaximumNumberOfEnums.m	(revision 15420)
@@ -9,3 +9,3 @@
 %      macro=MaximumNumberOfEnums()
 
-macro=561;
+macro=562;