Changeset 17353


Ignore:
Timestamp:
02/26/14 11:04:28 (11 years ago)
Author:
bdef
Message:

CHG: cleaning water transfer from tria and related stuff

Location:
issm/trunk-jpl/src/c
Files:
16 edited

Legend:

Unmodified
Added
Removed

  • issm/trunk-jpl/src/c/analyses/HydrologyDCEfficientAnalysis.cpp

    r17351 r17353  
    5656        iomodel->FetchDataToInput(elements,HydrologydcEplInitialThicknessEnum);
    5757        iomodel->FetchDataToInput(elements,HydrologydcSedimentTransmitivityEnum);
    58 
    59         iomodel->FetchDataToInput(elements,HydrologydcMaskEplactiveNodeEnum);
    6058       
    6159        elements->InputDuplicate(HydrologydcEplInitialThicknessEnum,HydrologydcEplThicknessEnum);
    … …  
    186184                        /*Transfer EPL part*/
    187185                        transfer=GetHydrologyKMatrixTransfer(basalelement,gauss);
    188                         basalelement->NodalFunctions(&basis[0],gauss);
    189186                        D_scalar=transfer*gauss->weight*Jdet*dt;
    190187                        TripleMultiply(basis,numnodes,1,0,
    … …  
    296293        return pe;
    297294}/*}}}*/
     295
    298296void HydrologyDCEfficientAnalysis::GetSolutionFromInputs(Vector<IssmDouble>* solution,Element* element){/*{{{*/
    299297        element->GetSolutionFromInputsOneDof(solution,EplHeadEnum);
    300298}/*}}}*/
     299
    301300void HydrologyDCEfficientAnalysis::InputUpdateFromSolution(IssmDouble* solution,Element* element){/*{{{*/
    302301
    … …  
    322321        /*Fetch dof list and allocate solution vector*/
    323322        basalelement->GetDofList(&doflist,NoneApproximationEnum,GsetEnum);
    324        
    325323        IssmDouble* eplHeads    = xNew<IssmDouble>(numnodes);
    326324
    … …  
    329327                eplHeads[i]=solution[doflist[i]];
    330328                if(xIsNan<IssmDouble>(eplHeads[i])) _error_("NaN found in solution vector");
    331 
    332329        }
    333330        /*Add input to the element: */
    … …  
    336333        /*Free ressources:*/
    337334        xDelete<IssmDouble>(eplHeads);
    338         //      xDelete<IssmDouble>(eplOldHeads);
    339335        xDelete<int>(doflist);
    340336        if(meshtype!=Mesh2DhorizontalEnum){basalelement->DeleteMaterials(); delete basalelement;};
    … …  
    355351        return rho_freshwater*g*epl_porosity*(water_compressibility+(epl_compressibility/epl_porosity));                 
    356352}/*}}}*/
     353
    357354IssmDouble HydrologyDCEfficientAnalysis::SedimentStoring(Element* element){/*{{{*/
    358355        IssmDouble rho_freshwater           = element->GetMaterialParameter(MaterialsRhoFreshwaterEnum);
    … …  
    364361        return rho_freshwater*g*sediment_porosity*sediment_thickness*(water_compressibility+(sediment_compressibility/sediment_porosity));               
    365362}/*}}}*/
     363
    366364IssmDouble HydrologyDCEfficientAnalysis::GetHydrologyDCInefficientHmax(Element* element, Gauss* gauss){/*{{{*/
    367365        int        hmax_flag;
    … …  
    369367        IssmDouble rho_ice,rho_water;
    370368        IssmDouble thickness,bed;
     369
    371370        /*Get the flag to the limitation method*/
    372371        element->FindParam(&hmax_flag,HydrologydcSedimentlimitFlagEnum);

  • issm/trunk-jpl/src/c/analyses/HydrologyDCInefficientAnalysis.cpp

    r17349 r17353  
    8181
    8282        iomodel->FetchDataToInput(elements,ThicknessEnum);
    83         iomodel->FetchDataToInput(elements,SurfaceEnum);
    8483        iomodel->FetchDataToInput(elements,BedEnum);
    8584        iomodel->FetchDataToInput(elements,MaskIceLevelsetEnum);
    8685        iomodel->FetchDataToInput(elements,MeshElementonbedEnum);
    87         iomodel->FetchDataToInput(elements,MeshElementonsurfaceEnum);
    8886        iomodel->FetchDataToInput(elements,BasalforcingsMeltingRateEnum);
    8987        iomodel->FetchDataToInput(elements,SedimentHeadEnum);
    … …  
    235233                                if(active_element){
    236234                                        transfer=GetHydrologyKMatrixTransfer(basalelement,gauss);
    237                                        
    238235                                        basalelement->NodalFunctions(&basis[0],gauss);
    239236                                        D_scalar=transfer*gauss->weight*Jdet*dt;

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

    r17338 r17353  
    494494                                name==HydrologydcEplThicknessEnum ||
    495495                                name==HydrologydcMaskEplactiveNodeEnum ||
    496                                 name==MeshVertexonbedEnum ||
    497                                 name==WaterTransferEnum
     496                                name==MeshVertexonbedEnum
    498497
    499498                                ) {

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

    r17349 r17353  
    283283                virtual void UpdateConstraintsExtrudeFromTop(void)=0;
    284284
    285                 //              virtual void GetHydrologyDCInefficientHmax(IssmDouble* ph_max,int index)=0;
    286285                virtual void GetHydrologyDCInefficientHmax(IssmDouble* ph_max, Node* innode)=0;
    287                 virtual void GetHydrologyTransfer(Vector<IssmDouble>* transfer)=0;
    288286                virtual void HydrologyEPLGetMask(Vector<IssmDouble>* mask)=0;
    289287                virtual void HydrologyEPLGetActive(Vector<IssmDouble>* active)=0;

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

    r17345 r17353  
    41794179}
    41804180/*}}}*/
    4181 /*FUNCTION Penta::GetHydrologyTransfer{{{*/
    4182 void  Penta::GetHydrologyTransfer(Vector<IssmDouble>* transfer){
    4183 
    4184         if (!IsOnBed()) return;
    4185 
    4186         Tria* tria=(Tria*)SpawnTria(0); //lower face is 0, upper face is 1.
    4187         tria->GetHydrologyTransfer(transfer);
    4188         delete tria->material; delete tria;
    4189 }
    4190 /*}}}*/
    41914181/*FUNCTION Penta::GetSolutionFromInputsOneDof {{{*/
    41924182void Penta::GetSolutionFromInputsOneDof(Vector<IssmDouble>* solution, int enum_type){

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

    r17349 r17353  
    243243                ElementMatrix* CreateEPLDomainMassMatrix(void);
    244244                void           GetHydrologyDCInefficientHmax(IssmDouble* ph_max, Node* innode);
    245                 //              void           GetHydrologyDCInefficientHmax(IssmDouble* ph_max,int index){_error_("not implemented yet");};
    246                 void           GetHydrologyTransfer(Vector<IssmDouble>* transfer);
    247245                void           HydrologyEPLGetActive(Vector<IssmDouble>* active_vec);
    248246                void           HydrologyEPLGetMask(Vector<IssmDouble>* vec_mask);

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

    r17349 r17353  
    151151
    152152                void    GetHydrologyDCInefficientHmax(IssmDouble* ph_max, Node* innode){_error_("not implemented yet");};
    153                 //              void    GetHydrologyDCInefficientHmax(IssmDouble* ph_max,int index){_error_("not implemented yet");};
    154                 void    GetHydrologyTransfer(Vector<IssmDouble>* transfer){_error_("not implemented yet");};
    155153                void    HydrologyEPLGetActive(Vector<IssmDouble>* active_vec){_error_("not implemented yet");};
    156154                void    HydrologyEPLGetMask(Vector<IssmDouble>* vec_mask){_error_("not implemented yet");};

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

    r17350 r17353  
    44994499}
    45004500/*}}}*/
    4501 /*FUNCTION Tria::GetHydrologyTransfer{{{*/
    4502 /*This thing should be useless and deleted soon*/
    4503 
    4504  void  Tria::GetHydrologyTransfer(Vector<IssmDouble>* transfer){
    4505 
    4506 /*      const int  numdof   = NDOF1 *NUMVERTICES; */
    4507 /*      int        *doflist = NULL; */
    4508 /*      int        analysis_type; */
    4509 /*      bool       isefficientlayer; */
    4510 /*      bool       active_element; */
    4511 /*      int        transfermethod; */
    4512 /*      IssmDouble leakage,h_max,dt,test; */
    4513 /*      IssmDouble wh_trans,sed_thick,relaxed; */
    4514 /*      IssmDouble epl_specificstoring,sedstoring; */
    4515 /*      IssmDouble activeEpl[numdof],epl_thickness[numdof],old_epl_thickness[numdof]; */
    4516 /*      IssmDouble epl_head[numdof],sed_head[numdof]; */
    4517 /*      IssmDouble epl_head_old[numdof],sed_head_old[numdof]; */
    4518 /*      IssmDouble preceding_transfer[numdof],sed_trans[numdof]; */
    4519 
    4520 /*      Input* active_element_input=NULL; */
    4521 
    4522 
    4523 /*      parameters->FindParam(&analysis_type,AnalysisTypeEnum); */
    4524                
    4525 /*      GetDofList(&doflist,NoneApproximationEnum,GsetEnum); */
    4526 
    4527 /*      /\*Get the flag to know if the efficient layer is present*\/ */
    4528 /*      this->parameters->FindParam(&isefficientlayer,HydrologydcIsefficientlayerEnum); */
    4529 
    4530 /*      if(isefficientlayer){ */
    4531 /*              /\*Also get the flag to the transfer method*\/ */
    4532 /*              this->parameters->FindParam(&transfermethod,HydrologydcTransferFlagEnum); */
    4533 /*              this->parameters->FindParam(&dt,TimesteppingTimeStepEnum);               */
    4534 /*              /\*Switch between the different transfer methods cases*\/ */
    4535 /*              switch(transfermethod){ */
    4536 /*              case 0: */
    4537 /*                      /\*Just keepping the transfer to zero, should be OK with the initial value of transfer*\/ */
    4538 /*                      break; */
    4539 /*              case 1: */
    4540        
    4541 /*                      active_element_input=inputs->GetInput(HydrologydcMaskEplactiveEltEnum); _assert_(active_element_input); */
    4542 /*                      active_element_input->GetInputValue(&active_element); */
    4543 
    4544 /*                      GetInputListOnVertices(&sed_head[0],SedimentHeadEnum); */
    4545 /*                      GetInputListOnVertices(&sed_head_old[0],SedimentHeadOldEnum); */
    4546 /*                      GetInputListOnVertices(&sed_trans[0],HydrologydcSedimentTransmitivityEnum); */
    4547 /*                      GetInputListOnVertices(&epl_head[0],EplHeadEnum); */
    4548 /*                      GetInputListOnVertices(&epl_head_old[0],EplHeadOldEnum); */
    4549 /*                      GetInputListOnVertices(&epl_thickness[0],HydrologydcEplThicknessEnum);                   */
    4550 
    4551 /*                      this->parameters->FindParam(&leakage,HydrologydcLeakageFactorEnum); */
    4552 
    4553 /*                      sed_thick = matpar->GetSedimentThickness(); */
    4554 
    4555 /*                      if(!active_element){ */
    4556 
    4557                                                
    4558 /*                              /\*No transfer if the EPL is not active*\/ */
    4559 /*                      } */
    4560 /*                      else{ */
    4561 /*                              GetInputListOnVertices(&preceding_transfer[0],WaterTransferEnum); */
    4562 /*                              sedstoring=matpar->GetSedimentStoring(); */
    4563 /*                              epl_specificstoring=matpar->GetEplSpecificStoring(); */
    4564                                        
    4565 /*                              for(int i=0;i<numdof;i++){ */
    4566 /*                                      this->GetHydrologyDCInefficientHmax(&h_max,this->nodes[i]); */
    4567 /*                                      /\*EPL head higher than sediment head, transfer from the epl to the sediment*\/ */
    4568 /*                                      if(epl_head[i]>sed_head[i]){ */
    4569 /*                                              if(sed_head[i]>=h_max){ */
    4570 /*                                                      wh_trans=0.0; */
    4571 /*                                              } */
    4572 /*                                              else{ */
    4573 /*                                                      wh_trans=epl_specificstoring*epl_thickness[i]*sed_trans[i]*(epl_head[i]-sed_head[i])/(leakage*sed_thick); */
    4574 /*                                              }                                                */
    4575 /*                                      } */
    4576 /*                                      /\* EPL head lower than sediment head, transfer from the sediment to the epl *\/ */
    4577 /*                                      else if(epl_head[i]<=sed_head[i]){ */
    4578 /*                                              wh_trans=sedstoring*sed_trans[i]*(epl_head[i]-sed_head[i])/(leakage*sed_thick); */
    4579 /*                                      }  */
    4580                                        
    4581 /*                                      /\*Assign output pointer*\/ */
    4582 /*                                      transfer->SetValue(doflist[i],wh_trans,INS_VAL); */
    4583 /*                              } */
    4584 /*                      } */
    4585 /*                      break; */
    4586 /*              default: */
    4587 /*                      _error_("no case higher than 1 for the Transfer method"); */
    4588 /*              } */
    4589 /*      } */
    4590 /*      /\*Free ressources:*\/ */
    4591 /*      xDelete<int>(doflist); */
    4592  }
    4593 /*}}}*/
    45944501/*FUNCTION Tria::HydrologyEPLGetActive {{{*/
    45954502void Tria::HydrologyEPLGetActive(Vector<IssmDouble>* active_vec){

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

    r17349 r17353  
    254254                ElementMatrix* CreateEPLDomainMassMatrix(void);
    255255                void           CreateHydrologyWaterVelocityInput(void);
    256                 //              void           GetHydrologyDCInefficientHmax(IssmDouble* ph_max,int index);
    257256                void           GetHydrologyDCInefficientHmax(IssmDouble* ph_max, Node* innode);
    258                 void           GetHydrologyTransfer(Vector<IssmDouble>* transfer);
    259257                void           HydrologyEPLGetActive(Vector<IssmDouble>* active_vec);
    260258                void           HydrologyEPLGetMask(Vector<IssmDouble>* vec_mask);

  • issm/trunk-jpl/src/c/classes/FemModel.cpp

    r17294 r17353  
    14071407}
    14081408/*}}}*/
    1409 void FemModel::HydrologyTransferx(void){ /*{{{*/
    1410 
    1411         Vector<IssmDouble>* transferg=NULL;
    1412 
    1413         /*Vector allocation*/
    1414         transferg=new Vector<IssmDouble>(nodes->NumberOfNodes());
    1415 
    1416         for (int i=0;i<elements->Size();i++){
    1417 
    1418                 Element* element=dynamic_cast<Element*>(elements->GetObjectByOffset(i));
    1419                 element->GetHydrologyTransfer(transferg);
    1420         }
    1421         /*Assemble*/
    1422         transferg->Assemble();
    1423 
    1424         /*Update Inputs*/
    1425         InputUpdateFromVectorx(this,transferg,WaterTransferEnum,NodesEnum);
    1426         delete transferg;
    1427 }
    1428 /*}}}*/
    14291409void FemModel::HydrologyEPLThicknessx(void){ /*{{{*/
    14301410

  • issm/trunk-jpl/src/c/classes/FemModel.h

    r17236 r17353  
    9595                void UpdateConstraintsExtrudeFromBasex();
    9696                void UpdateConstraintsExtrudeFromTopx();
    97                 void HydrologyTransferx(void);
    9897                void HydrologyEPLupdateDomainx(void);
    9998                void HydrologyEPLThicknessx(void);

  • issm/trunk-jpl/src/c/cores/hydrology_core.cpp

    r17012 r17353  
    9393                                if(VerboseSolution()) _printf0_("   saving results \n");
    9494                                if(isefficientlayer){
    95                                         int outputs[9] = {SedimentHeadEnum,SedimentHeadResidualEnum,EplHeadEnum,HydrologydcMaskEplactiveNodeEnum,HydrologydcMaskEplactiveEltEnum,EplHeadSlopeXEnum,EplHeadSlopeYEnum,HydrologydcEplThicknessEnum,WaterTransferEnum};
    96                                         femmodel->RequestedOutputsx(&femmodel->results,&outputs[0],9);
     95                                        int outputs[8] = {SedimentHeadEnum,SedimentHeadResidualEnum,EplHeadEnum,HydrologydcMaskEplactiveNodeEnum,HydrologydcMaskEplactiveEltEnum,EplHeadSlopeXEnum,EplHeadSlopeYEnum,HydrologydcEplThicknessEnum};
     96                                        femmodel->RequestedOutputsx(&femmodel->results,&outputs[0],8);
    9797                                }
    9898                                else{

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

    r17334 r17353  
    130130        HydrologyEfficientEnum,
    131131        HydrologySedimentKmaxEnum,
    132         WaterTransferEnum,
    133132        IndependentObjectEnum,
    134133        InversionControlParametersEnum,

  • issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp

    r17334 r17353  
    138138                case HydrologyEfficientEnum : return "HydrologyEfficient";
    139139                case HydrologySedimentKmaxEnum : return "HydrologySedimentKmax";
    140                 case WaterTransferEnum : return "WaterTransfer";
    141140                case IndependentObjectEnum : return "IndependentObject";
    142141                case InversionControlParametersEnum : return "InversionControlParameters";

  • issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp

    r17334 r17353  
    141141   if(stage==2){
    142142              if (strcmp(name,"HydrologySedimentKmax")==0) return HydrologySedimentKmaxEnum;
    143               else if (strcmp(name,"WaterTransfer")==0) return WaterTransferEnum;
    144143              else if (strcmp(name,"IndependentObject")==0) return IndependentObjectEnum;
    145144              else if (strcmp(name,"InversionControlParameters")==0) return InversionControlParametersEnum;
    … …  
    260259              else if (strcmp(name,"SteadystateRequestedOutputs")==0) return SteadystateRequestedOutputsEnum;
    261260              else if (strcmp(name,"Surface")==0) return SurfaceEnum;
     261              else if (strcmp(name,"ThermalIsenthalpy")==0) return ThermalIsenthalpyEnum;
    262262         else stage=3;
    263263   }
    264264   if(stage==3){
    265               if (strcmp(name,"ThermalIsenthalpy")==0) return ThermalIsenthalpyEnum;
    266               else if (strcmp(name,"ThermalIsdynamicbasalspc")==0) return ThermalIsdynamicbasalspcEnum;
     265              if (strcmp(name,"ThermalIsdynamicbasalspc")==0) return ThermalIsdynamicbasalspcEnum;
    267266              else if (strcmp(name,"ThermalMaxiter")==0) return ThermalMaxiterEnum;
    268267              else if (strcmp(name,"ThermalPenaltyFactor")==0) return ThermalPenaltyFactorEnum;
    … …  
    383382              else if (strcmp(name,"Results")==0) return ResultsEnum;
    384383              else if (strcmp(name,"GenericParam")==0) return GenericParamEnum;
     384              else if (strcmp(name,"AdolcParam")==0) return AdolcParamEnum;
    385385         else stage=4;
    386386   }
    387387   if(stage==4){
    388               if (strcmp(name,"AdolcParam")==0) return AdolcParamEnum;
    389               else if (strcmp(name,"BoolInput")==0) return BoolInputEnum;
     388              if (strcmp(name,"BoolInput")==0) return BoolInputEnum;
    390389              else if (strcmp(name,"BoolParam")==0) return BoolParamEnum;
    391390              else if (strcmp(name,"Contour")==0) return ContourEnum;
    … …  
    506505              else if (strcmp(name,"VzMesh")==0) return VzMeshEnum;
    507506              else if (strcmp(name,"Enthalpy")==0) return EnthalpyEnum;
     507              else if (strcmp(name,"EnthalpyPicard")==0) return EnthalpyPicardEnum;
    508508         else stage=5;
    509509   }
    510510   if(stage==5){
    511               if (strcmp(name,"EnthalpyPicard")==0) return EnthalpyPicardEnum;
    512               else if (strcmp(name,"ThicknessAbsGradient")==0) return ThicknessAbsGradientEnum;
     511              if (strcmp(name,"ThicknessAbsGradient")==0) return ThicknessAbsGradientEnum;
    513512              else if (strcmp(name,"ThicknessAlongGradient")==0) return ThicknessAlongGradientEnum;
    514513              else if (strcmp(name,"ThicknessAcrossGradient")==0) return ThicknessAcrossGradientEnum;
    … …  
    629628              else if (strcmp(name,"Sset")==0) return SsetEnum;
    630629              else if (strcmp(name,"Verbose")==0) return VerboseEnum;
     630              else if (strcmp(name,"TriangleInterp")==0) return TriangleInterpEnum;
    631631         else stage=6;
    632632   }
    633633   if(stage==6){
    634               if (strcmp(name,"TriangleInterp")==0) return TriangleInterpEnum;
    635               else if (strcmp(name,"BilinearInterp")==0) return BilinearInterpEnum;
     634              if (strcmp(name,"BilinearInterp")==0) return BilinearInterpEnum;
    636635              else if (strcmp(name,"NearestInterp")==0) return NearestInterpEnum;
    637636              else if (strcmp(name,"XY")==0) return XYEnum;

  • issm/trunk-jpl/src/c/solutionsequences/solutionsequence_hydro_nonlinear.cpp

    r17351 r17353  
    7171                delete analysis;
    7272        }
    73         /*For the initialization we compute the transfer without the mask if the EPL is not present*/
    74         femmodel->HydrologyTransferx();
    75 
    7673        /*The real computation starts here, outermost loop is on the two layer system*/
    7774        for(;;){
Note: See TracChangeset for help on using the changeset viewer.