Changeset 27246 for issm/trunk-jpl/src/c/modules/FloatingiceMeltingRatex/FloatingiceMeltingRatex.cpp

Timestamp:

08/29/22 08:52:51 (3 years ago)

Author:

vverjans

Message:

CHG: clean-up of autoregression method

File:

• issm/trunk-jpl/src/c/modules/FloatingiceMeltingRatex/FloatingiceMeltingRatex.cpp (modified) (5 diffs)

issm/trunk-jpl/src/c/modules/FloatingiceMeltingRatex/FloatingiceMeltingRatex.cpp

@@ -210 +210 @@
 }
 /*}}}*/
-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){/*{{{*/
 
@@ -247 +219 @@
    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;
@@ -261 +231 @@
    bool isstochastic;
    bool isdeepmeltingstochastic = false;
-   int M,N,Nphi,arorder,numbasins,my_rank;
+   int M,N,Nlags,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* termconstant   = NULL;
+   IssmDouble* trend          = NULL;
+   IssmDouble* lagcoefs       = NULL;
    IssmDouble* deepwaterel    = NULL;
    IssmDouble* upperwaterel   = NULL;
@@ -275 +245 @@
         /*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);
+   femmodel->parameters->FindParam(&termconstant,&M,BasalforcingsautoregressionconstEnum);               _assert_(M==numbasins);
+   femmodel->parameters->FindParam(&trend,&M,BasalforcingsautoregressiontrendEnum);               _assert_(M==numbasins);
+   femmodel->parameters->FindParam(&lagcoefs,&M,&Nlags,BasalforcingsarlagcoefsEnum);             _assert_(M==numbasins); _assert_(Nlags==arorder);
 
         /*Get basin-specific parameters*/
@@ -303 +273 @@
       Element* element = xDynamicCast<Element*>(object);
       /*Compute autoregression*/
-      element->Autoregression(isstepforar,arorder,telapsed_ar,beta0,beta1,phi,isdeepmeltingstochastic,BasalforcingsDeepwaterMeltingRateAutoregressionEnum);
+      element->Autoregression(isstepforar,arorder,telapsed_ar,tstep_ar,termconstant,trend,lagcoefs,isdeepmeltingstochastic,BasalforcingsDeepwaterMeltingRateAutoregressionEnum);
                 element->BasinLinearFloatingiceMeltingRate(deepwaterel,upperwatermelt,upperwaterel,perturbation);
         }
 
         /*Cleanup*/
-        xDelete<IssmDouble>(beta0);
-        xDelete<IssmDouble>(beta1);
-        xDelete<IssmDouble>(phi);
+        xDelete<IssmDouble>(termconstant);
+        xDelete<IssmDouble>(trend);
+        xDelete<IssmDouble>(lagcoefs);
         xDelete<IssmDouble>(deepwaterel);
         xDelete<IssmDouble>(upperwaterel);