Changeset 14175


Ignore:
Timestamp:
12/17/12 16:24:36 (12 years ago)
Author:
helsen
Message:

update SMB gradients code

Location:
issm/trunk-jpl/src
Files:
8 edited

Legend:

Unmodified
Added
Removed

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

    r14161 r14175  
    196196        SurfaceforcingsIssmbgradientsEnum,
    197197        SurfaceforcingsMonthlytemperaturesEnum,
    198         SurfaceforcingsHcEnum,
    199198        SurfaceforcingsHrefEnum,
    200199        SurfaceforcingsSmbrefEnum,
    201         SurfaceforcingsSmbPosMaxEnum,
    202         SurfaceforcingsSmbPosMinEnum,
    203         SurfaceforcingsAPosEnum,
    204200        SurfaceforcingsBPosEnum,
    205         SurfaceforcingsANegEnum,
    206201        SurfaceforcingsBNegEnum,
    207202        ThermalMaxiterEnum,

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

    r14004 r14175  
    27272727
    27282728        /*Recover SmbGradients*/
    2729         GetInputListOnVertices(&Hc[0],SurfaceforcingsHcEnum);
    27302729        GetInputListOnVertices(&Href[0],SurfaceforcingsHrefEnum);
    27312730        GetInputListOnVertices(&Smbref[0],SurfaceforcingsSmbrefEnum);
    2732         GetInputListOnVertices(&smb_pos_max[0],SurfaceforcingsSmbPosMaxEnum);
    2733         GetInputListOnVertices(&smb_pos_min[0],SurfaceforcingsSmbPosMinEnum);
    2734         GetInputListOnVertices(&a_pos[0],SurfaceforcingsAPosEnum);
    27352731        GetInputListOnVertices(&b_pos[0],SurfaceforcingsBPosEnum);
    2736         GetInputListOnVertices(&a_neg[0],SurfaceforcingsANegEnum);
    27372732        GetInputListOnVertices(&b_neg[0],SurfaceforcingsBNegEnum);
    27382733
    … …  
    27472742   // loop over all vertices
    27482743 for(i=0;i<NUMVERTICES;i++){
    2749      if(s[i]>Hc[i]){
    2750                   if(Href[i]>Hc[i]){smb[i]=Smbref[i]+b_pos[i]*(s[i]-Href[i]);}
    2751                   if(Href[i]<=Hc[i]){smb[i]=a_pos[i]+b_pos[i]*s[i];}
    2752                   if(smb[i]>smb_pos_max[i]){smb[i]=smb_pos_max[i];}
    2753                   if(smb[i]<smb_pos_min[i]){smb[i]=smb_pos_min[i];}
     2744     if(Smbref[i]>0){
     2745                  smb[i]=Smbref[i]+b_pos[i]*(s[i]-Href[i]);
    27542746          }
    27552747          else{
    2756                   if(Href[i]>Hc[i]){smb[i]=a_neg[i]+b_neg[i]*s[i];}
    2757                   if(Href[i]<=Hc[i]){smb[i]=Smbref[i]+b_neg[i]*(s[i]-Href[i]);}
     2748                  smb[i]=Smbref[i]+b_neg[i]*(s[i]-Href[i]);
    27582749          }
    27592750          smb[i]=smb[i]/rho_ice;      // SMB in m/y ice         

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

    r14159 r14175  
    23252325   IssmDouble h[NUMVERTICES];                                   // ice thickness (m)           
    23262326        IssmDouble s[NUMVERTICES];                                      // surface elevation (m)
    2327         IssmDouble a_pos[NUMVERTICES];                          // Hs-SMB relation parameter
    23282327        IssmDouble b_pos[NUMVERTICES];                          // Hs-SMB relation parameter
    2329         IssmDouble a_neg[NUMVERTICES];                          // Hs-SMB relation parameter
    23302328        IssmDouble b_neg[NUMVERTICES];                          // Hs-SMB relation paremeter
    2331         IssmDouble Hc[NUMVERTICES];                                     // elevation of transition between accumulation regime and ablation regime
    23322329        IssmDouble Href[NUMVERTICES];                                   // reference elevation from which deviations are used to calculate the SMB adjustment
    23332330        IssmDouble Smbref[NUMVERTICES];                         // reference SMB to which deviations are added
    2334         IssmDouble smb_pos_max[NUMVERTICES];            // maximum SMB value in the accumulation regime
    2335         IssmDouble smb_pos_min[NUMVERTICES];            // minimum SMB value in the accumulation regime
    23362331   IssmDouble rho_water;                   // density of fresh water
    23372332        IssmDouble rho_ice;                     // density of ice
    … …  
    23412336
    23422337        /*Recover SmbGradients*/
    2343         GetInputListOnVertices(&Hc[0],SurfaceforcingsHcEnum);
    23442338        GetInputListOnVertices(&Href[0],SurfaceforcingsHrefEnum);
    23452339        GetInputListOnVertices(&Smbref[0],SurfaceforcingsSmbrefEnum);
    2346         GetInputListOnVertices(&smb_pos_max[0],SurfaceforcingsSmbPosMaxEnum);
    2347         GetInputListOnVertices(&smb_pos_min[0],SurfaceforcingsSmbPosMinEnum);
    2348         GetInputListOnVertices(&a_pos[0],SurfaceforcingsAPosEnum);
    23492340        GetInputListOnVertices(&b_pos[0],SurfaceforcingsBPosEnum);
    2350         GetInputListOnVertices(&a_neg[0],SurfaceforcingsANegEnum);
    23512341        GetInputListOnVertices(&b_neg[0],SurfaceforcingsBNegEnum);
    23522342
    … …  
    23612351   // loop over all vertices
    23622352   for(i=0;i<NUMVERTICES;i++){
    2363      if(s[i]>Hc[i]){
    2364                   if(Href[i]>Hc[i]){smb[i]=Smbref[i]+b_pos[i]*(s[i]-Href[i]);}
    2365                   if(Href[i]<=Hc[i]){smb[i]=a_pos[i]+b_pos[i]*s[i];}
    2366                   if(smb[i]>smb_pos_max[i]){smb[i]=smb_pos_max[i];}
    2367                   if(smb[i]<smb_pos_min[i]){smb[i]=smb_pos_min[i];}
     2353     if(Smbref[i]>0){
     2354                  smb[i]=Smbref[i]+b_pos[i]*(s[i]-Href[i]);
    23682355          }
    23692356          else{
    2370                   if(Href[i]>Hc[i]){smb[i]=a_neg[i]+b_neg[i]*s[i];}
    2371                   if(Href[i]<=Hc[i]){smb[i]=Smbref[i]+b_neg[i]*(s[i]-Href[i]);}
     2357                  smb[i]=Smbref[i]+b_neg[i]*(s[i]-Href[i]);
    23722358          }
    23732359          smb[i]=smb[i]/rho_ice;      // SMB in m/y ice         
    2374   /*   printf("s %e \n",s[i]);
    2375      printf("Hsref %e \n",Href[i]);
    2376      printf("Hc %e \n",Hc[i]);
    2377      printf("Smbref %e \n",Smbref[i]);
    2378      printf("b_neg %e \n",b_neg[i]);
    2379      printf("smb %e \n",smb[i]);
    2380           _error_("stop-in-code"); */
    23812360                }  //end of the loop over the vertices
    23822361          /*Update inputs*/

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

    r14161 r14175  
    201201                case SurfaceforcingsIssmbgradientsEnum : return "SurfaceforcingsIssmbgradients";
    202202                case SurfaceforcingsMonthlytemperaturesEnum : return "SurfaceforcingsMonthlytemperatures";
    203                 case SurfaceforcingsHcEnum : return "SurfaceforcingsHc";
    204203                case SurfaceforcingsHrefEnum : return "SurfaceforcingsHref";
    205204                case SurfaceforcingsSmbrefEnum : return "SurfaceforcingsSmbref";
    206                 case SurfaceforcingsSmbPosMaxEnum : return "SurfaceforcingsSmbPosMax";
    207                 case SurfaceforcingsSmbPosMinEnum : return "SurfaceforcingsSmbPosMin";
    208                 case SurfaceforcingsAPosEnum : return "SurfaceforcingsAPos";
    209205                case SurfaceforcingsBPosEnum : return "SurfaceforcingsBPos";
    210                 case SurfaceforcingsANegEnum : return "SurfaceforcingsANeg";
    211206                case SurfaceforcingsBNegEnum : return "SurfaceforcingsBNeg";
    212207                case ThermalMaxiterEnum : return "ThermalMaxiter";

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

    r13521 r14175  
    7272        }
    7373        if(issmbgradients){
    74                 iomodel->FetchDataToInput(elements,SurfaceforcingsHcEnum);
    7574                iomodel->FetchDataToInput(elements,SurfaceforcingsHrefEnum);
    7675                iomodel->FetchDataToInput(elements,SurfaceforcingsSmbrefEnum);
    77                 iomodel->FetchDataToInput(elements,SurfaceforcingsSmbPosMaxEnum);
    78                 iomodel->FetchDataToInput(elements,SurfaceforcingsSmbPosMinEnum);
    79                 iomodel->FetchDataToInput(elements,SurfaceforcingsAPosEnum);
    8076                iomodel->FetchDataToInput(elements,SurfaceforcingsBPosEnum);
    81                 iomodel->FetchDataToInput(elements,SurfaceforcingsANegEnum);
    8277                iomodel->FetchDataToInput(elements,SurfaceforcingsBNegEnum);
    8378        }

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

    r14161 r14175  
    205205              else if (strcmp(name,"SurfaceforcingsIssmbgradients")==0) return SurfaceforcingsIssmbgradientsEnum;
    206206              else if (strcmp(name,"SurfaceforcingsMonthlytemperatures")==0) return SurfaceforcingsMonthlytemperaturesEnum;
    207               else if (strcmp(name,"SurfaceforcingsHc")==0) return SurfaceforcingsHcEnum;
    208207              else if (strcmp(name,"SurfaceforcingsHref")==0) return SurfaceforcingsHrefEnum;
    209208              else if (strcmp(name,"SurfaceforcingsSmbref")==0) return SurfaceforcingsSmbrefEnum;
    210               else if (strcmp(name,"SurfaceforcingsSmbPosMax")==0) return SurfaceforcingsSmbPosMaxEnum;
    211               else if (strcmp(name,"SurfaceforcingsSmbPosMin")==0) return SurfaceforcingsSmbPosMinEnum;
    212               else if (strcmp(name,"SurfaceforcingsAPos")==0) return SurfaceforcingsAPosEnum;
    213209              else if (strcmp(name,"SurfaceforcingsBPos")==0) return SurfaceforcingsBPosEnum;
    214               else if (strcmp(name,"SurfaceforcingsANeg")==0) return SurfaceforcingsANegEnum;
    215210              else if (strcmp(name,"SurfaceforcingsBNeg")==0) return SurfaceforcingsBNegEnum;
    216211              else if (strcmp(name,"ThermalMaxiter")==0) return ThermalMaxiterEnum;
    … …  
    261256              else if (strcmp(name,"NoneAnalysis")==0) return NoneAnalysisEnum;
    262257              else if (strcmp(name,"PrognosticAnalysis")==0) return PrognosticAnalysisEnum;
    263          else stage=3;
    264    }
    265    if(stage==3){
    266               if (strcmp(name,"PrognosticSolution")==0) return PrognosticSolutionEnum;
     258              else if (strcmp(name,"PrognosticSolution")==0) return PrognosticSolutionEnum;
    267259              else if (strcmp(name,"SteadystateSolution")==0) return SteadystateSolutionEnum;
    268260              else if (strcmp(name,"SurfaceSlopeAnalysis")==0) return SurfaceSlopeAnalysisEnum;
    269261              else if (strcmp(name,"SurfaceSlopeSolution")==0) return SurfaceSlopeSolutionEnum;
    270262              else if (strcmp(name,"SurfaceSlopeXAnalysis")==0) return SurfaceSlopeXAnalysisEnum;
    271               else if (strcmp(name,"SurfaceSlopeYAnalysis")==0) return SurfaceSlopeYAnalysisEnum;
     263         else stage=3;
     264   }
     265   if(stage==3){
     266              if (strcmp(name,"SurfaceSlopeYAnalysis")==0) return SurfaceSlopeYAnalysisEnum;
    272267              else if (strcmp(name,"ThermalAnalysis")==0) return ThermalAnalysisEnum;
    273268              else if (strcmp(name,"ThermalSolution")==0) return ThermalSolutionEnum;
    … …  
    384379              else if (strcmp(name,"QmuMelting")==0) return QmuMeltingEnum;
    385380              else if (strcmp(name,"ResetPenalties")==0) return ResetPenaltiesEnum;
    386          else stage=4;
    387    }
    388    if(stage==4){
    389               if (strcmp(name,"SegmentOnIceShelf")==0) return SegmentOnIceShelfEnum;
     381              else if (strcmp(name,"SegmentOnIceShelf")==0) return SegmentOnIceShelfEnum;
    390382              else if (strcmp(name,"SurfaceAbsVelMisfit")==0) return SurfaceAbsVelMisfitEnum;
    391383              else if (strcmp(name,"SurfaceArea")==0) return SurfaceAreaEnum;
    392384              else if (strcmp(name,"SurfaceAverageVelMisfit")==0) return SurfaceAverageVelMisfitEnum;
    393385              else if (strcmp(name,"SurfaceLogVelMisfit")==0) return SurfaceLogVelMisfitEnum;
    394               else if (strcmp(name,"SurfaceLogVxVyMisfit")==0) return SurfaceLogVxVyMisfitEnum;
     386         else stage=4;
     387   }
     388   if(stage==4){
     389              if (strcmp(name,"SurfaceLogVxVyMisfit")==0) return SurfaceLogVxVyMisfitEnum;
    395390              else if (strcmp(name,"SurfaceRelVelMisfit")==0) return SurfaceRelVelMisfitEnum;
    396391              else if (strcmp(name,"SurfaceSlopeX")==0) return SurfaceSlopeXEnum;
    … …  
    507502              else if (strcmp(name,"QmuErrName")==0) return QmuErrNameEnum;
    508503              else if (strcmp(name,"QmuInName")==0) return QmuInNameEnum;
    509          else stage=5;
    510    }
    511    if(stage==5){
    512               if (strcmp(name,"QmuOutName")==0) return QmuOutNameEnum;
     504              else if (strcmp(name,"QmuOutName")==0) return QmuOutNameEnum;
    513505              else if (strcmp(name,"Regular")==0) return RegularEnum;
    514506              else if (strcmp(name,"Scaled")==0) return ScaledEnum;
    515507              else if (strcmp(name,"Separate")==0) return SeparateEnum;
    516508              else if (strcmp(name,"Sset")==0) return SsetEnum;
    517               else if (strcmp(name,"Verbose")==0) return VerboseEnum;
     509         else stage=5;
     510   }
     511   if(stage==5){
     512              if (strcmp(name,"Verbose")==0) return VerboseEnum;
    518513              else if (strcmp(name,"TriangleInterp")==0) return TriangleInterpEnum;
    519514              else if (strcmp(name,"BilinearInterp")==0) return BilinearInterpEnum;

  • issm/trunk-jpl/src/c/shared/Numerics/UnitConversion.cpp

    r13622 r14175  
    6565                case SurfaceforcingsPrecipitationEnum:       scale=yts;break;     //m/yr
    6666                case SurfaceforcingsMassBalanceEnum:         scale=yts;break;     //m/yr
    67                 case SurfaceforcingsSmbPosMaxEnum:                              scale=yts;break;     //m/yr
    68                 case SurfaceforcingsSmbPosMinEnum:                              scale=yts;break;     //m/yr
    69                 case SurfaceforcingsAPosEnum:                                           scale=yts;break;     //m/yr
    7067                case SurfaceforcingsBPosEnum:                                           scale=yts;break;     //m/yr
    71                 case SurfaceforcingsANegEnum:                                           scale=yts;break;     //m/yr
    7268                case SurfaceforcingsBNegEnum:                                           scale=yts;break;     //m/yr
    7369                case SurfaceforcingsSmbrefEnum:                                 scale=yts;break;     //m/yr

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

    r13646 r14175  
    1111                issmbgradients = 0;
    1212                isdelta18o = 0;
    13                 hc = NaN;
    1413                href = NaN;
    1514                smbref = NaN;
    16                 smb_pos_max = NaN;
    17                 smb_pos_min = NaN;
    18                 a_pos = NaN;
    1915                b_pos = NaN;
    20                 a_neg = NaN;
    2116                b_neg = NaN;
    2217                monthlytemperatures = NaN;
    … …  
    6156                                        end
    6257                                elseif(obj.issmbgradients)
    63                                         md = checkfield(md,'surfaceforcings.hc','forcing',1,'NaN',1);
    6458                                        md = checkfield(md,'surfaceforcings.href','forcing',1,'NaN',1);
    6559                                        md = checkfield(md,'surfaceforcings.smbref','forcing',1,'NaN',1);
    66                                         md = checkfield(md,'surfaceforcings.smb_pos_max','forcing',1,'NaN',1);
    67                                         md = checkfield(md,'surfaceforcings.smb_pos_min','forcing',1,'NaN',1);
    68                                         md = checkfield(md,'surfaceforcings.a_pos','forcing',1,'NaN',1);
    6960                                        md = checkfield(md,'surfaceforcings.b_pos','forcing',1,'NaN',1);
    70                                         md = checkfield(md,'surfaceforcings.a_neg','forcing',1,'NaN',1);
    7161                                        md = checkfield(md,'surfaceforcings.b_neg','forcing',1,'NaN',1);
    7262                                else
    … …  
    9282                        fielddisplay(obj,'delta18o_surface','surface elevation of the delta18o site, required if pdd is activated and delta18o activated');
    9383                        fielddisplay(obj,'issmbgradients','is smb gradients method activated (0 or 1, default is 0)');
    94                         fielddisplay(obj,'hc',' elevation of intersection between accumulation and ablation regime required if smb gradients is activated');
    9584                        fielddisplay(obj,'href',' reference elevation from which deviation is used to calculate SMB adjustment in smb gradients method');
    9685                        fielddisplay(obj,'smbref',' reference smb from which deviation is calculated in smb gradients method');
    97                         fielddisplay(obj,'smb_pos_max',' maximum value of positive smb required if smb gradients is activated');
    98                         fielddisplay(obj,'smb_pos_min',' minimum value of positive smb required if smb gradients is activated');
    99                         fielddisplay(obj,'a_pos',' intercept of hs - smb regression line for accumulation regime required if smb gradients is activated');
    10086                        fielddisplay(obj,'b_pos',' slope of hs - smb regression line for accumulation regime required if smb gradients is activated');
    101                         fielddisplay(obj,'a_neg',' intercept of hs - smb regression line for ablation regime required if smb gradients is activated');
    10287                        fielddisplay(obj,'b_neg',' slope of hs - smb regression line for ablation regime required if smb gradients is activated');
    10388
    … …  
    122107                        WriteData(fid,'object',obj,'fieldname','issmbgradients','format','Boolean');
    123108                        if obj.issmbgradients,
    124                                 WriteData(fid,'object',obj,'fieldname','hc','format','DoubleMat','mattype',1);
    125109                                WriteData(fid,'object',obj,'fieldname','href','format','DoubleMat','mattype',1);
    126110                                WriteData(fid,'object',obj,'fieldname','smbref','format','DoubleMat','mattype',1);
    127                                 WriteData(fid,'object',obj,'fieldname','smb_pos_max','format','DoubleMat','mattype',1);
    128                                 WriteData(fid,'object',obj,'fieldname','smb_pos_min','format','DoubleMat','mattype',1);
    129                                 WriteData(fid,'object',obj,'fieldname','a_pos','format','DoubleMat','mattype',1);
    130111                                WriteData(fid,'object',obj,'fieldname','b_pos','format','DoubleMat','mattype',1);
    131                                 WriteData(fid,'object',obj,'fieldname','a_neg','format','DoubleMat','mattype',1);
    132112                                WriteData(fid,'object',obj,'fieldname','b_neg','format','DoubleMat','mattype',1);
    133113                        end
Note: See TracChangeset for help on using the changeset viewer.