Changeset 18855


Ignore:
Timestamp:
11/25/14 12:01:45 (10 years ago)
Author:
seroussi
Message:

CHG: moved Migrate Grounding line in Element.cpp

Location:
issm/trunk-jpl/src/c/classes/Elements
Files:
8 edited

Legend:

Unmodified
Added
Removed

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

    r18849 r18855  
    12011201
    12021202}/*}}}*/
     1203void       Element::MigrateGroundingLine(IssmDouble* phi_ungrounding){/*{{{*/
     1204
     1205        int         numvertices = this->GetNumberOfVertices();
     1206        int        i,migration_style;
     1207        IssmDouble bed_hydro,yts;
     1208        IssmDouble rho_water,rho_ice,density;
     1209        IssmDouble* melting = xNew<IssmDouble>(numvertices);
     1210        IssmDouble* phi     = xNew<IssmDouble>(numvertices);
     1211        IssmDouble* h       = xNew<IssmDouble>(numvertices);
     1212        IssmDouble* s       = xNew<IssmDouble>(numvertices);
     1213        IssmDouble* b       = xNew<IssmDouble>(numvertices);
     1214        IssmDouble* r       = xNew<IssmDouble>(numvertices);
     1215
     1216        /*Recover info at the vertices: */
     1217        parameters->FindParam(&migration_style,GroundinglineMigrationEnum);
     1218        parameters->FindParam(&yts,ConstantsYtsEnum);
     1219        GetInputListOnVertices(&h[0],ThicknessEnum);
     1220        GetInputListOnVertices(&s[0],SurfaceEnum);
     1221        GetInputListOnVertices(&b[0],BaseEnum);
     1222        GetInputListOnVertices(&r[0],BedEnum);
     1223        GetInputListOnVertices(&phi[0],MaskGroundediceLevelsetEnum);
     1224        rho_water   = matpar->GetRhoWater();
     1225        rho_ice     = matpar->GetRhoIce();
     1226        density     = rho_ice/rho_water;
     1227
     1228        if(migration_style == ContactEnum){
     1229                for(i = 0;i < numvertices;i++) phi[i] = phi_ungrounding[vertices[i]->Pid()];
     1230                this->AddInput(MaskGroundediceLevelsetEnum,&phi[0],P1Enum);
     1231
     1232                /*go through vertices, and update inputs, considering them to be TriaVertex type: */
     1233                for(i = 0;i < numvertices;i++){
     1234                        /*Ice shelf: if bed below bathymetry, impose it at the bathymetry and update surface, elso do nothing */
     1235                        if(phi[i] >= 0.){
     1236                                        b[i]  = r[i];
     1237                        }
     1238                }
     1239
     1240                /*Update inputs*/
     1241                this->AddInput(BaseEnum,&b[0],P1Enum);
     1242                return;
     1243        }
     1244
     1245        this->AddInput(MaskGroundediceLevelsetEnum,&phi[0],P1Enum);
     1246
     1247        /*go through vertices, and update inputs, considering them to be TriaVertex type: */
     1248        for(i=0;i<numvertices;i++){
     1249                /*Ice shelf: if bed below bathymetry, impose it at the bathymetry and update surface, elso do nothing */
     1250                if(phi[i]<=0.){
     1251                        if(b[i]<=r[i]){
     1252                                b[i]        = r[i];
     1253                                s[i]        = b[i]+h[i];
     1254                        }
     1255                }
     1256                /*Ice sheet: if hydrostatic bed above bathymetry, ice sheet starts to unground, elso do nothing */
     1257                /*Change only if AggressiveMigration or if the ice sheet is in contact with the ocean*/
     1258                else{ // phi>0
     1259                        bed_hydro=-density*h[i];
     1260                        if (bed_hydro>r[i]){
     1261                                /*Unground only if the element is connected to the ice shelf*/
     1262                                if(migration_style==AggressiveMigrationEnum || migration_style==SubelementMigrationEnum || migration_style==SubelementMigration2Enum){
     1263                                        s[i]        = (1-density)*h[i];
     1264                                        b[i]        = -density*h[i];
     1265                                }
     1266                                else if(migration_style==SoftMigrationEnum && phi_ungrounding[vertices[i]->Pid()]<0.){
     1267                                        s[i]        = (1-density)*h[i];
     1268                                        b[i]        = -density*h[i];
     1269                                }
     1270                                else{
     1271                                        if(migration_style!=SoftMigrationEnum) _error_("Error: migration should be Aggressive, Soft or Subelement");
     1272                                }
     1273                        }
     1274                }
     1275        }
     1276
     1277        /*Recalculate phi*/
     1278        for(i=0;i<numvertices;i++){
     1279                if(migration_style==SoftMigrationEnum){
     1280                        bed_hydro=-density*h[i];
     1281                        if(phi[i]<0. || bed_hydro<=r[i] || phi_ungrounding[vertices[i]->Pid()]<0.){
     1282                                phi[i]=h[i]+r[i]/density;
     1283                        }
     1284                }
     1285                else phi[i]=h[i]+r[i]/density;
     1286        }
     1287        this->AddInput(MaskGroundediceLevelsetEnum,&phi[0],P1Enum);
     1288
     1289        /*Update inputs*/
     1290        this->AddInput(SurfaceEnum,&s[0],P1Enum);
     1291        this->AddInput(BaseEnum,&b[0],P1Enum);
     1292
     1293        /*Delete*/
     1294        xDelete<IssmDouble>(melting);
     1295        xDelete<IssmDouble>(phi);
     1296        xDelete<IssmDouble>(r);
     1297        xDelete<IssmDouble>(b);
     1298        xDelete<IssmDouble>(s);
     1299        xDelete<IssmDouble>(h);
     1300
     1301}
     1302/*}}}*/
    12031303ElementVector* Element::NewElementVector(int approximation_enum){/*{{{*/
    12041304        return new ElementVector(nodes,this->GetNumberOfNodes(),this->parameters,approximation_enum);

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

    r18833 r18855  
    120120                bool       IsFloating();
    121121                void       LinearFloatingiceMeltingRate();
     122                void       MigrateGroundingLine(IssmDouble* sheet_ungrounding);
    122123                ElementVector*  NewElementVector(int approximation_enum=NoneApproximationEnum);
    123124                ElementMatrix*  NewElementMatrix(int approximation_enum=NoneApproximationEnum);
    … …  
    294295                virtual void UpdateConstraintsExtrudeFromTop(void)=0;
    295296
    296                 virtual void   MigrateGroundingLine(IssmDouble* sheet_ungrounding)=0;
    297297                virtual void   FSContactMigration(Vector<IssmDouble>* vertexgrounded,Vector<IssmDouble>* vertexfloating)=0;
    298298                virtual void   PotentialUngrounding(Vector<IssmDouble>* potential_sheet_ungrounding)=0;

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

    r18851 r18855  
    34283428}
    34293429/*}}}*/
    3430 void       Penta::MigrateGroundingLine(IssmDouble* phi_ungrounding){/*{{{*/
    3431 
    3432         int        i,migration_style;
    3433         IssmDouble bed_hydro,yts;
    3434         IssmDouble rho_water,rho_ice,density;
    3435         IssmDouble h[NUMVERTICES],s[NUMVERTICES],b[NUMVERTICES],r[NUMVERTICES];
    3436         IssmDouble melting[NUMVERTICES],phi[NUMVERTICES];
    3437 
    3438         if(!IsOnBase()) return;
    3439 
    3440         /*Recover info at the vertices: */
    3441         parameters->FindParam(&migration_style,GroundinglineMigrationEnum);
    3442         parameters->FindParam(&yts,ConstantsYtsEnum);
    3443         GetInputListOnVertices(&h[0],ThicknessEnum);
    3444         GetInputListOnVertices(&s[0],SurfaceEnum);
    3445         GetInputListOnVertices(&b[0],BaseEnum);
    3446         GetInputListOnVertices(&r[0],BedEnum);
    3447         GetInputListOnVertices(&phi[0],MaskGroundediceLevelsetEnum);
    3448         rho_water   = matpar->GetRhoWater();
    3449         rho_ice     = matpar->GetRhoIce();
    3450         density     = rho_ice/rho_water;
    3451 
    3452         /*go through vertices, and update inputs, considering them to be PentaVertex type: */
    3453         for(i=0;i<NUMVERTICES;i++){
    3454                 /*Ice shelf: if bed below bathymetry, impose it at the bathymetry and update surface, elso do nothing */
    3455                 if(phi[i]<=0){
    3456                         if(b[i]<=r[i]){
    3457                                 b[i]        = r[i];
    3458                                 s[i]        = b[i]+h[i];
    3459                         }
    3460                 }
    3461                 /*Ice sheet: if hydrostatic bed above bathymetry, ice sheet starts to unground, elso do nothing */
    3462                 /*Change only if AggressiveMigration or if the ice sheet is in contact with the ocean*/
    3463                 else{
    3464                         bed_hydro=-density*h[i];
    3465                         if(bed_hydro>r[i]){
    3466                                 /*Unground only if the element is connected to the ice shelf*/
    3467                                 if(migration_style==AggressiveMigrationEnum || migration_style==SubelementMigrationEnum || migration_style==SubelementMigration2Enum){
    3468                                         s[i]        = (1-density)*h[i];
    3469                                         b[i]        = -density*h[i];
    3470                                 }
    3471                                 else if(migration_style==SoftMigrationEnum && phi_ungrounding[vertices[i]->Pid()]<0.){
    3472                                         s[i]        = (1-density)*h[i];
    3473                                         b[i]        = -density*h[i];
    3474                                 }
    3475                                 else{
    3476                                         if(migration_style!=SoftMigrationEnum) _error_("Error: migration should be Aggressive, Soft or Subelement");
    3477                                 }
    3478                         }
    3479                 }
    3480         }
    3481 
    3482         /*Recalculate phi*/
    3483         for(i=0;i<NUMVERTICES;i++){
    3484                 if(migration_style==SoftMigrationEnum){
    3485                         bed_hydro=-density*h[i];
    3486                         if(phi[i]<0. || bed_hydro<=r[i] || phi_ungrounding[vertices[i]->Pid()]<0.){
    3487                                 phi[i]=h[i]+r[i]/density;
    3488                         }
    3489                 }
    3490                 else phi[i]=h[i]+r[i]/density;
    3491         }
    3492         this->inputs->AddInput(new PentaInput(MaskGroundediceLevelsetEnum,&phi[0],P1Enum));
    3493 
    3494         /*Update inputs*/
    3495         this->inputs->AddInput(new PentaInput(SurfaceEnum,&s[0],P1Enum));
    3496         this->inputs->AddInput(new PentaInput(BaseEnum,&b[0],P1Enum));
    3497 }
    3498 /*}}}*/
    34993430void       Penta::PotentialUngrounding(Vector<IssmDouble>* potential_ungrounding){/*{{{*/
    35003431

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

    r18842 r18855  
    135135                IssmDouble Masscon(IssmDouble* levelset){_error_("not implemented yet");};
    136136
    137                 void   MigrateGroundingLine(IssmDouble* sheet_ungrounding);
    138137                void   PotentialUngrounding(Vector<IssmDouble>* potential_sheet_ungrounding);
    139138                int    UpdatePotentialUngrounding(IssmDouble* potential_sheet_ungrounding,Vector<IssmDouble>* vec_nodes_on_iceshelf,IssmDouble* nodes_on_iceshelf);

  • issm/trunk-jpl/src/c/classes/Elements/Seg.h

    r18786 r18855  
    171171
    172172                void   PotentialUngrounding(Vector<IssmDouble>* potential_sheet_ungrounding){_error_("not implemented yet");};
    173                 void   MigrateGroundingLine(IssmDouble* sheet_ungrounding){_error_("not implemented yet");};
    174173                int    UpdatePotentialUngrounding(IssmDouble* vertices_potentially_ungrounding,Vector<IssmDouble>* vec_nodes_on_iceshelf,IssmDouble* nodes_on_iceshelf){_error_("not implemented yet");};
    175174                /*}}}*/

  • issm/trunk-jpl/src/c/classes/Elements/Tetra.h

    r18786 r18855  
    177177
    178178                void   PotentialUngrounding(Vector<IssmDouble>* potential_sheet_ungrounding){_error_("not implemented yet");};
    179                 void   MigrateGroundingLine(IssmDouble* sheet_ungrounding){_error_("not implemented yet");};
    180179                int    UpdatePotentialUngrounding(IssmDouble* vertices_potentially_ungrounding,Vector<IssmDouble>* vec_nodes_on_iceshelf,IssmDouble* nodes_on_iceshelf){_error_("not implemented yet");};
    181180                /*}}}*/

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

    r18836 r18855  
    32563256#endif
    32573257
    3258 void       Tria::MigrateGroundingLine(IssmDouble* phi_ungrounding){/*{{{*/
    3259 
    3260         int        i,migration_style;
    3261         bool       groundedelement = false;
    3262         IssmDouble bed_hydro,yts;
    3263         IssmDouble rho_water,rho_ice,density;
    3264         IssmDouble melting[NUMVERTICES],phi[NUMVERTICES];;
    3265         IssmDouble h[NUMVERTICES],s[NUMVERTICES],b[NUMVERTICES],r[NUMVERTICES];
    3266 
    3267         /*Recover info at the vertices: */
    3268         parameters->FindParam(&migration_style,GroundinglineMigrationEnum);
    3269         parameters->FindParam(&yts,ConstantsYtsEnum);
    3270         GetInputListOnVertices(&h[0],ThicknessEnum);
    3271         GetInputListOnVertices(&s[0],SurfaceEnum);
    3272         GetInputListOnVertices(&b[0],BaseEnum);
    3273         GetInputListOnVertices(&r[0],BedEnum);
    3274         GetInputListOnVertices(&phi[0],MaskGroundediceLevelsetEnum);
    3275         rho_water   = matpar->GetRhoWater();
    3276         rho_ice     = matpar->GetRhoIce();
    3277         density     = rho_ice/rho_water;
    3278 
    3279         if(migration_style == ContactEnum){
    3280                 for(i = 0;i < NUMVERTICES;i++) phi[i] = phi_ungrounding[vertices[i]->Pid()];
    3281                 this->inputs->AddInput(new TriaInput(MaskGroundediceLevelsetEnum,&phi[0],P1Enum));
    3282 
    3283                 /*go through vertices, and update inputs, considering them to be TriaVertex type: */
    3284                 for(i = 0;i < NUMVERTICES;i++){
    3285                         /*Ice shelf: if bed below bathymetry, impose it at the bathymetry and update surface, elso do nothing */
    3286                         if(phi[i] >= 0.){
    3287                                         b[i]  = r[i];
    3288                         }
    3289                 }
    3290 
    3291                 /*Update inputs*/
    3292                 this->inputs->AddInput(new TriaInput(BaseEnum,&b[0],P1Enum));
    3293                 return;
    3294         }
    3295 
    3296         this->inputs->AddInput(new TriaInput(MaskGroundediceLevelsetEnum,&phi[0],P1Enum));
    3297 
    3298         /*go through vertices, and update inputs, considering them to be TriaVertex type: */
    3299         for(i=0;i<NUMVERTICES;i++){
    3300                 /*Ice shelf: if bed below bathymetry, impose it at the bathymetry and update surface, elso do nothing */
    3301                 if(phi[i]<=0.){
    3302                         if(b[i]<=r[i]){
    3303                                 b[i]        = r[i];
    3304                                 s[i]        = b[i]+h[i];
    3305                         }
    3306                 }
    3307                 /*Ice sheet: if hydrostatic bed above bathymetry, ice sheet starts to unground, elso do nothing */
    3308                 /*Change only if AggressiveMigration or if the ice sheet is in contact with the ocean*/
    3309                 else{ // phi>0
    3310                         bed_hydro=-density*h[i];
    3311                         if (bed_hydro>r[i]){
    3312                                 /*Unground only if the element is connected to the ice shelf*/
    3313                                 if(migration_style==AggressiveMigrationEnum || migration_style==SubelementMigrationEnum || migration_style==SubelementMigration2Enum){
    3314                                         s[i]        = (1-density)*h[i];
    3315                                         b[i]        = -density*h[i];
    3316                                 }
    3317                                 else if(migration_style==SoftMigrationEnum && phi_ungrounding[vertices[i]->Pid()]<0.){
    3318                                         s[i]        = (1-density)*h[i];
    3319                                         b[i]        = -density*h[i];
    3320                                 }
    3321                                 else{
    3322                                         if(migration_style!=SoftMigrationEnum) _error_("Error: migration should be Aggressive, Soft or Subelement");
    3323                                 }
    3324                         }
    3325                 }
    3326         }
    3327 
    3328         /*Recalculate phi*/
    3329         for(i=0;i<NUMVERTICES;i++){
    3330                 if(migration_style==SoftMigrationEnum){
    3331                         bed_hydro=-density*h[i];
    3332                         if(phi[i]<0. || bed_hydro<=r[i] || phi_ungrounding[vertices[i]->Pid()]<0.){
    3333                                 phi[i]=h[i]+r[i]/density;
    3334                         }
    3335                 }
    3336                 else phi[i]=h[i]+r[i]/density;
    3337         }
    3338         this->inputs->AddInput(new TriaInput(MaskGroundediceLevelsetEnum,&phi[0],P1Enum));
    3339 
    3340         /*Update inputs*/
    3341         this->inputs->AddInput(new TriaInput(SurfaceEnum,&s[0],P1Enum));
    3342         this->inputs->AddInput(new TriaInput(BaseEnum,&b[0],P1Enum));
    3343 
    3344 }
    3345 /*}}}*/
    33463258void       Tria::PotentialUngrounding(Vector<IssmDouble>* potential_ungrounding){/*{{{*/
    33473259

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

    r18786 r18855  
    139139
    140140                void   PotentialUngrounding(Vector<IssmDouble>* potential_sheet_ungrounding);
    141                 void   MigrateGroundingLine(IssmDouble* sheet_ungrounding);
    142141                int    UpdatePotentialUngrounding(IssmDouble* vertices_potentially_ungrounding,Vector<IssmDouble>* vec_nodes_on_iceshelf,IssmDouble* nodes_on_iceshelf);
    143142
Note: See TracChangeset for help on using the changeset viewer.