Changeset 12870

Timestamp:

08/02/12 09:33:23 (13 years ago)

Author:

Mathieu Morlighem

Message:

BUG: Delta18O was not using the right array for temp and precip

Location:

issm/trunk-jpl/src/c

Files:

• classes/objects/Elements/Tria.cpp (modified) (4 diffs)
• solutions/transient_core.cpp (modified) (1 diff)

issm/trunk-jpl/src/c/classes/objects/Elements/Tria.cpp

@@ -910 +910 @@
 /*FUNCTION Tria::Delta18oParameterization{{{*/
 void  Tria::Delta18oParameterization(void){
+
+        int        i;
         IssmDouble monthlytemperatures[NUMVERTICES][12],monthlyprec[NUMVERTICES][12];
         IssmDouble TemperaturesPresentday[NUMVERTICES][12],TemperaturesLgm[NUMVERTICES][12];
         IssmDouble PrecipitationsPresentday[NUMVERTICES][12];
+        IssmDouble temp[NUMVERTICES];
         IssmDouble Delta18oPresent,Delta18oLgm,Delta18oTime;
         IssmDouble Delta18oSurfacePresent,Delta18oSurfaceLgm,Delta18oSurfaceTime;
         IssmDouble time,yts,finaltime;
+
+        /*Recover parameters*/
         this->parameters->FindParam(&time,TimeEnum);
         this->parameters->FindParam(&yts,ConstantsYtsEnum);
@@ -921 +926 @@
 
         /*Recover present day temperature and precipitation*/
-        Input*     input=inputs->GetInput(SurfaceforcingsTemperaturesPresentdayEnum); _assert_(input);
-        Input*     input2=inputs->GetInput(SurfaceforcingsTemperaturesLgmEnum); _assert_(input2);
-        Input*     input3=inputs->GetInput(SurfaceforcingsPrecipitationsPresentdayEnum); _assert_(input3);
+        Input* input=inputs->GetInput(SurfaceforcingsTemperaturesPresentdayEnum);    _assert_(input);
+        Input* input2=inputs->GetInput(SurfaceforcingsTemperaturesLgmEnum);          _assert_(input2);
+        Input* input3=inputs->GetInput(SurfaceforcingsPrecipitationsPresentdayEnum); _assert_(input3);
         GaussTria* gauss=new GaussTria();
-        for(int month=0;month<12;month++) {
-                for(int iv=0;iv<NUMVERTICES;iv++) {
+        for(int month=0;month<12;month++){
+                for(int iv=0;iv<NUMVERTICES;iv++){
                         gauss->GaussVertex(iv);
                         input->GetInputValue(&TemperaturesPresentday[iv][month],gauss,month/12.*yts);
@@ -934 +939 @@
                 }
         }
+
         /*Recover delta18o and Delta18oSurface at present day, lgm and at time t*/
         this->parameters->FindParam(&Delta18oPresent,SurfaceforcingsDelta18oEnum,finaltime);
-        this->parameters->FindParam(&Delta18oLgm,SurfaceforcingsDelta18oEnum,(finaltime-(21000*yts)));
+        this->parameters->FindParam(&Delta18oLgm,SurfaceforcingsDelta18oEnum,(finaltime-(21000*yts)));
         this->parameters->FindParam(&Delta18oTime,SurfaceforcingsDelta18oEnum,time);
         this->parameters->FindParam(&Delta18oSurfacePresent,SurfaceforcingsDelta18oSurfaceEnum,finaltime);
@@ -950 +956 @@
                                         &monthlytemperatures[iv][0], &monthlyprec[iv][0]);
         }
-        //if(id==1)     printf(" monthlytemperatures1 %f\n",monthlytemperatures[iv][0]);
-        //if(id==1)     printf(" monthlyprec1 %f\n",monthlyprec[iv][0]*yts);
+
         /*Update inputs*/ 
         TransientInput* NewTemperatureInput = new TransientInput(SurfaceforcingsMonthlytemperaturesEnum);
         TransientInput* NewPrecipitationInput = new TransientInput(SurfaceforcingsPrecipitationEnum);
         for (int imonth=0;imonth<12;imonth++) {
-                        TriaP1Input* newmonthinput1 = new TriaP1Input(SurfaceforcingsMonthlytemperaturesEnum,&TemperaturesPresentday[0][imonth]);
-                        //TriaP1Input* newmonthinput1 = new TriaP1Input(SurfaceforcingsMonthlytemperaturesEnum,&monthlytemperatures[0][imonth]);
-                        NewTemperatureInput->AddTimeInput(newmonthinput1,time+imonth/12.*yts);
-                        TriaP1Input* newmonthinput2 = new TriaP1Input(SurfaceforcingsPrecipitationEnum,&monthlyprec[0][imonth]);
-                        NewPrecipitationInput->AddTimeInput(newmonthinput2,time+imonth/12.*yts);
+                for(i=0;i<NUMVERTICES;i++) temp[i]=TemperaturesPresentday[i][imonth];
+                TriaP1Input* newmonthinput1 = new TriaP1Input(SurfaceforcingsMonthlytemperaturesEnum,&temp[0]);
+                //TriaP1Input* newmonthinput1 = new TriaP1Input(SurfaceforcingsMonthlytemperaturesEnum,&monthlytemperatures[0][imonth]);
+                NewTemperatureInput->AddTimeInput(newmonthinput1,time+imonth/12.*yts);
+                for(i=0;i<NUMVERTICES;i++) temp[i]=monthlyprec[i][imonth];
+                TriaP1Input* newmonthinput2 = new TriaP1Input(SurfaceforcingsPrecipitationEnum,&temp[0]);
+                NewPrecipitationInput->AddTimeInput(newmonthinput2,time+imonth/12.*yts);
         }
         this->inputs->AddInput(NewTemperatureInput);

issm/trunk-jpl/src/c/solutions/transient_core.cpp

@@ -139 +139 @@
                         InputToResultx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,SurfaceforcingsMassBalanceEnum);
                         InputToResultx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,MaskElementonfloatingiceEnum);
+                        bool isdelta18o;
+                        femmodel->parameters->FindParam(&isdelta18o,SurfaceforcingsIsdelta18oEnum);
+                        if(isdelta18o){
+                                InputToResultx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,SurfaceforcingsMonthlytemperaturesEnum);
+                                InputToResultx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,SurfaceforcingsPrecipitationEnum);
+                        }
                         RequestedOutputsx(femmodel->results,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,requested_outputs,numoutputs);