Index: issm/trunk-jpl/src/c/analyses/MasstransportAnalysis.cpp
===================================================================
--- issm/trunk-jpl/src/c/analyses/MasstransportAnalysis.cpp	(revision 26835)
+++ issm/trunk-jpl/src/c/analyses/MasstransportAnalysis.cpp	(revision 26836)
@@ -219,4 +219,7 @@
 			iomodel->FetchDataToInput(inputs,elements,"md.basalforcings.ocean_temp",BasalforcingsOceanTempEnum);
 			break;
+		case AutoregressionLinearFloatingMeltRateEnum:
+			iomodel->FetchDataToInput(inputs,elements,"md.basalforcings.basin_id",BasalforcingsLinearBasinIdEnum);
+			break;
 		default:
 			_error_("Basal forcing model "<<EnumToStringx(basalforcing_model)<<" not supported yet");
Index: issm/trunk-jpl/src/c/classes/Elements/Element.cpp
===================================================================
--- issm/trunk-jpl/src/c/classes/Elements/Element.cpp	(revision 26835)
+++ issm/trunk-jpl/src/c/classes/Elements/Element.cpp	(revision 26836)
@@ -69,6 +69,7 @@
 
    /*Get Basin ID and Basin coefficients*/
-   if(enum_type==SMBautoregressionEnum)                   this->GetInputValue(&basinid,SmbBasinsIdEnum);
-   if(enum_type==FrontalForcingsRignotAutoregressionEnum) this->GetInputValue(&basinid,FrontalForcingsBasinIdEnum);
+   if(enum_type==SMBautoregressionEnum)                               this->GetInputValue(&basinid,SmbBasinsIdEnum);
+   if(enum_type==FrontalForcingsRignotAutoregressionEnum)             this->GetInputValue(&basinid,FrontalForcingsBasinIdEnum);
+   if(enum_type==BasalforcingsDeepwaterMeltingRateAutoregressionEnum) this->GetInputValue(&basinid,BasalforcingsLinearBasinIdEnum);
    for(int i=0;i<arorder;i++) phi_basin[i] = phi[basinid*arorder+i];
    beta0_basin   = beta0[basinid];
@@ -92,6 +93,7 @@
       xDelete<IssmDouble>(oldvarspin);
    }
-   if(enum_type==SMBautoregressionEnum)                   this->inputs->SetArrayInput(SmbValuesAutoregressionEnum,this->lid,varspin,numvertices*arorder);
-   if(enum_type==FrontalForcingsRignotAutoregressionEnum) this->inputs->SetArrayInput(ThermalforcingValuesAutoregressionEnum,this->lid,varspin,numvertices*arorder);
+   if(enum_type==SMBautoregressionEnum)                               this->inputs->SetArrayInput(SmbValuesAutoregressionEnum,this->lid,varspin,numvertices*arorder);
+   if(enum_type==FrontalForcingsRignotAutoregressionEnum)             this->inputs->SetArrayInput(ThermalforcingValuesAutoregressionEnum,this->lid,varspin,numvertices*arorder);
+   if(enum_type==BasalforcingsDeepwaterMeltingRateAutoregressionEnum) this->inputs->SetArrayInput(BasalforcingsDeepwaterMeltingRateValuesAutoregressionEnum,this->lid,varspin,numvertices*arorder);
 
    /*Cleanup and return*/
@@ -123,4 +125,10 @@
          outenum_type   = FrontalForcingsThermalForcingEnum;
          break;
+		case(BasalforcingsDeepwaterMeltingRateAutoregressionEnum):
+         arenum_type    = BasalforcingsDeepwaterMeltingRateValuesAutoregressionEnum;
+         basinenum_type = BasalforcingsLinearBasinIdEnum;
+         noiseenum_type = BasalforcingsDeepwaterMeltingRateNoiseEnum;
+         outenum_type   = BasalforcingsSpatialDeepwaterMeltingRateEnum;
+         break;
    }
 
@@ -172,4 +180,35 @@
    xDelete<IssmDouble>(varlist);
    xDelete<IssmDouble>(valuesautoregression);
+}/*}}}*/
+void       Element::BasinLinearFloatingiceMeltingRate(IssmDouble* deepwaterel,IssmDouble* upperwatermelt,IssmDouble* upperwaterel,IssmDouble* perturbation){/*{{{*/
+
+	const int NUM_VERTICES = this->GetNumberOfVertices();
+
+	int basinid;
+   IssmDouble alpha;
+   IssmDouble base[MAXVERTICES];
+   IssmDouble values[MAXVERTICES];
+   IssmDouble deepwatermelt[MAXVERTICES];
+
+	/*Get element-specific values*/
+	this->GetInputValue(&basinid,BasalforcingsLinearBasinIdEnum);
+	this->GetInputListOnVertices(&base[0],BaseEnum);
+   this->GetInputListOnVertices(&deepwatermelt[0],BasalforcingsSpatialDeepwaterMeltingRateEnum);
+
+	/*Compute melt rate at vertices according to basin-specific values of input arguments*/
+   for(int i=0;i<NUM_VERTICES;i++){
+		if(base[i]>=upperwaterel[basinid]){
+         values[i]=upperwatermelt[basinid];
+      }
+      else if (base[i]<deepwaterel[basinid]){
+         values[i]=deepwatermelt[i];
+      }
+      else{
+         alpha = (base[i]-upperwaterel[basinid])/(deepwaterel[basinid]-upperwaterel[basinid]);
+         values[i]=deepwatermelt[i]*alpha+(1.-alpha)*upperwatermelt[basinid];
+      }
+   }
+
+   this->AddInput(BasalforcingsFloatingiceMeltingRateEnum,&values[0],P1Enum);
 }/*}}}*/
 void       Element::ComputeLambdaS(){/*{{{*/
@@ -2301,5 +2340,4 @@
 	IssmDouble base[MAXVERTICES];
 	IssmDouble values[MAXVERTICES];
-	IssmDouble perturbation[MAXVERTICES];
 	IssmDouble time;
 
@@ -2312,5 +2350,4 @@
 
 	this->GetInputListOnVertices(&base[0],BaseEnum);
-   this->GetInputListOnVertices(&perturbation[0],BasalforcingsPerturbationMeltingRateEnum);
 	for(int i=0;i<NUM_VERTICES;i++){
 		if(base[i]>=upperwaterel){
@@ -2324,6 +2361,4 @@
 			values[i]=deepwatermelt*alpha+(1.-alpha)*upperwatermelt;
 		}
-
-      values[i]+=perturbation[i];
 	}
 
Index: issm/trunk-jpl/src/c/classes/Elements/Element.h
===================================================================
--- issm/trunk-jpl/src/c/classes/Elements/Element.h	(revision 26835)
+++ issm/trunk-jpl/src/c/classes/Elements/Element.h	(revision 26836)
@@ -70,4 +70,5 @@
 		void					 AutoregressionInit(int numbasins,int arorder,int nspin,IssmDouble starttime,IssmDouble tstep_ar,IssmDouble tinit_ar,IssmDouble* beta0,IssmDouble* beta1,IssmDouble* phi,int enum_type);
       void               Autoregression(bool isstepforar,int arorder,IssmDouble telapsed_ar,IssmDouble* beta0,IssmDouble* beta1,IssmDouble* phi,bool isfieldstochastic,int enum_type);
+		void               BasinLinearFloatingiceMeltingRate(IssmDouble* deepwaterel,IssmDouble* upperwatermelt,IssmDouble* upperwaterel,IssmDouble* perturbation);
 		void               ComputeLambdaS(void);
 		void               ComputeNewDamage();
Index: issm/trunk-jpl/src/c/modules/FloatingiceMeltingRatex/FloatingiceMeltingRatex.cpp
===================================================================
--- issm/trunk-jpl/src/c/modules/FloatingiceMeltingRatex/FloatingiceMeltingRatex.cpp	(revision 26835)
+++ issm/trunk-jpl/src/c/modules/FloatingiceMeltingRatex/FloatingiceMeltingRatex.cpp	(revision 26836)
@@ -45,4 +45,8 @@
 			if(VerboseSolution())_printf0_("        call BeckmannGoosse Floating melting rate module\n");
 			BeckmannGoosseFloatingiceMeltingRatex(femmodel);
+			break;
+		case AutoregressionLinearFloatingMeltRateEnum:
+			if(VerboseSolution())_printf0_("        call Autoregression Linear Floating melting rate module\n");
+			AutoregressionLinearFloatingiceMeltingRatex(femmodel);
 			break;
 		default:
@@ -206,2 +210,99 @@
 }
 /*}}}*/
