Index: /issm/trunk-jpl/src/c/classes/Elements/Element.h
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Element.h	(revision 20870)
+++ /issm/trunk-jpl/src/c/classes/Elements/Element.h	(revision 20871)
@@ -298,4 +298,5 @@
 		#ifdef _HAVE_SEALEVELRISE_
 		virtual IssmDouble    GetArea3D(void)=0;
+		virtual IssmDouble    GetAreaSpherical(void)=0;
 		virtual IssmDouble    OceanAverage(IssmDouble* Sg)=0;
 		virtual IssmDouble    OceanArea(void)=0;
Index: /issm/trunk-jpl/src/c/classes/Elements/Penta.h
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Penta.h	(revision 20870)
+++ /issm/trunk-jpl/src/c/classes/Elements/Penta.h	(revision 20871)
@@ -64,4 +64,5 @@
 		void           FSContactMigration(Vector<IssmDouble>* vertexgrounded,Vector<IssmDouble>* vertexfloating);
 		IssmDouble     GetArea3D(void){_error_("not implemented yet!");};
+		IssmDouble     GetAreaSpherical(void){_error_("not implemented yet!");};
 		void           GetAreaCoordinates(IssmDouble *area_coordinates,IssmDouble* xyz_zero,IssmDouble* xyz_list,int numpoints);
 		Element*       GetBasalElement(void);
Index: /issm/trunk-jpl/src/c/classes/Elements/Seg.h
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Seg.h	(revision 20870)
+++ /issm/trunk-jpl/src/c/classes/Elements/Seg.h	(revision 20871)
@@ -163,4 +163,5 @@
 		void        ZeroLevelsetCoordinates(IssmDouble** pxyz_zero,IssmDouble* xyz_list,int levelsetenum){_error_("not implemented");};
 		IssmDouble     GetArea3D(void){_error_("not implemented yet!");};
+		IssmDouble     GetAreaSpherical(void){_error_("not implemented yet!");};
 
 #ifdef _HAVE_GIA_
Index: /issm/trunk-jpl/src/c/classes/Elements/Tetra.h
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Tetra.h	(revision 20870)
+++ /issm/trunk-jpl/src/c/classes/Elements/Tetra.h	(revision 20871)
@@ -63,4 +63,5 @@
 		void        FSContactMigration(Vector<IssmDouble>* vertexgrounded,Vector<IssmDouble>* vertexfloating){_error_("not implemented yet");};
 		IssmDouble     GetArea3D(void){_error_("not implemented yet!");};
+		IssmDouble     GetAreaSpherical(void){_error_("not implemented yet!");};
 		Element*    GetBasalElement(void){_error_("not implemented yet");};
 		int         GetElementType(void);
Index: /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp	(revision 20870)
+++ /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp	(revision 20871)
@@ -912,4 +912,32 @@
 	return area_fraction*this->GetArea();
 }/*}}}*/
