source: issm/oecreview/Archive/24684-25833/ISSM-25465-25466.diff

../trunk-jpl/src/c/cores/transient_core.cpp

@@ -72 +72 @@
         if(isoceancoupling) OceanExchangeDatax(femmodel,true);
         #endif
 
-        IssmDouble  output_value;
-        int         num_dependents;
-        IssmPDouble *dependents;
-        DataSet*    dependent_objects=NULL;
-        IssmDouble  J=0.;
-
+        DataSet* dependent_objects=NULL;
         if(iscontrol && isautodiff){
-
-                femmodel->parameters->SetParam(starttime,TimeEnum);
-                femmodel->parameters->FindParam(&num_dependents,AutodiffNumDependentsEnum);
+                femmodel->parameters->SetParam(starttime,TimeEnum); /*FIXME: DO WE REALLY NEED THIS?*/
                 femmodel->parameters->FindParam(&dependent_objects,AutodiffDependentObjectsEnum);
-
-                /*Go through our dependent variables, and compute the response:*/
-                dependents=xNew<IssmPDouble>(num_dependents);
         }
 
         if(isslr) sealevelrise_core_geometry(femmodel);
@@ -157 +147 @@
                         /*Go through our dependent variables, and compute the response:*/
                         for(int i=0;i<dependent_objects->Size();i++){
                                 DependentObject* dep=(DependentObject*)dependent_objects->GetObjectByOffset(i);
-                                 dep->Responsex(&output_value,femmodel);
-                                 dep->AddValue(output_value);
-                                }
+                                IssmDouble  output_value;
+                                dep->Responsex(&output_value,femmodel);
+                                dep->AddValue(output_value);
                         }
+                }
         }
 
         if(!iscontrol || !isautodiff) femmodel->RequestedDependentsx();
@@ -283 +274 @@
         }
 
 }/*}}}*/
+
+#ifdef _HAVE_CODIPACK_
+void transient_ad(FemModel* femmodel){/*{{{*/
+
+        /*parameters: */
+        IssmDouble finaltime,dt,yts,starttime;
+        bool       isoceancoupling,iscontrol,isautodiff,isgroundingline,isslr;
+        bool       save_results,dakota_analysis;
+        int        timestepping;
+        int        output_frequency;
+        int        sb_coupling_frequency;
+        int        recording_frequency;
+        int        groundingline_migration,smb_model,amr_frequency,amr_restart;
+        int        numoutputs;
+        char     **requested_outputs = NULL;
+
+        /*intermediary: */
+        int        step;
+        IssmDouble time;
+
+        /*first, figure out if there was a check point, if so, do a reset of the FemModel* femmodel structure. */
+        femmodel->parameters->FindParam(&recording_frequency,SettingsRecordingFrequencyEnum);
+        if(recording_frequency) femmodel->Restart();
+
+        /*then recover parameters common to all solutions*/
+        femmodel->parameters->FindParam(&step,StepEnum);
+        femmodel->parameters->FindParam(&time,TimeEnum);
+        femmodel->parameters->FindParam(&finaltime,TimesteppingFinalTimeEnum);
+        femmodel->parameters->FindParam(&starttime,TimesteppingStartTimeEnum);
+        femmodel->parameters->FindParam(&yts,ConstantsYtsEnum);
+        femmodel->parameters->FindParam(&dakota_analysis,QmuIsdakotaEnum);
+        femmodel->parameters->FindParam(&output_frequency,SettingsOutputFrequencyEnum);
+        femmodel->parameters->FindParam(&sb_coupling_frequency,SettingsSbCouplingFrequencyEnum);
+        femmodel->parameters->FindParam(&timestepping,TimesteppingTypeEnum);
+        femmodel->parameters->FindParam(&isslr,TransientIsslrEnum);
+        femmodel->parameters->FindParam(&isgroundingline,TransientIsgroundinglineEnum);
+        femmodel->parameters->FindParam(&amr_frequency,TransientAmrFrequencyEnum);
+        femmodel->parameters->FindParam(&iscontrol,InversionIscontrolEnum);
+        femmodel->parameters->FindParam(&isautodiff,AutodiffIsautodiffEnum);
+        femmodel->parameters->FindParam(&numoutputs,TransientNumRequestedOutputsEnum);
+        if(isgroundingline) femmodel->parameters->FindParam(&groundingline_migration,GroundinglineMigrationEnum);
+        if(numoutputs)      femmodel->parameters->FindParam(&requested_outputs,&numoutputs,TransientRequestedOutputsEnum);
+
+        DataSet*    dependent_objects=NULL;
+        if(iscontrol){
+                femmodel->parameters->SetParam(starttime,TimeEnum);
+                femmodel->parameters->FindParam(&dependent_objects,AutodiffDependentObjectsEnum);
+        }
+
+        if(isslr) sealevelrise_core_geometry(femmodel);
+
+        while(time < finaltime - (yts*DBL_EPSILON)){ //make sure we run up to finaltime.
+
+                /*Time Increment*/
+                switch(timestepping){
+                        case AdaptiveTimesteppingEnum:
+                                femmodel->TimeAdaptx(&dt);
+                                if(time+dt>finaltime) dt=finaltime-time;
+                                femmodel->parameters->SetParam(dt,TimesteppingTimeStepEnum);
+                                break;
+                        case FixedTimesteppingEnum:
+                                femmodel->parameters->FindParam(&dt,TimesteppingTimeStepEnum);
+                                break;
+                        default:
+                                _error_("Time stepping \""<<EnumToStringx(timestepping)<<"\" not supported yet");
+                }
+                step+=1;
+                time+=dt;
+                femmodel->parameters->SetParam(time,TimeEnum);
+                femmodel->parameters->SetParam(step,StepEnum);
+                femmodel->parameters->SetParam(false,SaveResultsEnum);
+
+                /*Run transient step!*/
+                transient_step(femmodel);
+
+                /*unload results*/
+                if(VerboseSolution()) _printf0_("   computing requested outputs\n");
+                femmodel->RequestedOutputsx(&femmodel->results,requested_outputs,numoutputs,false);
+
+                if(recording_frequency && (step%recording_frequency==0)){
+                        if(VerboseSolution()) _printf0_("   checkpointing model \n");
+                        femmodel->CheckPoint();
+                }
+
+                /*Go through our dependent variables, and compute the response:*/
+                for(int i=0;i<dependent_objects->Size();i++){
+                        DependentObject* dep=(DependentObject*)dependent_objects->GetObjectByOffset(i);
+                        IssmDouble  output_value;
+                        dep->Responsex(&output_value,femmodel);
+                        dep->AddValue(output_value);
+                }
+        }
+
+        if(!iscontrol) femmodel->RequestedDependentsx();
+        if(iscontrol) femmodel->parameters->SetParam(dependent_objects,AutodiffDependentObjectsEnum);
+
+        /*Free ressources:*/
+        if(numoutputs){for(int i=0;i<numoutputs;i++){xDelete<char>(requested_outputs[i]);} xDelete<char*>(requested_outputs);}
+
+}/*}}}*/
+#endif