+void AutoregressionLinearFloatingiceMeltingRateInitx(FemModel* femmodel){/*{{{*/
+
+	/*Initialization step of AutoregressionLinearFloatingiceMeltingRatex*/
+   int M,N,Nphi,arorder,numbasins,my_rank;
+   IssmDouble starttime,tstep_ar,tinit_ar;
+   femmodel->parameters->FindParam(&numbasins,BasalforcingsLinearNumBasinsEnum);
+   femmodel->parameters->FindParam(&arorder,BasalforcingsAutoregressiveOrderEnum);
+   IssmDouble* beta0    = NULL;
+   IssmDouble* beta1    = NULL;
+   IssmDouble* phi      = NULL;
+   femmodel->parameters->FindParam(&starttime,TimesteppingStartTimeEnum);
+   femmodel->parameters->FindParam(&tstep_ar,BasalforcingsAutoregressionTimestepEnum);
+   femmodel->parameters->FindParam(&tinit_ar,BasalforcingsAutoregressionInitialTimeEnum);
+   femmodel->parameters->FindParam(&beta0,&M,BasalforcingsBeta0Enum);    _assert_(M==numbasins);
+   femmodel->parameters->FindParam(&beta1,&M,BasalforcingsBeta1Enum);    _assert_(M==numbasins);
+   femmodel->parameters->FindParam(&phi,&M,&Nphi,BasalforcingsPhiEnum);  _assert_(M==numbasins); _assert_(Nphi==arorder);
+	
+	/*AR model spin-up with 0 noise to initialize BasalforcingsDeepwaterMeltingRateValuesAutoregressionEnum (688 = log(0.001)/log(0.99): decaying time of inluence of phi[0]=0.99 to 0.001 of beta_0*/
+   int nspin = 688;
+   for(Object* &object:femmodel->elements->objects){
+      Element* element = xDynamicCast<Element*>(object); //generate element object
+      element->AutoregressionInit(numbasins,arorder,nspin,starttime,tstep_ar,tinit_ar,beta0,beta1,phi,BasalforcingsDeepwaterMeltingRateAutoregressionEnum);
+	}
+	/*Cleanup*/
+   xDelete<IssmDouble>(beta0);
+   xDelete<IssmDouble>(beta1);
+   xDelete<IssmDouble>(phi);
+}/*}}}*/
+void AutoregressionLinearFloatingiceMeltingRatex(FemModel* femmodel){/*{{{*/
+
+	/*Get time parameters*/
+   IssmDouble time,dt,starttime,tstep_ar;
+   femmodel->parameters->FindParam(&time,TimeEnum);
+   femmodel->parameters->FindParam(&dt,TimesteppingTimeStepEnum);
+   femmodel->parameters->FindParam(&starttime,TimesteppingStartTimeEnum);
+   femmodel->parameters->FindParam(&tstep_ar,BasalforcingsAutoregressionTimestepEnum);
+
+   /*Initialize module at first time step*/
+   if(time<=starttime+dt){AutoregressionLinearFloatingiceMeltingRateInitx(femmodel);}
+   /*Determine if this is a time step for the AR model*/
+   bool isstepforar = false;
+
+   #ifndef _HAVE_AD_
+   if((fmod(time,tstep_ar)<fmod((time-dt),tstep_ar)) || (time<=starttime+dt) || tstep_ar==dt) isstepforar = true;
+   #else
+   _error_("not implemented yet");
+   #endif
+
+   /*Load parameters*/
+   bool isstochastic;
+   bool isdeepmeltingstochastic = false;
+   int M,N,Nphi,arorder,numbasins,my_rank;
+   femmodel->parameters->FindParam(&numbasins,BasalforcingsLinearNumBasinsEnum);
+	femmodel->parameters->FindParam(&arorder,BasalforcingsAutoregressiveOrderEnum);
+   IssmDouble tinit_ar;
+   IssmDouble* beta0          = NULL;
+   IssmDouble* beta1          = NULL;
+   IssmDouble* phi            = NULL;
+   IssmDouble* deepwaterel    = NULL;
+   IssmDouble* upperwaterel   = NULL;
+   IssmDouble* upperwatermelt = NULL;
+   IssmDouble* perturbation   = NULL;
+
+	/*Get autoregressive parameters*/
+   femmodel->parameters->FindParam(&tinit_ar,BasalforcingsAutoregressionInitialTimeEnum);
+   femmodel->parameters->FindParam(&beta0,&M,BasalforcingsBeta0Enum);               _assert_(M==numbasins);
+   femmodel->parameters->FindParam(&beta1,&M,BasalforcingsBeta1Enum);               _assert_(M==numbasins);
+   femmodel->parameters->FindParam(&phi,&M,&Nphi,BasalforcingsPhiEnum);             _assert_(M==numbasins); _assert_(Nphi==arorder);
+
+	/*Get basin-specific parameters*/
+   femmodel->parameters->FindParam(&deepwaterel,&M,BasalforcingsDeepwaterElevationEnum);            _assert_(M==numbasins);
+   femmodel->parameters->FindParam(&upperwaterel,&M,BasalforcingsUpperwaterElevationEnum);          _assert_(M==numbasins);
+   femmodel->parameters->FindParam(&upperwatermelt,&M,BasalforcingsUpperwaterMeltingRateEnum);      _assert_(M==numbasins);
+
+	/*Evaluate whether stochasticity on DeepwaterMeltingRate is requested*/
+	femmodel->parameters->FindParam(&isstochastic,StochasticForcingIsStochasticForcingEnum);
+   if(isstochastic){
+      int  numstochasticfields;
+      int* stochasticfields;
+      femmodel->parameters->FindParam(&numstochasticfields,StochasticForcingNumFieldsEnum);
+      femmodel->parameters->FindParam(&stochasticfields,&N,StochasticForcingFieldsEnum); _assert_(N==numstochasticfields);
+      for(int i=0;i<numstochasticfields;i++){
+         if(stochasticfields[i]==BasalforcingsDeepwaterMeltingRateAutoregressionEnum) isdeepmeltingstochastic = true;
+      }
+      xDelete<int>(stochasticfields);
+   }
+   /*Time elapsed with respect to AR model initial time*/
+   IssmDouble telapsed_ar = time-tinit_ar;
+
+	/*Loop over each element to compute FloatingiceMeltingRate at vertices*/
+   for(Object* &object:femmodel->elements->objects){
+      Element* element = xDynamicCast<Element*>(object);
+      /*Compute autoregression*/
+      element->Autoregression(isstepforar,arorder,telapsed_ar,beta0,beta1,phi,isdeepmeltingstochastic,BasalforcingsDeepwaterMeltingRateAutoregressionEnum);
+		element->BasinLinearFloatingiceMeltingRate(deepwaterel,upperwatermelt,upperwaterel,perturbation);
+	}
+}/*}}}*/
Index: issm/trunk-jpl/src/c/modules/FloatingiceMeltingRatex/FloatingiceMeltingRatex.h
===================================================================
--- issm/trunk-jpl/src/c/modules/FloatingiceMeltingRatex/FloatingiceMeltingRatex.h	(revision 26835)
+++ issm/trunk-jpl/src/c/modules/FloatingiceMeltingRatex/FloatingiceMeltingRatex.h	(revision 26836)
@@ -16,4 +16,6 @@
 void FloatingiceMeltingRateIsmip6x(FemModel* femmodel);
 void BeckmannGoosseFloatingiceMeltingRatex(FemModel* femmodel);
