Changeset 12748


Ignore:
Timestamp:
07/26/12 10:08:15 (13 years ago)
Author:
lemorzad
Message:

addind delta18o temperature and precipitation methode

Location:
issm/trunk-jpl
Files:
10 added
21 edited

Legend:

Unmodified
Added
Removed

  • issm/trunk-jpl/externalpackages/matlab/install.sh

    r12504 r12748  
    55
    66#Select or create a new simlink
    7 ln -s /usr/local/pkgs/matlab-7.6/ install
     7#ln -s /usr/local/pkgs/matlab-7.6/ install
    88#ln -s /discover/vis/mathworks/matlab_r2011b/ install
    99#ln -s /usr/local/matlab704/ install
    1010#ln -s /usr/local/matlab711/ install
    1111#ln -s /usr/local/matlab712/ install
     12
     13ln -s /home/kevin/perso.dir/matlab10 install
    1214
    1315# Macintosh (OSX) simlink

  • issm/trunk-jpl/src/c/EnumDefinitions/EnumDefinitions.h

    r12744 r12748  
    149149        SettingsResultsAsPatchesEnum,
    150150        SettingsWaitonlockEnum,
     151        SurfaceforcingsDelta18oEnum,
     152        SurfaceforcingsDelta18oTemperaturesPresentdayEnum,
     153        SurfaceforcingsDelta18oTemperaturesLgmEnum,
     154        SurfaceforcingsDelta18oSurfaceEnum,
     155        SurfaceforcingsIsdelta18oEnum,
     156        SurfaceforcingsPrecipitationsPresentdayEnum,
    151157        DebugProfilingEnum,
    152158        ProfilingCurrentMemEnum,

  • issm/trunk-jpl/src/c/Makefile.am

    r12744 r12748  
    212212                                        ./shared/Elements/GetNumberOfDofs.cpp\
    213213                                        ./shared/Elements/PddSurfaceMassBalance.cpp\
     214                                        ./shared/Elements/ComputeDelta18oTemperaturePrecipitation.cpp\
    214215                                        ./shared/String/sharedstring.h\
    215216                                        ./shared/Wrapper/wrappershared.h\
    … …  
    294295                                        ./modules/PositiveDegreeDayx/PositiveDegreeDayx.h\
    295296                                        ./modules/PositiveDegreeDayx/PositiveDegreeDayx.cpp\
     297                                        ./modules/Delta18oParameterizationx/Delta18oParameterizationx.h\
     298                                        ./modules/Delta18oParameterizationx/Delta18oParameterizationx.cpp\
    296299                                        ./modules/SmbGradientsx/SmbGradientsx.h\
    297300                                        ./modules/SmbGradientsx/SmbGradientsx.cpp\
    … …  
    818821matlabpetsc_sources= ./matlab/io/MatlabMatrixToPetscMatrix.cpp\
    819822                                         ./matlab/io/MatlabVectorToPetscVector.cpp
    820        
     823
    821824#}}}
    822825#Modules sources{{{

  • issm/trunk-jpl/src/c/modules/EnumToStringx/EnumToStringx.cpp

    r12744 r12748  
    154154                case SettingsResultsAsPatchesEnum : return "SettingsResultsAsPatches";
    155155                case SettingsWaitonlockEnum : return "SettingsWaitonlock";
     156                case SurfaceforcingsDelta18oEnum : return "SurfaceforcingsDelta18o";
     157                case SurfaceforcingsDelta18oTemperaturesPresentdayEnum : return "SurfaceforcingsDelta18oTemperaturesPresentday";
     158                case SurfaceforcingsDelta18oTemperaturesLgmEnum : return "SurfaceforcingsDelta18oTemperaturesLgm";
     159                case SurfaceforcingsDelta18oSurfaceEnum : return "SurfaceforcingsDelta18oSurface";
     160                case SurfaceforcingsIsdelta18oEnum : return "SurfaceforcingsIsdelta18o";
     161                case SurfaceforcingsPrecipitationsPresentdayEnum : return "SurfaceforcingsPrecipitationsPresentday";
    156162                case DebugProfilingEnum : return "DebugProfiling";
    157163                case ProfilingCurrentMemEnum : return "ProfilingCurrentMem";

  • issm/trunk-jpl/src/c/modules/ModelProcessorx/CreateParameters.cpp

    r12696 r12748  
    2525        Parameters *parameters       = NULL;
    2626        IssmDouble     *requestedoutputs = NULL;
    27        
     27        bool   isdelta18o;
     28        IssmDouble *Delta18oT, *Delta18oSurfaceT;
     29        int nlinesD18o, ncolsD18o, nlinesD18osurf, ncolsD18osurf;
     30
    2831        if(*pparameters)return; //do not create parameters twice!
    2932
    … …  
    9093        parameters->AddObject(iomodel->CopyConstantObject(InversionTaoEnum));
    9194        parameters->AddObject(iomodel->CopyConstantObject(SurfaceforcingsIspddEnum));
     95        parameters->AddObject(iomodel->CopyConstantObject(SurfaceforcingsIsdelta18oEnum));
    9296        parameters->AddObject(iomodel->CopyConstantObject(SurfaceforcingsIssmbgradientsEnum));
     97
     98        iomodel->Constant(&isdelta18o,SurfaceforcingsIsdelta18oEnum);
     99        if(isdelta18o){
     100                iomodel->FetchData(&Delta18oT,&nlinesD18o,&ncolsD18o,SurfaceforcingsDelta18oEnum);
     101                _assert_(nlinesD18o==2);
     102                parameters->AddObject(new TransientParam(SurfaceforcingsDelta18oEnum,&Delta18oT[0],&Delta18oT[ncolsD18o],ncolsD18o));
     103                 
     104                iomodel->FetchData(&Delta18oSurfaceT,&nlinesD18osurf,&ncolsD18osurf,SurfaceforcingsDelta18oSurfaceEnum);
     105                _assert_(nlinesD18osurf==2);
     106                parameters->AddObject(new TransientParam(SurfaceforcingsDelta18oSurfaceEnum,&Delta18oSurfaceT[0],&Delta18oSurfaceT[ncolsD18osurf],ncolsD18osurf));
     107        }
    93108
    94109        /*some parameters that did not come with the iomodel: */

  • issm/trunk-jpl/src/c/modules/ModelProcessorx/Prognostic/UpdateElementsPrognostic.cpp

    r12677 r12748  
    2020        int    stabilization;
    2121        bool   dakota_analysis;
    22         bool   ispdd;
    2322        bool   issmbgradients;
    2423
    … …  
    2928        iomodel->Constant(&dakota_analysis,QmuIsdakotaEnum);
    3029        iomodel->FetchData(1,MeshElementsEnum);
    31         iomodel->Constant(&ispdd,SurfaceforcingsIspddEnum);
    3230        iomodel->Constant(&issmbgradients,SurfaceforcingsIssmbgradientsEnum);
    3331
    … …  
    6967                iomodel->FetchDataToInput(elements,TemperatureEnum);
    7068        }
    71         if(ispdd){
    72           iomodel->FetchDataToInput(elements,VyEnum);
    73           iomodel->FetchDataToInput(elements,ThermalSpctemperatureEnum);
    74           iomodel->FetchDataToInput(elements,SurfaceforcingsPrecipitationEnum);
    75           iomodel->FetchDataToInput(elements,SurfaceforcingsMonthlytemperaturesEnum);
     69        if(issmbgradients){
     70                iomodel->FetchDataToInput(elements,SurfaceforcingsHcEnum);
     71                iomodel->FetchDataToInput(elements,SurfaceforcingsSmbPosMaxEnum);
     72                iomodel->FetchDataToInput(elements,SurfaceforcingsSmbPosMinEnum);
     73                iomodel->FetchDataToInput(elements,SurfaceforcingsAPosEnum);
     74                iomodel->FetchDataToInput(elements,SurfaceforcingsBPosEnum);
     75                iomodel->FetchDataToInput(elements,SurfaceforcingsANegEnum);
     76                iomodel->FetchDataToInput(elements,SurfaceforcingsBNegEnum);
    7677        }
    77         if(issmbgradients){
    78           iomodel->FetchDataToInput(elements,SurfaceforcingsHcEnum);
    79           iomodel->FetchDataToInput(elements,SurfaceforcingsSmbPosMaxEnum);
    80           iomodel->FetchDataToInput(elements,SurfaceforcingsSmbPosMinEnum);
    81           iomodel->FetchDataToInput(elements,SurfaceforcingsAPosEnum);
    82           iomodel->FetchDataToInput(elements,SurfaceforcingsBPosEnum);
    83           iomodel->FetchDataToInput(elements,SurfaceforcingsANegEnum);
    84           iomodel->FetchDataToInput(elements,SurfaceforcingsBNegEnum);
    85         }
    86         else{
    87                 iomodel->FetchDataToInput(elements,SurfaceforcingsMassBalanceEnum);
    88         }
     78        //else{
     79        //      iomodel->FetchDataToInput(elements,SurfaceforcingsMassBalanceEnum);
     80        //}
    8981
    9082        /*Free data: */

  • issm/trunk-jpl/src/c/modules/ModelProcessorx/Transient/UpdateElementsTransient.cpp

    r9340 r12748  
    1616void    UpdateElementsTransient(Elements* elements, Parameters* parameters,IoModel* iomodel,int analysis_counter,int analysis_type){
    1717
    18         /*nothing for now: */
     18        bool   ispdd;
     19        bool   isdelta18o;
     20        IssmDouble *size, Delta18oTimeSerie,Delta18oSurfaceTimeSerie ;
    1921
     22        /*Fetch data needed: */
     23        iomodel->Constant(&ispdd,SurfaceforcingsIspddEnum);
     24        iomodel->Constant(&isdelta18o,SurfaceforcingsIsdelta18oEnum);
     25
     26        if(ispdd){
     27                iomodel->FetchDataToInput(elements,ThermalSpctemperatureEnum);
     28                if(isdelta18o){
     29                          iomodel->FetchDataToInput(elements,SurfaceforcingsDelta18oEnum);
     30                          iomodel->FetchDataToInput(elements,SurfaceforcingsDelta18oSurfaceEnum);
     31 
     32                          iomodel->FetchDataToInput(elements,SurfaceforcingsDelta18oTemperaturesLgmEnum);
     33                          iomodel->FetchDataToInput(elements,SurfaceforcingsDelta18oTemperaturesPresentdayEnum);
     34                          iomodel->FetchDataToInput(elements,SurfaceforcingsPrecipitationsPresentdayEnum);
     35                }
     36                else{
     37                          iomodel->FetchDataToInput(elements,SurfaceforcingsPrecipitationEnum);
     38                          iomodel->FetchDataToInput(elements,SurfaceforcingsMonthlytemperaturesEnum);
     39                }
     40        }
     41        else{
     42                iomodel->FetchDataToInput(elements,SurfaceforcingsMassBalanceEnum);
     43        }
    2044}

  • issm/trunk-jpl/src/c/modules/StringToEnumx/StringToEnumx.cpp

    r12744 r12748  
    158158              else if (strcmp(name,"SettingsResultsAsPatches")==0) return SettingsResultsAsPatchesEnum;
    159159              else if (strcmp(name,"SettingsWaitonlock")==0) return SettingsWaitonlockEnum;
     160              else if (strcmp(name,"SurfaceforcingsDelta18o")==0) return SurfaceforcingsDelta18oEnum;
     161              else if (strcmp(name,"SurfaceforcingsDelta18oTemperaturesPresentday")==0) return SurfaceforcingsDelta18oTemperaturesPresentdayEnum;
     162              else if (strcmp(name,"SurfaceforcingsDelta18oTemperaturesLgm")==0) return SurfaceforcingsDelta18oTemperaturesLgmEnum;
     163              else if (strcmp(name,"SurfaceforcingsDelta18oSurface")==0) return SurfaceforcingsDelta18oSurfaceEnum;
     164              else if (strcmp(name,"SurfaceforcingsIsdelta18o")==0) return SurfaceforcingsIsdelta18oEnum;
     165              else if (strcmp(name,"SurfaceforcingsPrecipitationsPresentday")==0) return SurfaceforcingsPrecipitationsPresentdayEnum;
    160166              else if (strcmp(name,"DebugProfiling")==0) return DebugProfilingEnum;
    161167              else if (strcmp(name,"ProfilingCurrentMem")==0) return ProfilingCurrentMemEnum;
    … …  
    255261              else if (strcmp(name,"ControlInput")==0) return ControlInputEnum;
    256262              else if (strcmp(name,"DatasetInput")==0) return DatasetInputEnum;
    257               else if (strcmp(name,"DofIndexing")==0) return DofIndexingEnum;
     263         else stage=3;
     264   }
     265   if(stage==3){
     266              if (strcmp(name,"DofIndexing")==0) return DofIndexingEnum;
    258267              else if (strcmp(name,"DoubleInput")==0) return DoubleInputEnum;
    259268              else if (strcmp(name,"DoubleMatArrayParam")==0) return DoubleMatArrayParamEnum;
    … …  
    261270              else if (strcmp(name,"DoubleParam")==0) return DoubleParamEnum;
    262271              else if (strcmp(name,"DoubleVecParam")==0) return DoubleVecParamEnum;
    263          else stage=3;
    264    }
    265    if(stage==3){
    266               if (strcmp(name,"Element")==0) return ElementEnum;
     272              else if (strcmp(name,"Element")==0) return ElementEnum;
    267273              else if (strcmp(name,"ElementResult")==0) return ElementResultEnum;
    268274              else if (strcmp(name,"ExternalResult")==0) return ExternalResultEnum;
    … …  
    378384              else if (strcmp(name,"DragCoefficientAbsGradient")==0) return DragCoefficientAbsGradientEnum;
    379385              else if (strcmp(name,"TransientInput")==0) return TransientInputEnum;
    380               else if (strcmp(name,"Outputfilename")==0) return OutputfilenameEnum;
     386         else stage=4;
     387   }
     388   if(stage==4){
     389              if (strcmp(name,"Outputfilename")==0) return OutputfilenameEnum;
    381390              else if (strcmp(name,"Waterfraction")==0) return WaterfractionEnum;
    382391              else if (strcmp(name,"Watercolumn")==0) return WatercolumnEnum;
    … …  
    384393              else if (strcmp(name,"ViscousHeating")==0) return ViscousHeatingEnum;
    385394              else if (strcmp(name,"QmuTemperature")==0) return QmuTemperatureEnum;
    386          else stage=4;
    387    }
    388    if(stage==4){
    389               if (strcmp(name,"HydrologyWaterVx")==0) return HydrologyWaterVxEnum;
     395              else if (strcmp(name,"HydrologyWaterVx")==0) return HydrologyWaterVxEnum;
    390396              else if (strcmp(name,"HydrologyWaterVy")==0) return HydrologyWaterVyEnum;
    391397              else if (strcmp(name,"StressTensor")==0) return StressTensorEnum;

  • issm/trunk-jpl/src/c/modules/modules.h

    r12677 r12748  
    2525#include "./CreateNodalConstraintsx/CreateNodalConstraintsx.h"
    2626#include "./DakotaResponsesx/DakotaResponsesx.h"
     27#include "./Delta18oParameterizationx/Delta18oParameterizationx.h"
    2728#include "./DragCoefficientAbsGradientx/DragCoefficientAbsGradientx.h"
    2829#include "./ElementConnectivityx/ElementConnectivityx.h"

  • issm/trunk-jpl/src/c/objects/Elements/Element.h

    r12704 r12748  
    7272                virtual void   PotentialSheetUngrounding(Vector* potential_sheet_ungrounding)=0;
    7373                virtual void   PositiveDegreeDay(IssmDouble* pdds,IssmDouble* pds,IssmDouble signorm)=0;
     74                virtual void   Delta18oParameterization(void)=0;
    7475                virtual void   SmbGradients()=0;
    7576                virtual int    UpdatePotentialSheetUngrounding(IssmDouble* potential_sheet_ungrounding,Vector* vec_nodes_on_iceshelf,IssmDouble* nodes_on_iceshelf)=0;

  • issm/trunk-jpl/src/c/objects/Elements/Penta.cpp

    r12704 r12748  
    683683        this->results=new Results();
    684684
     685}
     686/*}}}*/
     687/*FUNCTION Penta::Delta18oParameterization{{{*/
     688void  Penta::Delta18oParameterization(void){
     689
     690        IssmDouble monthlytemperatures[NUMVERTICES][12],monthlyprec[NUMVERTICES][12];
     691        IssmDouble Delta18oTemperaturesPresentday[NUMVERTICES][12],Delta18oTemperaturesLgm[NUMVERTICES][12];
     692        IssmDouble PrecipitationsPresentday[NUMVERTICES][12];
     693        IssmDouble Delta18oPresent,Delta18oLgm,Delta18oTime;
     694        IssmDouble Delta18oSurfacePresent,Delta18oSurfaceLgm,Delta18oSurfaceTime;
     695        IssmDouble time,yts,finaltime;
     696        this->parameters->FindParam(&time,TimeEnum);
     697        this->parameters->FindParam(&yts,ConstantsYtsEnum);
     698        this->parameters->FindParam(&finaltime,TimesteppingFinalTimeEnum);
     699
     700        /*Recover present day temperature and precipitation*/
     701        Input*     input=inputs->GetInput(SurfaceforcingsDelta18oTemperaturesPresentdayEnum); _assert_(input);
     702        Input*     input2=inputs->GetInput(SurfaceforcingsDelta18oTemperaturesLgmEnum); _assert_(input2);
     703        Input*     input3=inputs->GetInput(SurfaceforcingsPrecipitationsPresentdayEnum); _assert_(input3);
     704        GaussPenta* gauss=new GaussPenta();
     705        for(int month=0;month<12;month++) {
     706          for(int iv=0;iv<NUMVERTICES;iv++) {
     707                gauss->GaussVertex(iv);
     708                input->GetInputValue(&Delta18oTemperaturesPresentday[iv][month],gauss,month/12.*yts);
     709                input2->GetInputValue(&Delta18oTemperaturesLgm[iv][month],gauss,month/12.*yts);
     710                input3->GetInputValue(&PrecipitationsPresentday[iv][month],gauss,month/12.*yts);
     711                monthlyprec[iv][month]=monthlyprec[iv][month]*yts; // convertion to m/yr
     712          }
     713        }
     714       
     715        /*Recover delta18o and Delta18oSurface at present day, lgm and at time t*/
     716        this->parameters->FindParam(&Delta18oPresent,SurfaceforcingsDelta18oEnum,finaltime*yts);
     717        this->parameters->FindParam(&Delta18oLgm,SurfaceforcingsDelta18oEnum,(finaltime-21000)*yts);
     718        this->parameters->FindParam(&Delta18oTime,SurfaceforcingsDelta18oEnum,time*yts);
     719        this->parameters->FindParam(&Delta18oSurfacePresent,SurfaceforcingsDelta18oSurfaceEnum,finaltime*yts);
     720        this->parameters->FindParam(&Delta18oSurfaceLgm,SurfaceforcingsDelta18oSurfaceEnum,(finaltime-21000)*yts);
     721        this->parameters->FindParam(&Delta18oSurfaceTime,SurfaceforcingsDelta18oSurfaceEnum,time*yts);
     722
     723        /*Compute the temperature and precipitation*/
     724        for(int iv=0;iv<NUMVERTICES;iv++){
     725          ComputeDelta18oTemperaturePrecipitation(Delta18oSurfacePresent, Delta18oSurfaceLgm, Delta18oSurfaceTime,
     726                                              Delta18oPresent, Delta18oLgm, Delta18oTime,
     727                                              &PrecipitationsPresentday[iv][0],
     728                                              &Delta18oTemperaturesLgm[iv][0], &Delta18oTemperaturesPresentday[iv][0],
     729                                              &monthlytemperatures[iv][0], &monthlyprec[iv][0]);
     730        }
     731 
     732        /*Update inputs*/
     733        TransientInput* NewTemperatureInput = new TransientInput(SurfaceforcingsMonthlytemperaturesEnum);
     734        TransientInput* NewPrecipitationInput = new TransientInput(SurfaceforcingsPrecipitationEnum);
     735        for (int imonth=0;imonth<12;imonth++) {
     736                for(int iv=0;iv<NUMVERTICES;iv++) {
     737                        PentaP1Input* newmonthinput1 = new PentaP1Input(SurfaceforcingsMonthlytemperaturesEnum,&monthlytemperatures[iv][imonth]);
     738                        NewTemperatureInput->AddTimeInput(newmonthinput1,imonth/12.*yts);
     739                        PentaP1Input* newmonthinput2 = new PentaP1Input(SurfaceforcingsPrecipitationEnum,&monthlyprec[iv][imonth]);
     740                        NewPrecipitationInput->AddTimeInput(newmonthinput2,imonth/12.*yts);
     741                }
     742        }
     743        this->inputs->AddInput(NewTemperatureInput);
     744        this->inputs->AddInput(NewPrecipitationInput);
    685745}
    686746/*}}}*/

  • issm/trunk-jpl/src/c/objects/Elements/Penta.h

    r12704 r12748  
    8484                void   CreatePVector(Vector* pf);
    8585                void   CreateJacobianMatrix(Matrix* Jff);
     86                void   Delta18oParameterization(void);
    8687                void   DeleteResults(void);
    8788                int    GetNodeIndex(Node* node);

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

    r12704 r12748  
    906906        this->results=new Results();
    907907
     908}
     909/*}}}*/
     910/*FUNCTION Tria::Delta18oParameterization{{{*/
     911void  Tria::Delta18oParameterization(void){
     912
     913        IssmDouble monthlytemperatures[NUMVERTICES][12],monthlyprec[NUMVERTICES][12];
     914        IssmDouble Delta18oTemperaturesPresentday[NUMVERTICES][12],Delta18oTemperaturesLgm[NUMVERTICES][12];
     915        IssmDouble PrecipitationsPresentday[NUMVERTICES][12];
     916        IssmDouble Delta18oPresent,Delta18oLgm,Delta18oTime;
     917        IssmDouble Delta18oSurfacePresent,Delta18oSurfaceLgm,Delta18oSurfaceTime;
     918        IssmDouble time,yts,finaltime;
     919        this->parameters->FindParam(&time,TimeEnum);
     920        this->parameters->FindParam(&yts,ConstantsYtsEnum);
     921        this->parameters->FindParam(&finaltime,TimesteppingFinalTimeEnum);
     922
     923        /*Recover present day temperature and precipitation*/
     924        Input*     input=inputs->GetInput(SurfaceforcingsDelta18oTemperaturesPresentdayEnum); _assert_(input);
     925        Input*     input2=inputs->GetInput(SurfaceforcingsDelta18oTemperaturesLgmEnum); _assert_(input2);
     926        Input*     input3=inputs->GetInput(SurfaceforcingsPrecipitationsPresentdayEnum); _assert_(input3);
     927        GaussTria* gauss=new GaussTria();
     928        for(int month=0;month<12;month++) {
     929                for(int iv=0;iv<NUMVERTICES;iv++) {
     930                        gauss->GaussVertex(iv);
     931                        input->GetInputValue(&Delta18oTemperaturesPresentday[iv][month],gauss,month/12.*yts);
     932                        input2->GetInputValue(&Delta18oTemperaturesLgm[iv][month],gauss,month/12.*yts);
     933                        input3->GetInputValue(&PrecipitationsPresentday[iv][month],gauss,month/12.*yts);
     934                        monthlyprec[iv][month]=monthlyprec[iv][month]*yts; // convertion to m/yr
     935                }
     936        }
     937
     938        /*Recover delta18o and Delta18oSurface at present day, lgm and at time t*/
     939        this->parameters->FindParam(&Delta18oPresent,SurfaceforcingsDelta18oEnum,finaltime*yts);
     940        this->parameters->FindParam(&Delta18oLgm,SurfaceforcingsDelta18oEnum,(finaltime-21000)*yts);
     941        this->parameters->FindParam(&Delta18oTime,SurfaceforcingsDelta18oEnum,time*yts);
     942        this->parameters->FindParam(&Delta18oSurfacePresent,SurfaceforcingsDelta18oSurfaceEnum,finaltime*yts);
     943        this->parameters->FindParam(&Delta18oSurfaceLgm,SurfaceforcingsDelta18oSurfaceEnum,(finaltime-21000)*yts);
     944        this->parameters->FindParam(&Delta18oSurfaceTime,SurfaceforcingsDelta18oSurfaceEnum,time*yts);
     945   
     946        /*Compute the temperature and precipitation*/
     947        for(int iv=0;iv<NUMVERTICES;iv++){
     948          ComputeDelta18oTemperaturePrecipitation(Delta18oSurfacePresent, Delta18oSurfaceLgm, Delta18oSurfaceTime,
     949                                              Delta18oPresent, Delta18oLgm, Delta18oTime,
     950                                              &PrecipitationsPresentday[iv][0],
     951                                              &Delta18oTemperaturesLgm[iv][0], &Delta18oTemperaturesPresentday[iv][0],
     952                                              &monthlytemperatures[iv][0], &monthlyprec[iv][0]);
     953        }
     954 
     955        /*Update inputs*/
     956        TransientInput* NewTemperatureInput = new TransientInput(SurfaceforcingsMonthlytemperaturesEnum);
     957        TransientInput* NewPrecipitationInput = new TransientInput(SurfaceforcingsPrecipitationEnum);
     958        for (int imonth=0;imonth<12;imonth++) {
     959                for(int iv=0;iv<NUMVERTICES;iv++) {
     960                        TriaP1Input* newmonthinput1 = new TriaP1Input(SurfaceforcingsMonthlytemperaturesEnum,&monthlytemperatures[iv][imonth]);
     961                        NewTemperatureInput->AddTimeInput(newmonthinput1,imonth/12.*yts);
     962                        TriaP1Input* newmonthinput2 = new TriaP1Input(SurfaceforcingsPrecipitationEnum,&monthlyprec[iv][imonth]);
     963                        NewPrecipitationInput->AddTimeInput(newmonthinput2,imonth/12.*yts);
     964                }
     965        }
     966        this->inputs->AddInput(NewTemperatureInput);
     967        this->inputs->AddInput(NewPrecipitationInput);
    908968}
    909969/*}}}*/

  • issm/trunk-jpl/src/c/objects/Elements/Tria.h

    r12704 r12748  
    8080                void   CreatePVector(Vector* pf);
    8181                void   CreateJacobianMatrix(Matrix* Jff);
     82                void   Delta18oParameterization(void);
    8283                int    GetNodeIndex(Node* node);
    8384                int    Sid();

  • issm/trunk-jpl/src/c/shared/Elements/PddSurfaceMassBalance.cpp

    r12705 r12748  
    4343  IssmDouble frzndd = 0.; 
    4444 
    45   IssmDouble tstar;                        // monthly mean surface temp
     45  IssmDouble tstar;          // monthly mean surface temp
    4646  IssmDouble Tsum= 0.;       // average summer (JJA) temperature
    4747  IssmDouble Tsurf = 0.;     // average annual temperature   

  • issm/trunk-jpl/src/c/shared/Elements/elements.h

    r12704 r12748  
    1313IssmDouble Paterson(IssmDouble temperature);
    1414IssmDouble Arrhenius(IssmDouble temperature,IssmDouble depth,IssmDouble n);
    15 IssmDouble PddSurfaceMassBlance(IssmDouble* monthlytemperatures,  IssmDouble* monthlyprec, IssmDouble* pdds, IssmDouble* pds, IssmDouble signorm, IssmDouble yts, IssmDouble h, IssmDouble s, IssmDouble rho_ice, IssmDouble rho_water);
     15IssmDouble PddSurfaceMassBlance(IssmDouble* monthlytemperatures,  IssmDouble* monthlyprec, IssmDouble* pdds, IssmDouble* pds,
     16                                IssmDouble signorm, IssmDouble yts, IssmDouble h, IssmDouble s, IssmDouble rho_ice, IssmDouble rho_water);
     17void ComputeDelta18oTemperaturePrecipitation(IssmDouble Delta18oSurfacePresent, IssmDouble Delta18oSurfaceLgm, IssmDouble Delta18oSurfaceTime,
     18                                     IssmDouble Delta18oPresent, IssmDouble Delta18oLgm, IssmDouble Delta18oTime,
     19                                     IssmDouble* PrecipitationsPresentday,
     20                                     IssmDouble* Delta18oTemperaturesLgm, IssmDouble* Delta18oTemperaturesPresentday,
     21                                             IssmDouble* monthlytemperaturesout, IssmDouble* monthlyprecout);
    1622void   GetVerticesCoordinates(IssmDouble* xyz,  Node** nodes, int numvertices);
    1723int    GetNumberOfDofs( Node** nodes,int numnodes,int setenum,int approximation_enum);

  • issm/trunk-jpl/src/c/solutions/prognostic_core.cpp

    r12677 r12748  
    1616        /*parameters: */
    1717        bool save_results;
    18         bool ispdd;
    1918        bool issmbgradients;
    2019
    … …  
    2423        /*recover parameters: */
    2524        femmodel->parameters->FindParam(&save_results,SaveResultsEnum);
    26         femmodel->parameters->FindParam(&ispdd,SurfaceforcingsIspddEnum);
    2725        femmodel->parameters->FindParam(&issmbgradients,SurfaceforcingsIssmbgradientsEnum);
    2826
    29         if(ispdd){
    30           if(VerboseSolution()) _pprintLine_("   call positive degree day module");
    31           PositiveDegreeDayx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters);
    32         }
    3327
    3428        if(issmbgradients){

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

    r12635 r12748  
    2424        /*parameters: */
    2525        IssmDouble starttime,finaltime,dt,yts;
    26         bool   isdiagnostic,isprognostic,isthermal,isgroundingline,isenthalpy;
     26        bool   isdiagnostic,isprognostic,isthermal,isgroundingline,isenthalpy,ispdd,isdelta18o;
    2727        bool   save_results,dakota_analysis;
    2828        bool   time_adapt=false;
    … …  
    5353        femmodel->parameters->FindParam(&numoutputs,TransientNumRequestedOutputsEnum);
    5454        if(numoutputs) femmodel->parameters->FindParam(&requested_outputs,&numoutputs,TransientRequestedOutputsEnum);
     55        femmodel->parameters->FindParam(&ispdd,SurfaceforcingsIspddEnum);
     56        femmodel->parameters->FindParam(&isdelta18o,SurfaceforcingsIsdelta18oEnum);
    5557
    5658        /*initialize: */
    … …  
    106108                        #endif
    107109                }
     110               
     111                if(ispdd){
     112                        if(isdelta18o){
     113                                if(VerboseSolution()) _pprintLine_("   call Delta18oParametrization module");
     114                                Delta18oParameterizationx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters);
     115                        }
     116                        if(VerboseSolution()) _pprintLine_("   call positive degree day module");
     117                        PositiveDegreeDayx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters);
     118                }
    108119
    109120                if(isdiagnostic){

  • issm/trunk-jpl/src/m/classes/surfaceforcings.m

    r12677 r12748  
    1010                ispdd = 0;
    1111                issmbgradients = 0;
     12                isdelta18o = 0;
    1213           hc = NaN;
    1314                smb_pos_max = NaN;
    … …  
    1819                b_neg = NaN;
    1920                monthlytemperatures = NaN;
     21                delta18o = NaN;
     22                delta18o_surface = NaN;
     23                temperatures_presentday = NaN;
     24                temperatures_lgm = NaN;
     25                precipitations_presentday = NaN;
    2026        end
    2127        methods
    … …  
    3339                  obj.ispdd=0;
    3440                  obj.issmbgradients=0;
     41                  obj.isdelta18o=0;
    3542
    3643                end % }}}
    … …  
    4148                                checkfield(md,'surfaceforcings.issmbgradients','numel',1,'values',[0 1]);
    4249                                if(obj.ispdd)
    43                                         md = checkfield(md,'surfaceforcings.monthlytemperatures','forcing',1,'NaN',1);
     50                                        if(obj.isdelta18o==0)
     51                                                md = checkfield(md,'surfaceforcings.monthlytemperatures','forcing',1,'NaN',1);
     52                                                md = checkfield(md,'surfaceforcings.precipitation','forcing',1,'NaN',1);
     53                                        else
     54                                                md = checkfield(md,'surfaceforcings.delta18o','NaN',1);
     55                                                md = checkfield(md,'surfaceforcings.delta18o_surface','NaN',1);
     56                                                md = checkfield(md,'surfaceforcings.temperatures_presentday','size',[md.mesh.numberofvertices+1 12],'NaN',1);
     57                                                md = checkfield(md,'surfaceforcings.temperatures_lgm','size',[md.mesh.numberofvertices+1 12],'NaN',1);
     58                                                md = checkfield(md,'surfaceforcings.precipitations_presentday','size',[md.mesh.numberofvertices+1 12],'NaN',1);
     59                                        end
    4460                                elseif(obj.issmbgradients)
    4561                                        checkfield(md,'surfaceforcings.hc','forcing',1,'NaN',1);
    … …  
    6177                        disp(sprintf('   surface forcings parameters:'));
    6278
    63                         fielddisplay(obj,'precipitation','surface precipitation [m/yr water eq]');
    6479                        fielddisplay(obj,'mass_balance','surface mass balance [m/yr ice eq]');
    6580                        fielddisplay(obj,'ispdd','is pdd activated (0 or 1, default is 0)');
    66                         fielddisplay(obj,'monthlytemperatures','monthly surface temperatures required if pdd is activated');
     81                        fielddisplay(obj,'isdelta18o','is temperature and precipitation delta18o parametrisation activated (0 or 1, default is 0)');
     82                        fielddisplay(obj,'monthlytemperatures','monthly surface temperatures [Kelvin], required if pdd is activated and delta18o not activated');
     83                        fielddisplay(obj,'precipitation','surface precipitation [m/yr water eq]');
     84                        fielddisplay(obj,'temperatures_presentday','monthly present day surface temperatures [Kelvin], required if pdd is activated and delta18o activated');
     85                        fielddisplay(obj,'temperatures_lgm','monthly LGM surface temperatures [Kelvin], required if pdd is activated and delta18o activated');
     86                        fielddisplay(obj,'precipitations_presentday','monthly surface precipitation [m/yr water eq], required if pdd is activated and delta18o activated');
     87                        fielddisplay(obj,'delta18o_surface','surface elevation of the delta18o site, required if pdd is activated and delta18o activated');
    6788                        fielddisplay(obj,'issmbgradients','is smb gradients method activated (0 or 1, default is 0)');
    6889                        fielddisplay(obj,'hc',' elevation of intersection between accumulation and ablation regime required if smb gradients is activated');
    … …  
    79100                        WriteData(fid,'object',obj,'fieldname','mass_balance','format','DoubleMat','mattype',1);
    80101                        WriteData(fid,'object',obj,'fieldname','ispdd','format','Boolean');
     102                        WriteData(fid,'object',obj,'fieldname','isdelta18o','format','Boolean');
    81103                        if obj.ispdd,
     104                          if obj.isdelta18o
     105                                WriteData(fid,'object',obj,'fieldname','temperatures_presentday','format','DoubleMat','mattype',1);
     106                                WriteData(fid,'object',obj,'fieldname','temperatures_lgm','format','DoubleMat','mattype',1);
     107                                WriteData(fid,'object',obj,'fieldname','precipitations_presentday','format','DoubleMat','mattype',1);
     108                                WriteData(fid,'object',obj,'fieldname','delta18o_surface','format','DoubleMat','mattype',1);
     109                          else
    82110                                WriteData(fid,'object',obj,'fieldname','monthlytemperatures','format','DoubleMat','mattype',1);
     111                                WriteData(fid,'object',obj,'fieldname','precipitation','format','DoubleMat','mattype',1);
     112                          end
    83113                        end
    84114                        WriteData(fid,'object',obj,'fieldname','issmbgradients','format','Boolean');

  • issm/trunk-jpl/src/m/enum/EnumDefinitions.py

    r12745 r12748  
    13771377        return StringToEnum('SettingsWaitonlock')
    13781378
     1379def SurfaceforcingsDelta18oEnum():
     1380        """
     1381        SURFACEFORCINGSDELTA18OENUM - Enum of SurfaceforcingsDelta18o
     1382
     1383           Usage:
     1384              macro=SurfaceforcingsDelta18oEnum()
     1385        """
     1386
     1387        return StringToEnum('SurfaceforcingsDelta18o')
     1388
     1389def SurfaceforcingsDelta18oTemperaturesPresentdayEnum():
     1390        """
     1391        SURFACEFORCINGSDELTA18OTEMPERATURESPRESENTDAYENUM - Enum of SurfaceforcingsDelta18oTemperaturesPresentday
     1392
     1393           Usage:
     1394              macro=SurfaceforcingsDelta18oTemperaturesPresentdayEnum()
     1395        """
     1396
     1397        return StringToEnum('SurfaceforcingsDelta18oTemperaturesPresentday')
     1398
     1399def SurfaceforcingsDelta18oTemperaturesLgmEnum():
     1400        """
     1401        SURFACEFORCINGSDELTA18OTEMPERATURESLGMENUM - Enum of SurfaceforcingsDelta18oTemperaturesLgm
     1402
     1403           Usage:
     1404              macro=SurfaceforcingsDelta18oTemperaturesLgmEnum()
     1405        """
     1406
     1407        return StringToEnum('SurfaceforcingsDelta18oTemperaturesLgm')
     1408
     1409def SurfaceforcingsDelta18oSurfaceEnum():
     1410        """
     1411        SURFACEFORCINGSDELTA18OSURFACEENUM - Enum of SurfaceforcingsDelta18oSurface
     1412
     1413           Usage:
     1414              macro=SurfaceforcingsDelta18oSurfaceEnum()
     1415        """
     1416
     1417        return StringToEnum('SurfaceforcingsDelta18oSurface')
     1418
     1419def SurfaceforcingsIsdelta18oEnum():
     1420        """
     1421        SURFACEFORCINGSISDELTA18OENUM - Enum of SurfaceforcingsIsdelta18o
     1422
     1423           Usage:
     1424              macro=SurfaceforcingsIsdelta18oEnum()
     1425        """
     1426
     1427        return StringToEnum('SurfaceforcingsIsdelta18o')
     1428
     1429def SurfaceforcingsPrecipitationsPresentdayEnum():
     1430        """
     1431        SURFACEFORCINGSPRECIPITATIONSPRESENTDAYENUM - Enum of SurfaceforcingsPrecipitationsPresentday
     1432
     1433           Usage:
     1434              macro=SurfaceforcingsPrecipitationsPresentdayEnum()
     1435        """
     1436
     1437        return StringToEnum('SurfaceforcingsPrecipitationsPresentday')
     1438
    13791439def DebugProfilingEnum():
    13801440        """
    … …  
    45054565        """
    45064566
    4507         return 449
    4508 
     4567        return 455
     4568

  • issm/trunk-jpl/src/m/enum/MaximumNumberOfEnums.m

    r12745 r12748  
    99%      macro=MaximumNumberOfEnums()
    1010
    11 macro=449;
     11macro=455;
Note: See TracChangeset for help on using the changeset viewer.