Index: /issm/trunk-jpl/src/c/analyses/StressbalanceAnalysis.cpp
===================================================================
--- /issm/trunk-jpl/src/c/analyses/StressbalanceAnalysis.cpp	(revision 16873)
+++ /issm/trunk-jpl/src/c/analyses/StressbalanceAnalysis.cpp	(revision 16874)
@@ -975,16 +975,15 @@
 }/*}}}*/
 ElementMatrix* StressbalanceAnalysis::CreateKMatrixSSAFriction(Element* element){/*{{{*/
+
+	/*Return if element is inactive*/
+	if(element->IsFloating()) return NULL;
+
 	/*Intermediaries*/
 	bool       mainlyfloating;
-	int        analysis_type,migration_style,point1;
-	IssmDouble alpha2,Jdet,fraction1,fraction2;
-	IssmDouble phi=1.0;
-	IssmDouble D_scalar,gllevelset;
-	IssmDouble *xyz_list = NULL;
-	Friction   *friction = NULL;
-	Gauss*     gauss     = NULL;
-
-	/*Return if element is inactive*/
-	if(element->IsFloating()) return NULL;
+	int         migration_style,point1;
+	IssmDouble  alpha2,Jdet,fraction1,fraction2;
+	IssmDouble  gllevelset,phi=1.;
+	IssmDouble *xyz_list  = NULL;
+	Gauss*      gauss     = NULL;
 
 	/*Fetch number of nodes and dof for this finite element*/
@@ -993,19 +992,18 @@
 
 	/*Initialize Element matrix and vectors*/
-	ElementMatrix* Ke     = element->NewElementMatrix(SSAApproximationEnum);
-	IssmDouble*    B      = xNew<IssmDouble>(2*numdof);
-	IssmDouble*    D      = xNewZeroInit<IssmDouble>(2*2);
+	ElementMatrix* Ke      = element->NewElementMatrix(SSAApproximationEnum);
+	IssmDouble*    B       = xNew<IssmDouble>(2*numdof);
+	IssmDouble     D[2][2] = {0.};
 
 	/*Retrieve all inputs and parameters*/
 	element->GetVerticesCoordinates(&xyz_list);
 	element->FindParam(&migration_style,GroundinglineMigrationEnum);
-	Input* surface_input=element->GetInput(SurfaceEnum); _assert_(surface_input);
-	Input* vx_input=element->GetInput(VxEnum);           _assert_(vx_input);
-	Input* vy_input=element->GetInput(VyEnum);           _assert_(vy_input);
-	Input* gllevelset_input=NULL;
-	element->FindParam(&analysis_type,AnalysisTypeEnum);
+	Input* surface_input    = element->GetInput(SurfaceEnum); _assert_(surface_input);
+	Input* vx_input         = element->GetInput(VxEnum);      _assert_(vx_input);
+	Input* vy_input         = element->GetInput(VyEnum);      _assert_(vy_input);
+	Input* gllevelset_input = NULL;
 
 	/*build friction object, used later on: */
-	//friction=new Friction("2d",inputs,matpar,analysis_type);
+	Friction* friction=new Friction(element,2);
 
 	/*Recover portion of element that is grounded*/
@@ -1022,22 +1020,19 @@
 	/* Start  looping on the number of gaussian points: */
 	for(int ig=gauss->begin();ig<gauss->end();ig++){
-
 		gauss->GaussPoint(ig);
 
-		//friction->GetAlpha2(&alpha2, gauss,VxEnum,VyEnum,VzEnum);
-		alpha2=1.0;
+		friction->GetAlpha2(&alpha2, gauss,vx_input,vy_input,NULL);
 		if(migration_style==SubelementMigrationEnum) alpha2=phi*alpha2;
 		if(migration_style==SubelementMigration2Enum){
 			gllevelset_input->GetInputValue(&gllevelset, gauss);
-			if(gllevelset<0) alpha2=0;
+			if(gllevelset<0.) alpha2=0.;
 		}
 
 		this->GetBSSAFriction(B, element, xyz_list, gauss);
 		element->JacobianDeterminant(&Jdet, xyz_list,gauss);
-		D_scalar=alpha2*gauss->weight*Jdet;
-		for(int i=0;i<2;i++) D[i*2+i]=D_scalar;
+		for(int i=0;i<2;i++) D[i][i]=alpha2*gauss->weight*Jdet;
 
 		TripleMultiply(B,2,numdof,1,
-					D,2,2,0,
+					&D[0][0],2,2,0,
 					B,2,numdof,0,
 					&Ke->values[0],1);
@@ -1051,5 +1046,4 @@
 	delete friction;
 	xDelete<IssmDouble>(xyz_list);
-	xDelete<IssmDouble>(D);
 	xDelete<IssmDouble>(B);
 	return Ke;
@@ -1252,5 +1246,4 @@
 	delete gauss;
 	return pe;
-	return NULL;
 }/*}}}*/
 void StressbalanceAnalysis::GetBSSA(IssmDouble* B,Element* element,IssmDouble* xyz_list,Gauss* gauss){/*{{{*/