Changeset 12870


Ignore:
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:
2 edited

Legend:

Unmodified
Added
Removed
  • issm/trunk-jpl/src/c/classes/objects/Elements/Tria.cpp

    r12866 r12870  
    910910/*FUNCTION Tria::Delta18oParameterization{{{*/
    911911void  Tria::Delta18oParameterization(void){
     912
     913        int        i;
    912914        IssmDouble monthlytemperatures[NUMVERTICES][12],monthlyprec[NUMVERTICES][12];
    913915        IssmDouble TemperaturesPresentday[NUMVERTICES][12],TemperaturesLgm[NUMVERTICES][12];
    914916        IssmDouble PrecipitationsPresentday[NUMVERTICES][12];
     917        IssmDouble temp[NUMVERTICES];
    915918        IssmDouble Delta18oPresent,Delta18oLgm,Delta18oTime;
    916919        IssmDouble Delta18oSurfacePresent,Delta18oSurfaceLgm,Delta18oSurfaceTime;
    917920        IssmDouble time,yts,finaltime;
     921
     922        /*Recover parameters*/
    918923        this->parameters->FindParam(&time,TimeEnum);
    919924        this->parameters->FindParam(&yts,ConstantsYtsEnum);
     
    921926
    922927        /*Recover present day temperature and precipitation*/
    923         Input*     input=inputs->GetInput(SurfaceforcingsTemperaturesPresentdayEnum); _assert_(input);
    924         Input*     input2=inputs->GetInput(SurfaceforcingsTemperaturesLgmEnum); _assert_(input2);
    925         Input*     input3=inputs->GetInput(SurfaceforcingsPrecipitationsPresentdayEnum); _assert_(input3);
     928        Input* input=inputs->GetInput(SurfaceforcingsTemperaturesPresentdayEnum);    _assert_(input);
     929        Input* input2=inputs->GetInput(SurfaceforcingsTemperaturesLgmEnum);          _assert_(input2);
     930        Input* input3=inputs->GetInput(SurfaceforcingsPrecipitationsPresentdayEnum); _assert_(input3);
    926931        GaussTria* gauss=new GaussTria();
    927         for(int month=0;month<12;month++) {
    928                 for(int iv=0;iv<NUMVERTICES;iv++) {
     932        for(int month=0;month<12;month++){
     933                for(int iv=0;iv<NUMVERTICES;iv++){
    929934                        gauss->GaussVertex(iv);
    930935                        input->GetInputValue(&TemperaturesPresentday[iv][month],gauss,month/12.*yts);
     
    934939                }
    935940        }
     941
    936942        /*Recover delta18o and Delta18oSurface at present day, lgm and at time t*/
    937943        this->parameters->FindParam(&Delta18oPresent,SurfaceforcingsDelta18oEnum,finaltime);
    938         this->parameters->FindParam(&Delta18oLgm,SurfaceforcingsDelta18oEnum,(finaltime-(21000*yts)));
     944        this->parameters->FindParam(&Delta18oLgm,SurfaceforcingsDelta18oEnum,(finaltime-(21000*yts)));
    939945        this->parameters->FindParam(&Delta18oTime,SurfaceforcingsDelta18oEnum,time);
    940946        this->parameters->FindParam(&Delta18oSurfacePresent,SurfaceforcingsDelta18oSurfaceEnum,finaltime);
     
    950956                                        &monthlytemperatures[iv][0], &monthlyprec[iv][0]);
    951957        }
    952         //if(id==1)     printf(" monthlytemperatures1 %f\n",monthlytemperatures[iv][0]);
    953         //if(id==1)     printf(" monthlyprec1 %f\n",monthlyprec[iv][0]*yts);
     958
    954959        /*Update inputs*/
    955960        TransientInput* NewTemperatureInput = new TransientInput(SurfaceforcingsMonthlytemperaturesEnum);
    956961        TransientInput* NewPrecipitationInput = new TransientInput(SurfaceforcingsPrecipitationEnum);
    957962        for (int imonth=0;imonth<12;imonth++) {
    958                         TriaP1Input* newmonthinput1 = new TriaP1Input(SurfaceforcingsMonthlytemperaturesEnum,&TemperaturesPresentday[0][imonth]);
    959                         //TriaP1Input* newmonthinput1 = new TriaP1Input(SurfaceforcingsMonthlytemperaturesEnum,&monthlytemperatures[0][imonth]);
    960                         NewTemperatureInput->AddTimeInput(newmonthinput1,time+imonth/12.*yts);
    961                         TriaP1Input* newmonthinput2 = new TriaP1Input(SurfaceforcingsPrecipitationEnum,&monthlyprec[0][imonth]);
    962                         NewPrecipitationInput->AddTimeInput(newmonthinput2,time+imonth/12.*yts);
     963                for(i=0;i<NUMVERTICES;i++) temp[i]=TemperaturesPresentday[i][imonth];
     964                TriaP1Input* newmonthinput1 = new TriaP1Input(SurfaceforcingsMonthlytemperaturesEnum,&temp[0]);
     965                //TriaP1Input* newmonthinput1 = new TriaP1Input(SurfaceforcingsMonthlytemperaturesEnum,&monthlytemperatures[0][imonth]);
     966                NewTemperatureInput->AddTimeInput(newmonthinput1,time+imonth/12.*yts);
     967                for(i=0;i<NUMVERTICES;i++) temp[i]=monthlyprec[i][imonth];
     968                TriaP1Input* newmonthinput2 = new TriaP1Input(SurfaceforcingsPrecipitationEnum,&temp[0]);
     969                NewPrecipitationInput->AddTimeInput(newmonthinput2,time+imonth/12.*yts);
    963970        }
    964971        this->inputs->AddInput(NewTemperatureInput);
  • issm/trunk-jpl/src/c/solutions/transient_core.cpp

    r12832 r12870  
    139139                        InputToResultx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,SurfaceforcingsMassBalanceEnum);
    140140                        InputToResultx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,MaskElementonfloatingiceEnum);
     141                        bool isdelta18o;
     142                        femmodel->parameters->FindParam(&isdelta18o,SurfaceforcingsIsdelta18oEnum);
     143                        if(isdelta18o){
     144                                InputToResultx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,SurfaceforcingsMonthlytemperaturesEnum);
     145                                InputToResultx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,SurfaceforcingsPrecipitationEnum);
     146                        }
    141147                        RequestedOutputsx(femmodel->results,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,requested_outputs,numoutputs);
    142148
Note: See TracChangeset for help on using the changeset viewer.