+void AutoregressionLinearFloatingiceMeltingRateInitx(FemModel* femmodel);
+void AutoregressionLinearFloatingiceMeltingRatex(FemModel* femmodel);
 
 #endif  /* _FloatingiceMeltingRatex_H*/
Index: issm/trunk-jpl/src/c/modules/ModelProcessorx/CreateParameters.cpp
===================================================================
--- issm/trunk-jpl/src/c/modules/ModelProcessorx/CreateParameters.cpp	(revision 26835)
+++ issm/trunk-jpl/src/c/modules/ModelProcessorx/CreateParameters.cpp	(revision 26836)
@@ -252,4 +252,29 @@
 			break;
 		case BeckmannGoosseFloatingMeltRateEnum:
+			break;
+		case AutoregressionLinearFloatingMeltRateEnum:
+			/*Add parameters that are not in standard nbvertices format*/
+         parameters->AddObject(iomodel->CopyConstantObject("md.basalforcings.num_basins",BasalforcingsLinearNumBasinsEnum));
+         parameters->AddObject(iomodel->CopyConstantObject("md.basalforcings.ar_order",BasalforcingsAutoregressiveOrderEnum));
+         parameters->AddObject(iomodel->CopyConstantObject("md.basalforcings.ar_initialtime",BasalforcingsAutoregressionInitialTimeEnum));
+         parameters->AddObject(iomodel->CopyConstantObject("md.basalforcings.ar_timestep",BasalforcingsAutoregressionTimestepEnum));
+			iomodel->FetchData(&transparam,&M,&N,"md.basalforcings.beta0");
+         parameters->AddObject(new DoubleVecParam(BasalforcingsBeta0Enum,transparam,N));
+         xDelete<IssmDouble>(transparam);
+         iomodel->FetchData(&transparam,&M,&N,"md.basalforcings.beta1");
+         parameters->AddObject(new DoubleVecParam(BasalforcingsBeta1Enum,transparam,N));
+         xDelete<IssmDouble>(transparam);
+         iomodel->FetchData(&transparam,&M,&N,"md.basalforcings.phi");
+         parameters->AddObject(new DoubleMatParam(BasalforcingsPhiEnum,transparam,M,N));
+         xDelete<IssmDouble>(transparam);
+			iomodel->FetchData(&transparam,&M,&N,"md.basalforcings.upperwater_melting_rate");
+         parameters->AddObject(new DoubleVecParam(BasalforcingsUpperwaterMeltingRateEnum,transparam,N));
+         xDelete<IssmDouble>(transparam);   
+         iomodel->FetchData(&transparam,&M,&N,"md.basalforcings.upperwater_elevation");
+         parameters->AddObject(new DoubleVecParam(BasalforcingsUpperwaterElevationEnum,transparam,N));
+         xDelete<IssmDouble>(transparam); 
+			iomodel->FetchData(&transparam,&M,&N,"md.basalforcings.deepwater_elevation");
+         parameters->AddObject(new DoubleVecParam(BasalforcingsDeepwaterElevationEnum,transparam,N));
+         xDelete<IssmDouble>(transparam);   
 			break;
 		default:
Index: issm/trunk-jpl/src/c/shared/Enum/Enum.vim
===================================================================
--- issm/trunk-jpl/src/c/shared/Enum/Enum.vim	(revision 26835)
+++ issm/trunk-jpl/src/c/shared/Enum/Enum.vim	(revision 26836)
@@ -64,8 +64,14 @@
 syn keyword cConstant BalancethicknessStabilizationEnum
 syn keyword cConstant BarystaticContributionsEnum
