Index: /issm/trunk-jpl/src/c/analyses/GLheightadvectionAnalysis.cpp
===================================================================
--- /issm/trunk-jpl/src/c/analyses/GLheightadvectionAnalysis.cpp	(revision 23461)
+++ /issm/trunk-jpl/src/c/analyses/GLheightadvectionAnalysis.cpp	(revision 23462)
@@ -38,4 +38,5 @@
 		iomodel->FetchDataToInput(elements,"md.mesh.vertexonsurface",MeshVertexonsurfaceEnum);
 	}
+	iomodel->FetchDataToInput(elements,"md.mesh.vertexonboundary",MeshVertexonboundaryEnum);
 }/*}}}*/
 void GLheightadvectionAnalysis::UpdateParameters(Parameters* parameters,IoModel* iomodel,int solution_enum,int analysis_enum){/*{{{*/
@@ -61,5 +62,5 @@
 	/*Intermediaries */
 	int        domaintype,dim;
-	IssmDouble Jdet,D_scalar;
+	IssmDouble Jdet,D_scalar,onboundary;
 	IssmDouble vel,vx,vy;
 	IssmDouble* xyz_list = NULL;
@@ -89,4 +90,5 @@
    Input* vx_input=element->GetInput(VxEnum); _assert_(vx_input);
 	Input* vy_input=element->GetInput(VyEnum); _assert_(vy_input);
+	Input* bc_input=element->GetInput(MeshVertexonboundaryEnum); _assert_(bc_input);
 
 	IssmDouble h = element->CharacteristicLength();
@@ -104,7 +106,16 @@
 		D_scalar=gauss->weight*Jdet;
 
+		/*Get velocity*/
+		vx_input->GetInputValue(&vx,gauss);
+		vy_input->GetInputValue(&vy,gauss);
+		bc_input->GetInputValue(&onboundary,gauss);
+		if(onboundary>0.){
+			/*We do not want to advect garbage, make sure only diffusion is applied on boundary*/
+			vx = 0.; vy = 0.;
+		}
+
 		/*Diffusion */
 		if(sqrt(vx*vx+vy*vy)<1000./31536000.){
-			IssmPDouble kappa = -10;
+			IssmPDouble kappa = -10.;
 			for(int i=0;i<numnodes;i++){
 				for(int j=0;j<numnodes;j++){
@@ -113,7 +124,6 @@
 			}
 		}
+
 		/*Advection: */
-		vx_input->GetInputValue(&vx,gauss);
-		vy_input->GetInputValue(&vy,gauss);
 		for(int i=0;i<numnodes;i++){
 			for(int j=0;j<numnodes;j++){
@@ -122,5 +132,5 @@
 		}
 
-		/*Artifficial diffusivity*/
+		/*Artificial diffusivity*/
 		vel=sqrt(vx*vx + vy*vy)+1.e-14;
 		D[0][0]=D_scalar*h/(2.*vel)*fabs(vx*vx);  D[0][1]=D_scalar*h/(2.*vel)*fabs(vx*vy);