+IssmDouble Tria::GetAreaSpherical(void){/*{{{*/
+
+	bool spherical=true;
+	IssmDouble llr_list[NUMVERTICES][3];
+	IssmDouble x1,y1,z1,x2,y2,z2,x3,y3,z3;
+	IssmDouble arc12,arc23,arc31,semi_peri,excess; 
+
+	/*retrieve coordinates: lat,long,radius */
+	::GetVerticesCoordinates(&llr_list[0][0],vertices,NUMVERTICES,spherical);
+	x1=llr_list[0][0]/180*PI; y1=llr_list[0][1]/180*PI; z1=llr_list[0][2];
+	x2=llr_list[1][0]/180*PI; y2=llr_list[1][1]/180*PI; z2=llr_list[1][2];
+	x3=llr_list[2][0]/180*PI; y3=llr_list[2][1]/180*PI; z3=llr_list[2][2];
+
+	/*compute great circle distance between vertices */
+	arc12=2.*asin(sqrt(pow(sin((x2-x1)/2),2.0)+cos(x1)*cos(x2)*pow(sin((y2-y1)/2),2)));
+	arc23=2.*asin(sqrt(pow(sin((x3-x2)/2),2.0)+cos(x2)*cos(x3)*pow(sin((y3-y2)/2),2)));
+	arc31=2.*asin(sqrt(pow(sin((x1-x3)/2),2.0)+cos(x3)*cos(x1)*pow(sin((y1-y3)/2),2)));
+
+	/*semi parameter */ 
+	semi_peri=(arc12+arc23+arc31)/2; 
+
+	/*spherical excess */
+	excess=4.*atan(sqrt(tan(semi_peri/2)*tan((semi_peri-arc12)/2)*tan((semi_peri-arc23)/2)*tan((semi_peri-arc31)/2))); 
+
+	/*area = excess*radius^2 */
+	return excess*pow((z1+z2+z3)/3,2); 
+}
+/*}}}*/
 void       Tria::GetAreaCoordinates(IssmDouble* area_coordinates,IssmDouble* xyz_zero,IssmDouble* xyz_list,int numpoints){/*{{{*/
 	/*Computeportion of the element that is grounded*/ 
@@ -3543,5 +3571,5 @@
 IssmDouble    Tria::OceanArea(void){ /*{{{*/
 
-	if(IsWaterInElement()) return GetArea3D();
+	if(IsWaterInElement()) return GetAreaSpherical();
 	else return 0;
 
@@ -3555,5 +3583,5 @@
 
 		/*Compute area of element:*/
-		area=GetArea3D();
+		area=GetAreaSpherical();
 
 		/*Average Sg over vertices:*/
@@ -3658,5 +3686,5 @@
 
 	/*Compute area of element:*/
-	area=GetArea3D();
+	area=GetAreaSpherical();
 
 	/*Compute ice thickness change: */
@@ -3761,5 +3789,5 @@
 
 	/*Compute area of element:*/
-	area=GetArea3D();
+	area=GetAreaSpherical();
 
 	/* Where is the centroid of this element?:{{{*/
Index: /issm/trunk-jpl/src/c/classes/Elements/Tria.h
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Tria.h	(revision 20870)
+++ /issm/trunk-jpl/src/c/classes/Elements/Tria.h	(revision 20871)
@@ -157,4 +157,5 @@
 		IssmDouble 	   GetArea3D(void);
 		IssmDouble 	   GetAreaIce(void);
+		IssmDouble 	   GetAreaSpherical(void);
 		void           GetAreaCoordinates(IssmDouble *area_coordinates,IssmDouble* xyz_zero,IssmDouble* xyz_list,int numpoints);
 		int            GetElementType(void);
Index: /issm/trunk-jpl/src/c/classes/FemModel.cpp
===================================================================
--- /issm/trunk-jpl/src/c/classes/FemModel.cpp	(revision 20870)
+++ /issm/trunk-jpl/src/c/classes/FemModel.cpp	(revision 20871)
@@ -2322,5 +2322,5 @@
 		Element* element=xDynamicCast<Element*>(elements->GetObjectByOffset(i));
 		oceanarea_cpu += element->OceanArea();
-		eartharea_cpu+= element->GetArea3D();
+		eartharea_cpu+= element->GetAreaSpherical();
 	}
 	ISSM_MPI_Reduce (&oceanarea_cpu,&oceanarea,1,ISSM_MPI_DOUBLE,ISSM_MPI_SUM,0,IssmComm::GetComm() );
@@ -2382,5 +2382,5 @@
 		Element* element=xDynamicCast<Element*>(elements->GetObjectByOffset(i));
 		oceanarea_cpu += element->OceanArea();
-		eartharea_cpu+= element->GetArea3D();
+		eartharea_cpu+= element->GetAreaSpherical();
 	}
 	ISSM_MPI_Reduce (&oceanarea_cpu,&oceanarea,1,ISSM_MPI_DOUBLE,ISSM_MPI_SUM,0,IssmComm::GetComm() );