+syn keyword cConstant BasalforcingsAutoregressionInitialTimeEnum
+syn keyword cConstant BasalforcingsAutoregressionTimestepEnum
+syn keyword cConstant BasalforcingsAutoregressiveOrderEnum
+syn keyword cConstant BasalforcingsBeta0Enum
+syn keyword cConstant BasalforcingsBeta1Enum
 syn keyword cConstant BasalforcingsBottomplumedepthEnum
 syn keyword cConstant BasalforcingsCrustthicknessEnum
 syn keyword cConstant BasalforcingsDeepwaterElevationEnum
 syn keyword cConstant BasalforcingsDeepwaterMeltingRateEnum
+syn keyword cConstant BasalforcingsDeepwaterMeltingRateNoiseEnum
 syn keyword cConstant BasalforcingsDtbgEnum
 syn keyword cConstant BasalforcingsEnum
@@ -77,7 +83,9 @@
 syn keyword cConstant BasalforcingsIsmip6NumBasinsEnum
 syn keyword cConstant BasalforcingsIsmip6TfDepthsEnum
+syn keyword cConstant BasalforcingsLinearNumBasinsEnum
 syn keyword cConstant BasalforcingsLowercrustheatEnum
 syn keyword cConstant BasalforcingsMantleconductivityEnum
 syn keyword cConstant BasalforcingsNusseltEnum
+syn keyword cConstant BasalforcingsPhiEnum
 syn keyword cConstant BasalforcingsPicoAverageOverturningEnum
 syn keyword cConstant BasalforcingsPicoAverageSalinityEnum
@@ -594,7 +602,10 @@
 syn keyword cConstant BalancethicknessThickeningRateEnum
 syn keyword cConstant BasalCrevasseEnum
+syn keyword cConstant BasalforcingsDeepwaterMeltingRateAutoregressionEnum
+syn keyword cConstant BasalforcingsDeepwaterMeltingRateValuesAutoregressionEnum
 syn keyword cConstant BasalforcingsFloatingiceMeltingRateEnum
 syn keyword cConstant BasalforcingsGeothermalfluxEnum
 syn keyword cConstant BasalforcingsGroundediceMeltingRateEnum
+syn keyword cConstant BasalforcingsLinearBasinIdEnum
 syn keyword cConstant BasalforcingsPerturbationMeltingRateEnum
 syn keyword cConstant BasalforcingsSpatialDeepwaterElevationEnum
@@ -1212,4 +1223,5 @@
 syn keyword cConstant ArrheniusEnum
 syn keyword cConstant AutodiffJacobianEnum
+syn keyword cConstant AutoregressionLinearFloatingMeltRateEnum
 syn keyword cConstant Balancethickness2AnalysisEnum
 syn keyword cConstant Balancethickness2SolutionEnum
Index: issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h
===================================================================
--- issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h	(revision 26835)
+++ issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h	(revision 26836)
@@ -58,8 +58,14 @@
 	BalancethicknessStabilizationEnum,
 	BarystaticContributionsEnum,
+	BasalforcingsAutoregressionInitialTimeEnum,
+	BasalforcingsAutoregressionTimestepEnum,
+	BasalforcingsAutoregressiveOrderEnum,
+	BasalforcingsBeta0Enum,
+	BasalforcingsBeta1Enum,
 	BasalforcingsBottomplumedepthEnum,
 	BasalforcingsCrustthicknessEnum,
 	BasalforcingsDeepwaterElevationEnum,
 	BasalforcingsDeepwaterMeltingRateEnum,
+	BasalforcingsDeepwaterMeltingRateNoiseEnum,
 	BasalforcingsDtbgEnum,
 	BasalforcingsEnum,
@@ -71,7 +77,9 @@
 	BasalforcingsIsmip6NumBasinsEnum,
 	BasalforcingsIsmip6TfDepthsEnum,
+	BasalforcingsLinearNumBasinsEnum,
 	BasalforcingsLowercrustheatEnum,
 	BasalforcingsMantleconductivityEnum,
 	BasalforcingsNusseltEnum,
+	BasalforcingsPhiEnum,
 	BasalforcingsPicoAverageOverturningEnum,
 	BasalforcingsPicoAverageSalinityEnum,
@@ -590,7 +598,10 @@
 	BalancethicknessThickeningRateEnum,
 	BasalCrevasseEnum,
+	BasalforcingsDeepwaterMeltingRateAutoregressionEnum,
+	BasalforcingsDeepwaterMeltingRateValuesAutoregressionEnum,
 	BasalforcingsFloatingiceMeltingRateEnum,
 	BasalforcingsGeothermalfluxEnum,
 	BasalforcingsGroundediceMeltingRateEnum,
+	BasalforcingsLinearBasinIdEnum,
 	BasalforcingsPerturbationMeltingRateEnum,
 	BasalforcingsSpatialDeepwaterElevationEnum,
@@ -1211,4 +1222,5 @@
 	ArrheniusEnum,
 	AutodiffJacobianEnum,
+	AutoregressionLinearFloatingMeltRateEnum,
 	Balancethickness2AnalysisEnum,
 	Balancethickness2SolutionEnum,
Index: issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp
===================================================================
--- issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp	(revision 26835)
+++ issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp	(revision 26836)
@@ -66,8 +66,14 @@
 		case BalancethicknessStabilizationEnum : return "BalancethicknessStabilization";
 		case BarystaticContributionsEnum : return "BarystaticContributions";
+		case BasalforcingsAutoregressionInitialTimeEnum : return "BasalforcingsAutoregressionInitialTime";
+		case BasalforcingsAutoregressionTimestepEnum : return "BasalforcingsAutoregressionTimestep";
+		case BasalforcingsAutoregressiveOrderEnum : return "BasalforcingsAutoregressiveOrder";
+		case BasalforcingsBeta0Enum : return "BasalforcingsBeta0";
+		case BasalforcingsBeta1Enum : return "BasalforcingsBeta1";
 		case BasalforcingsBottomplumedepthEnum : return "BasalforcingsBottomplumedepth";
 		case BasalforcingsCrustthicknessEnum : return "BasalforcingsCrustthickness";
 		case BasalforcingsDeepwaterElevationEnum : return "BasalforcingsDeepwaterElevation";
 		case BasalforcingsDeepwaterMeltingRateEnum : return "BasalforcingsDeepwaterMeltingRate";
+		case BasalforcingsDeepwaterMeltingRateNoiseEnum : return "BasalforcingsDeepwaterMeltingRateNoise";
 		case BasalforcingsDtbgEnum : return "BasalforcingsDtbg";
 		case BasalforcingsEnum : return "Basalforcings";
@@ -79,7 +85,9 @@
 		case BasalforcingsIsmip6NumBasinsEnum : return "BasalforcingsIsmip6NumBasins";
 		case BasalforcingsIsmip6TfDepthsEnum : return "BasalforcingsIsmip6TfDepths";
+		case BasalforcingsLinearNumBasinsEnum : return "BasalforcingsLinearNumBasins";
 		case BasalforcingsLowercrustheatEnum : return "BasalforcingsLowercrustheat";
 		case BasalforcingsMantleconductivityEnum : return "BasalforcingsMantleconductivity";
 		case BasalforcingsNusseltEnum : return "BasalforcingsNusselt";
+		case BasalforcingsPhiEnum : return "BasalforcingsPhi";
 		case BasalforcingsPicoAverageOverturningEnum : return "BasalforcingsPicoAverageOverturning";
 		case BasalforcingsPicoAverageSalinityEnum : return "BasalforcingsPicoAverageSalinity";
@@ -596,7 +604,10 @@
 		case BalancethicknessThickeningRateEnum : return "BalancethicknessThickeningRate";
 		case BasalCrevasseEnum : return "BasalCrevasse";
+		case BasalforcingsDeepwaterMeltingRateAutoregressionEnum : return "BasalforcingsDeepwaterMeltingRateAutoregression";
+		case BasalforcingsDeepwaterMeltingRateValuesAutoregressionEnum : return "BasalforcingsDeepwaterMeltingRateValuesAutoregression";
 		case BasalforcingsFloatingiceMeltingRateEnum : return "BasalforcingsFloatingiceMeltingRate";
 		case BasalforcingsGeothermalfluxEnum : return "BasalforcingsGeothermalflux";
 		case BasalforcingsGroundediceMeltingRateEnum : return "BasalforcingsGroundediceMeltingRate";
+		case BasalforcingsLinearBasinIdEnum : return "BasalforcingsLinearBasinId";
 		case BasalforcingsPerturbationMeltingRateEnum : return "BasalforcingsPerturbationMeltingRate";
 		case BasalforcingsSpatialDeepwaterElevationEnum : return "BasalforcingsSpatialDeepwaterElevation";
@@ -1214,4 +1225,5 @@
 		case ArrheniusEnum : return "Arrhenius";
 		case AutodiffJacobianEnum : return "AutodiffJacobian";
+		case AutoregressionLinearFloatingMeltRateEnum : return "AutoregressionLinearFloatingMeltRate";
 		case Balancethickness2AnalysisEnum : return "Balancethickness2Analysis";
 		case Balancethickness2SolutionEnum : return "Balancethickness2Solution";
Index: issm/trunk-jpl/src/c/shared/Enum/Enumjl.vim
===================================================================
--- issm/trunk-jpl/src/c/shared/Enum/Enumjl.vim	(revision 26835)
+++ issm/trunk-jpl/src/c/shared/Enum/Enumjl.vim	(revision 26836)
@@ -57,8 +57,14 @@
 syn keyword juliaConstC BalancethicknessStabilizationEnum
 syn keyword juliaConstC BarystaticContributionsEnum
+syn keyword juliaConstC BasalforcingsAutoregressionInitialTimeEnum
+syn keyword juliaConstC BasalforcingsAutoregressionTimestepEnum
+syn keyword juliaConstC BasalforcingsAutoregressiveOrderEnum
+syn keyword juliaConstC BasalforcingsBeta0Enum
+syn keyword juliaConstC BasalforcingsBeta1Enum
 syn keyword juliaConstC BasalforcingsBottomplumedepthEnum
 syn keyword juliaConstC BasalforcingsCrustthicknessEnum
 syn keyword juliaConstC BasalforcingsDeepwaterElevationEnum
 syn keyword juliaConstC BasalforcingsDeepwaterMeltingRateEnum
+syn keyword juliaConstC BasalforcingsDeepwaterMeltingRateNoiseEnum
 syn keyword juliaConstC BasalforcingsDtbgEnum
 syn keyword juliaConstC BasalforcingsEnum
@@ -70,7 +76,9 @@
 syn keyword juliaConstC BasalforcingsIsmip6NumBasinsEnum
 syn keyword juliaConstC BasalforcingsIsmip6TfDepthsEnum
+syn keyword juliaConstC BasalforcingsLinearNumBasinsEnum
 syn keyword juliaConstC BasalforcingsLowercrustheatEnum
 syn keyword juliaConstC BasalforcingsMantleconductivityEnum
 syn keyword juliaConstC BasalforcingsNusseltEnum
+syn keyword juliaConstC BasalforcingsPhiEnum
 syn keyword juliaConstC BasalforcingsPicoAverageOverturningEnum
 syn keyword juliaConstC BasalforcingsPicoAverageSalinityEnum
@@ -587,7 +595,10 @@
 syn keyword juliaConstC BalancethicknessThickeningRateEnum
 syn keyword juliaConstC BasalCrevasseEnum
+syn keyword juliaConstC BasalforcingsDeepwaterMeltingRateAutoregressionEnum
+syn keyword juliaConstC BasalforcingsDeepwaterMeltingRateValuesAutoregressionEnum
 syn keyword juliaConstC BasalforcingsFloatingiceMeltingRateEnum
 syn keyword juliaConstC BasalforcingsGeothermalfluxEnum
 syn keyword juliaConstC BasalforcingsGroundediceMeltingRateEnum
+syn keyword juliaConstC BasalforcingsLinearBasinIdEnum
 syn keyword juliaConstC BasalforcingsPerturbationMeltingRateEnum
 syn keyword juliaConstC BasalforcingsSpatialDeepwaterElevationEnum
@@ -1205,4 +1216,5 @@
 syn keyword juliaConstC ArrheniusEnum
 syn keyword juliaConstC AutodiffJacobianEnum
+syn keyword juliaConstC AutoregressionLinearFloatingMeltRateEnum
 syn keyword juliaConstC Balancethickness2AnalysisEnum
 syn keyword juliaConstC Balancethickness2SolutionEnum
Index: issm/trunk-jpl/src/c/shared/io/Marshalling/IoCodeConversions.cpp
===================================================================
--- issm/trunk-jpl/src/c/shared/io/Marshalling/IoCodeConversions.cpp	(revision 26835)
+++ issm/trunk-jpl/src/c/shared/io/Marshalling/IoCodeConversions.cpp	(revision 26836)
@@ -255,4 +255,5 @@
 		case 7: return BasalforcingsIsmip6Enum;
 		case 8: return BeckmannGoosseFloatingMeltRateEnum;
+		case 9: return AutoregressionLinearFloatingMeltRateEnum;
 		default: _error_("Marshalled Basal Forcings code \""<<enum_in<<"\" not supported yet");
 	}
Index: issm/trunk-jpl/src/m/classes/autoregressionlinearbasalforcings.m
===================================================================
--- issm/trunk-jpl/src/m/classes/autoregressionlinearbasalforcings.m	(revision 26836)
+++ issm/trunk-jpl/src/m/classes/autoregressionlinearbasalforcings.m	(revision 26836)
@@ -0,0 +1,117 @@
+%AUTOREGRESSION LINEAR BASAL FORCINGS class definition
+%
+%   Usage:
+%      autoregressionlinearbasalforcings=spatiallinearbasalforcings();
+
+classdef autoregressionlinearbasalforcings
+	
+	%VV: not verified yet (27Jan2022)
+	
+	properties (SetAccess=public) 
+		num_basins                = 0;
+		beta0                     = NaN;
+      beta1                     = NaN;
+      ar_order                  = 0;
+      ar_initialtime            = 0;
+      ar_timestep               = 0;
+      phi                       = NaN;
+      basin_id                  = NaN;
+		groundedice_melting_rate  = NaN;
+		deepwater_elevation       = NaN;
+		upperwater_melting_rate   = NaN;
+		upperwater_elevation      = NaN;
+		geothermalflux            = NaN;
+	end
+	methods
+		function self = autoregressionlinearbasalforcings(varargin) % {{{
+			switch nargin
+				case 0
+					self=setdefaultparameters(self);
+				otherwise
+					error('constructor not supported');
+			end
+		end % }}}
+		function self = extrude(self,md) % {{{
+			self.groundedice_melting_rate=project3d(md,'vector',self.groundedice_melting_rate,'type','node','layer',1); 
+			self.deepwater_elevation=project3d(md,'vector',self.deepwater_elevation,'type','node','layer',1); 
+			self.upperwater_melting_rate=project3d(md,'vector',self.upperwater_melting_rate,'type','node','layer',1); 
+			self.upperwater_elevation=project3d(md,'vector',self.upperwater_elevation,'type','node','layer',1); 
+			self.geothermalflux=project3d(md,'vector',self.geothermalflux,'type','node','layer',1); %bedrock only gets geothermal flux
+		end % }}}
+		function self = initialize(self,md) % {{{
+
+			if isnan(self.groundedice_melting_rate),
+				self.groundedice_melting_rate=zeros(md.mesh.numberofvertices,1);
+				disp('      no basalforcings.groundedice_melting_rate specified: values set as zero');
+			end
+
+		end % }}}
+		function self = setdefaultparameters(self) % {{{
+			%Nothing for now
+		end % }}}
+		function md = checkconsistency(self,md,solution,analyses) % {{{
+
+			if ismember('MasstransportAnalysis',analyses) & ~(strcmp(solution,'TransientSolution') & md.transient.ismasstransport==0),
+				md = checkfield(md,'fieldname','basalforcings.num_basins','numel',1,'NaN',1,'Inf',1,'>',0);
+				md = checkfield(md,'fieldname','basalforcings.groundedice_melting_rate','NaN',1,'Inf',1,'timeseries',1);
+				md = checkfield(md,'fieldname','basalforcings.deepwater_elevation','NaN',1,'Inf',1,'size',[1,md.basalforcings.num_basins],'numel',md.basalforcings.num_basins);
+				md = checkfield(md,'fieldname','basalforcings.upperwater_melting_rate','NaN',1,'Inf',1,'>=',0,'size',[1,md.basalforcings.num_basins],'numel',md.basalforcings.num_basins);
+				md = checkfield(md,'fieldname','basalforcings.upperwater_elevation','NaN',1,'Inf',1,'<=',0,'size',[1,md.basalforcings.num_basins],'numel',md.basalforcings.num_basins);
+            md = checkfield(md,'fieldname','basalforcings.basin_id','Inf',1,'>=',0,'<=',md.basalforcings.num_basins,'size',[md.mesh.numberofelements,1]);
+            md = checkfield(md,'fieldname','basalforcings.beta0','NaN',1,'Inf',1,'size',[1,md.basalforcings.num_basins],'numel',md.basalforcings.num_basins); 
+				md = checkfield(md,'fieldname','basalforcings.beta1','NaN',1,'Inf',1,'size',[1,md.basalforcings.num_basins],'numel',md.basalforcings.num_basins); 
+				md = checkfield(md,'fieldname','basalforcings.ar_order','numel',1,'NaN',1,'Inf',1,'>=',0);
+            md = checkfield(md,'fieldname','basalforcings.ar_initialtime','numel',1,'NaN',1,'Inf',1);
+            md = checkfield(md,'fieldname','basalforcings.ar_timestep','numel',1,'NaN',1,'Inf',1,'>=',md.timestepping.time_step); %autoregression time step cannot be finer than ISSM timestep
+            md = checkfield(md,'fieldname','basalforcings.phi','NaN',1,'Inf',1,'size',[md.basalforcings.num_basins,md.basalforcings.ar_order]);
+			end
+			if ismember('BalancethicknessAnalysis',analyses),
+				error('not implemented yet!');
+			end
+			if ismember('ThermalAnalysis',analyses) & ~(strcmp(solution,'TransientSolution') & md.transient.isthermal==0),
+				error('not implemented yet!');
+			end
+			if numel(md.basalforcings.geothermalflux)>1
+            md = checkfield(md,'fieldname','basalforcings.geothermalflux','NaN',1,'Inf',1,'timeseries',1,'>=',0);
+         end
+		end % }}}
+		function disp(self) % {{{
+			disp(sprintf('   autoregression linear basal forcings parameters:'));
+			disp(sprintf('   autoregressive model is applied for deepwater_melting_rate'));
+
+			fielddisplay(self,'num_basins','number of different basins [unitless]');
+         fielddisplay(self,'basin_id','basin number assigned to each element [unitless]');
+         fielddisplay(self,'beta0','basin-specific intercept values [m/yr] (if beta_1==0 mean=beta_0/(1-sum(phi)))');
+         fielddisplay(self,'beta1','basin-specific trend values [m  yr^(-2)]');
+         fielddisplay(self,'ar_order','order of the autoregressive model [unitless]');
+         fielddisplay(self,'ar_initialtime','initial time assumed in the autoregressive model parameterization [yr]');
+         fielddisplay(self,'ar_timestep','time resolution of the autoregressive model [yr]');
+         fielddisplay(self,'phi','basin-specific vectors of lag coefficients [unitless]');
+			fielddisplay(self,'deepwater_elevation','basin-specific elevation of ocean deepwater [m]');
+			fielddisplay(self,'upperwater_melting_rate','basin-specific basal melting rate (positive if melting applied for floating ice whith base >= upperwater_elevation) [m/yr]');
+			fielddisplay(self,'upperwater_elevation','basin-specific elevation of ocean upperwater [m]');
+			fielddisplay(self,'groundedice_melting_rate','node-specific basal melting rate (positive if melting) [m/yr]');
+			fielddisplay(self,'geothermalflux','node-specific geothermal heat flux [W/m^2]');
+
+		end % }}}
+		function marshall(self,prefix,md,fid) % {{{
+
+			yts=md.constants.yts;
+
+			WriteData(fid,prefix,'name','md.basalforcings.model','data',9,'format','Integer');
+			WriteData(fid,prefix,'object',self,'fieldname','groundedice_melting_rate','format','DoubleMat','name','md.basalforcings.groundedice_melting_rate','mattype',1,'scale',1./yts,'timeserieslength',md.mesh.numberofvertices+1,'yts',md.constants.yts);
+			WriteData(fid,prefix,'object',self,'fieldname','geothermalflux','name','md.basalforcings.geothermalflux','format','DoubleMat','mattype',1,'timeserieslength',md.mesh.numberofvertices+1,'yts',md.constants.yts);
+			WriteData(fid,prefix,'object',self,'fieldname','num_basins','format','Integer');
+         WriteData(fid,prefix,'object',self,'fieldname','ar_order','format','Integer');
+         WriteData(fid,prefix,'object',self,'fieldname','ar_initialtime','format','Double','scale',yts);
+         WriteData(fid,prefix,'object',self,'fieldname','ar_timestep','format','Double','scale',yts);
+         WriteData(fid,prefix,'object',self,'fieldname','basin_id','data',self.basin_id-1,'name','md.basalforcings.basin_id','format','IntMat','mattype',2); %0-indexed
+         WriteData(fid,prefix,'object',self,'fieldname','beta0','format','DoubleMat','name','md.basalforcings.beta0','scale',1./yts,'yts',yts);
+         WriteData(fid,prefix,'object',self,'fieldname','beta1','format','DoubleMat','name','md.basalforcings.beta1','scale',1./(yts^2),'yts',yts);
+         WriteData(fid,prefix,'object',self,'fieldname','phi','format','DoubleMat','name','md.basalforcings.phi','yts',yts);	
+			WriteData(fid,prefix,'object',self,'fieldname','deepwater_elevation','format','DoubleMat','name','md.basalforcings.deepwater_elevation');
+			WriteData(fid,prefix,'object',self,'fieldname','upperwater_melting_rate','format','DoubleMat','name','md.basalforcings.upperwater_melting_rate','scale',1./yts);
+			WriteData(fid,prefix,'object',self,'fieldname','upperwater_elevation','format','DoubleMat','name','md.basalforcings.upperwater_elevation');
+		end % }}}
+	end
+end
Index: issm/trunk-jpl/src/m/classes/autoregressionlinearbasalforcings.py
===================================================================
--- issm/trunk-jpl/src/m/classes/autoregressionlinearbasalforcings.py	(revision 26836)
+++ issm/trunk-jpl/src/m/classes/autoregressionlinearbasalforcings.py	(revision 26836)
@@ -0,0 +1,112 @@
+import numpy as np
+
+from checkfield import *
+from fielddisplay import fielddisplay
+from project3d import *
+from WriteData import *
+
+class autoregressionlinearbasalforcings(object):
+    """autoregressionlinearbasalforcings class definition
+
+    Usage:
+        autoregressionlinearbasalforcings = autoregressionlinearbasalforcings()
+    """
+
+    def __init__(self, *args):  # {{{
+        self.num_basins = 0
+        self.beta0 = np.nan
+        self.beta1 = np.nan
+        self.ar_order = 0
+        self.ar_initialtime = 0
+        self.ar_timestep = 0
+        self.phi = np.nan
+        self.basin_id = np.nan
+        self.groundedice_melting_rate = np.nan
+        self.deepwater_elevation = np.nan
+        self.upperwater_melting_rate = np.nan
+        self.upperwater_elevation = np.nan
+        self.geothermalflux = np.nan
+
+        nargs = len(args)
+        if nargs == 0:
+            self.setdefaultparameters()
+        else:
+            raise Exception('constructor not supported')
+    # }}}
+
+    def __repr__(self):  # {{{
+        s = '   surface forcings parameters:\n'
+        s += '   autoregressive model is applied for deepwater_melting_rate\n'
+        s += '{}\n'.format(fielddisplay(self, 'num_basins', 'number of different basins [unitless]'))
+        s += '{}\n'.format(fielddisplay(self, 'basin_id', 'basin number assigned to each element [unitless]'))
+        s += '{}\n'.format(fielddisplay(self, 'beta0', 'basin-specific intercept values [m ice eq./yr] (if beta_1==0 mean=beta_0/(1-sum(phi)))'))
+        s += '{}\n'.format(fielddisplay(self, 'beta1', 'basin-specific trend values [m ice eq. yr^(-2)]'))
+        s += '{}\n'.format(fielddisplay(self, 'ar_order', 'order of the autoregressive model [unitless]'))
+        s += '{}\n'.format(fielddisplay(self, 'ar_initialtime', 'initial time assumed in the autoregressive model parameterization [yr]'))
+        s += '{}\n'.format(fielddisplay(self, 'ar_timestep', 'time resolution of the autoregressive model [yr]'))
+        s += '{}\n'.format(fielddisplay(self, 'phi', 'basin-specific vectors of lag coefficients [unitless]'))
+        s += '{}\n'.format(fielddisplay(self, 'deepwater_elevation', 'basin-specific elevation of ocean deepwater [m]'))
+        s += '{}\n'.format(fielddisplay(self, 'upperwater_melting_rate', 'basin-specic basal melting rate (positive if melting applied for floating ice whith base >= upperwater_elevation) [m/yr]'))
+        s += '{}\n'.format(fielddisplay(self, 'upperwater_elevation', 'basin-specific elevation of ocean upperwater [m]'))
+        s += '{}\n'.format(fielddisplay(self, 'groundedice_melting_rate','node-specific basal melting rate (positive if melting) [m/yr]'))
+        s += '{}\n'.format(fielddisplay(self, 'geothermalflux','node-specific geothermal heat flux [W/m^2]'))
+        return s
+    # }}}
+
+    def setdefaultparameters(self): #{{{
+        #Nothing for now
+        return self
+    # }}}
+
+    def extrude(self, md):  # {{{
+        return self # Nothing for now
+    # }}}
+
+    def initialize(self, md):  # {{{
+        if np.all(np.isnan(self.groundedice_melting_rate)):
+            self.groundedice_melting_rate = np.zeros((md.mesh.numberofvertices))
+            print("      no basalforcings.groundedice_melting_rate specified: values set as zero")
+        return self
+    # }}}
+
+    def checkconsistency(self, md, solution, analyses):  # {{{
+        if 'MasstransportAnalysis' in analyses:
+            md = checkfield(md, 'fieldname', 'basalforcings.num_basins', 'numel', 1, 'NaN', 1, 'Inf', 1, '>', 0)
+            md = checkfield(md,'fieldname','basalforcings.groundedice_melting_rate','NaN',1,'Inf',1,'timeseries',1)
+            md = checkfield(md, 'fieldname', 'basalforcings.deepwater_elevation', 'NaN', 1, 'Inf', 1, 'size', [1, md.basalforcings.num_basins], 'numel', md.basalforcings.num_basins) 
+            md = checkfield(md, 'fieldname', 'basalforcings.upperwater_elevation', 'NaN', 1, 'Inf', 1,'<=',0, 'size', [1, md.basalforcings.num_basins], 'numel', md.basalforcings.num_basins) 
+            md = checkfield(md, 'fieldname', 'basalforcings.upperwater_melting_rate', 'NaN', 1, 'Inf', 1,'>=',0, 'size', [1, md.basalforcings.num_basins], 'numel', md.basalforcings.num_basins) 
+            md = checkfield(md, 'fieldname', 'basalforcings.basin_id', 'Inf', 1, '>=', 0, '<=', md.basalforcings.num_basins, 'size', [md.mesh.numberofelements])
+            md = checkfield(md, 'fieldname', 'basalforcings.beta0', 'NaN', 1, 'Inf', 1, 'size', [1, md.basalforcings.num_basins], 'numel', md.basalforcings.num_basins) # Scheme fails if passed as column vector
+            md = checkfield(md, 'fieldname', 'basalforcings.beta1', 'NaN', 1, 'Inf', 1, 'size', [1, md.basalforcings.num_basins], 'numel', md.basalforcings.num_basins) # Scheme fails if passed as column vector; NOTE: As opposed to MATLAB implementation, pass list
+            md = checkfield(md, 'fieldname', 'basalforcings.ar_order', 'numel', 1, 'NaN', 1, 'Inf', 1, '>=', 0)
+            md = checkfield(md, 'fieldname', 'basalforcings.ar_initialtime', 'numel', 1, 'NaN', 1, 'Inf', 1)
+            md = checkfield(md, 'fieldname', 'basalforcings.ar_timestep', 'numel', 1, 'NaN', 1, 'Inf', 1, '>=', md.timestepping.time_step) # Autoregression time step cannot be finer than ISSM timestep
+            md = checkfield(md, 'fieldname', 'basalforcings.phi', 'NaN', 1, 'Inf', 1, 'size', [md.basalforcings.num_basins, md.basalforcings.ar_order])
+        if 'BalancethicknessAnalysis' in analyses:
+            raise Exception('not implemented yet!')
+        if 'ThermalAnalysis' in analyses and not solution == 'TransientSolution' and not md.transient.isthermal:
+            raise Exception('not implemented yet!')
+
+        return md
+    # }}}
+
+    def marshall(self, prefix, md, fid):  # {{{
+        yts = md.constants.yts
+
+        WriteData(fid, prefix, 'name', 'md.basalforcings.model', 'data', 9, 'format', 'Integer')
+        WriteData(fid, prefix, 'object', self, 'fieldname', 'groundedice_melting_rate', 'name', 'md.basalforcings.groundedice_melting_rate', 'format', 'DoubleMat', 'mattype', 1, 'scale', 1. / yts, 'timeserieslength', md.mesh.numberofvertices + 1, 'yts', yts)
+        WriteData(fid, prefix, 'object', self, 'fieldname', 'geothermalflux', 'name', 'md.basalforcings.geothermalflux', 'format', 'DoubleMat', 'mattype', 1, 'timeserieslength', md.mesh.numberofvertices + 1, 'yts', yts)
+        WriteData(fid, prefix, 'object', self, 'fieldname', 'num_basins', 'format', 'Integer')
+        WriteData(fid, prefix, 'object', self, 'fieldname', 'ar_order', 'format', 'Integer')
+        WriteData(fid, prefix, 'object', self, 'fieldname', 'ar_initialtime', 'format', 'Double', 'scale', yts)
+        WriteData(fid, prefix, 'object', self, 'fieldname', 'ar_timestep', 'format', 'Double', 'scale', yts)
+        WriteData(fid, prefix, 'object', self, 'fieldname', 'basin_id', 'data', self.basin_id, 'name', 'md.basalforcings.basin_id', 'format', 'IntMat', 'mattype', 2) # 0-indexed
+        WriteData(fid, prefix, 'object', self, 'fieldname', 'beta0', 'format', 'DoubleMat', 'name', 'md.basalforcings.beta0', 'scale', 1 / yts, 'yts', yts)
+        WriteData(fid, prefix, 'object', self, 'fieldname', 'beta1', 'format', 'DoubleMat', 'name', 'md.basalforcings.beta1', 'scale', 1 / (yts ** 2), 'yts', yts)
+        WriteData(fid, prefix, 'object', self, 'fieldname', 'phi', 'format', 'DoubleMat', 'name', 'md.basalforcings.phi', 'yts', yts)
+        WriteData(fid, prefix, 'object', self, 'fieldname', 'deepwater_elevation', 'format', 'DoubleMat', 'name', 'md.basalforcings.deepwater_elevation')
+        WriteData(fid, prefix, 'object', self, 'fieldname', 'upperwater_melting_rate', 'format', 'DoubleMat', 'name', 'md.basalforcings.upperwater_melting_rate', 'scale', 1 / yts, 'yts', yts)
+        WriteData(fid, prefix, 'object', self, 'fieldname', 'upperwater_elevation', 'format', 'DoubleMat', 'name', 'md.basalforcings.upperwater_elevation')
+
+    # }}}
