Index: /issm/trunk-jpl/src/c/analyses/hydrology_core.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/hydrology_core.cpp (revision 15014) +++ /issm/trunk-jpl/src/c/analyses/hydrology_core.cpp (revision 15015) @@ -79,7 +79,10 @@ else if (hydrology_model==HydrologydcEnum){ - InputDuplicatex(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,SedimentHeadEnum,SedimentHeadOldEnum); + InputDuplicatex(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,SedimentHeadEnum,SedimentHeadOldEnum); femmodel->parameters->FindParam(&isefficientlayer,HydrologydcIsefficientlayerEnum); + + if(VerboseSolution()) _pprintLine_(" computing water transfer"); + if(VerboseSolution()) _pprintLine_(" computing water head"); solutionsequence_hydro_nonlinear(femmodel); Index: /issm/trunk-jpl/src/c/classes/Elements/Element.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Element.h (revision 15014) +++ /issm/trunk-jpl/src/c/classes/Elements/Element.h (revision 15015) @@ -136,4 +136,5 @@ virtual void GetHydrologyDCInefficientHmax(IssmDouble* ph_max, Node* innode)=0; virtual void BasisIntegral(Vector* basisg)=0; + virtual void GetHydrologyTransfer(Vector* transfer)=0; #endif }; Index: /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp (revision 15014) +++ /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp (revision 15015) @@ -9344,5 +9344,5 @@ } /*}}}*/ -/*FUNCTION Tria::BasisIntegral {{{*/ +/*FUNCTION Penta::BasisIntegral {{{*/ void Penta::BasisIntegral(Vector* basisg){ _error_("Hydrological stuff not suported in Penta"); @@ -9350,3 +9350,11 @@ /*}}}*/ +/*}}}*/ +/*FUNCTION Penta::GetHydrologyTransfer{{{*/ +void Penta::GetHydrologyTransfer(Vector* transfer){ + _error_("Hydrological stuff not suported in Penta"); +} +/*}}}*/ + + #endif Index: /issm/trunk-jpl/src/c/classes/Elements/Penta.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Penta.h (revision 15014) +++ /issm/trunk-jpl/src/c/classes/Elements/Penta.h (revision 15015) @@ -307,4 +307,5 @@ void GetHydrologyDCInefficientHmax(IssmDouble* ph_max, Node* innode); void BasisIntegral(Vector* gbasis); + void GetHydrologyTransfer(Vector* transfer); #endif #ifdef _HAVE_THERMAL_ Index: /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp (revision 15014) +++ /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp (revision 15015) @@ -6466,4 +6466,63 @@ } /*}}}*/ +/*FUNCTION Tria::GetHydrologyTransfer{{{*/ +void Tria::GetHydrologyTransfer(Vector* transfer){ + + const int numdof = NDOF1 *NUMVERTICES; + int *doflist = NULL; + bool isefficientlayer; + int transfermethod; + IssmDouble sed_trans,sed_thick; + IssmDouble leakage; + IssmDouble wh_trans[numdof]={0.0}; + IssmDouble storing[numdof]; + IssmDouble epl_head[numdof],sed_head[numdof]; + + + GetDofList(&doflist,NoneApproximationEnum,GsetEnum); + + /*Get the flag to know if the efficient layer is present*/ + this->parameters->FindParam(&isefficientlayer,HydrologydcIsefficientlayerEnum); + + if(!isefficientlayer){ + transfer->SetValues(numdof,doflist,&wh_trans[0],ADD_VAL); + return; + } + /*Also get the flag to the transfer method (TO DO)*/ + this->parameters->FindParam(&transfermethod,HydrologydcTransferFlagEnum); + + /*Switch between the different transfer methods cases(TO DO)*/ + switch(transfermethod){ + case 0: + /*Just kipping the transfer to zero*/ + break; + case 1: + this->parameters->FindParam(&leakage,HydrologydcLeakageFactorEnum); + + sed_trans = matpar->GetSedimentTransmitivity(); + sed_thick = matpar->GetSedimentThickness(); + + GetInputListOnVertices(&sed_head[0],SedimentHeadEnum); + GetInputListOnVertices(&epl_head[0],EplHeadEnum); + + for(int i=0;iepl_head[i]){ + storing[i]=matpar->GetSedimentStoring(); + wh_trans[i]=sed_trans*storing[i]*(sed_head[i]-epl_head[i])/(leakage*sed_thick); + } + else{ + storing[i]=matpar->GetEplStoring(); + wh_trans[i]=sed_trans*storing[i]*(sed_head[i]-epl_head[i])/(leakage*sed_thick); + } + } + default: + _error_("no case higher than 1 for the Transfer method"); + } + + /*Assign output pointer*/ + transfer->SetValues(numdof,doflist,&wh_trans[0],ADD_VAL); +} + +/*}}}*/ #endif Index: /issm/trunk-jpl/src/c/classes/Elements/Tria.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Tria.h (revision 15014) +++ /issm/trunk-jpl/src/c/classes/Elements/Tria.h (revision 15015) @@ -254,4 +254,5 @@ void GetHydrologyDCInefficientHmax(IssmDouble* ph_max, Node* innode); void BasisIntegral(Vector* gbasis); + void GetHydrologyTransfer(Vector* transfer); #endif #ifdef _HAVE_BALANCED_ Index: /issm/trunk-jpl/src/c/classes/FemModel.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/FemModel.cpp (revision 15014) +++ /issm/trunk-jpl/src/c/classes/FemModel.cpp (revision 15015) @@ -1594,2 +1594,22 @@ } /*}}}*/ + +void FemModel::HydrologyTransferx(void){ /*{{{*/ + + Vector* transferg=NULL; + + /*Vector allocation*/ + transferg=new Vector(nodes->NumberOfNodes()); + + for (int i=0;iSize();i++){ + Element* element=dynamic_cast(elements->GetObjectByOffset(i)); + element->GetHydrologyTransfer(transferg); + } + /*Assemble*/ + transferg->Assemble(); + + /*Update Inputs*/ + InputUpdateFromVectorx(elements,nodes,vertices,loads,materials,parameters,transferg,WaterTransferEnum,NodesEnum); + +} +/*}}}*/ Index: /issm/trunk-jpl/src/c/classes/FemModel.h =================================================================== --- /issm/trunk-jpl/src/c/classes/FemModel.h (revision 15014) +++ /issm/trunk-jpl/src/c/classes/FemModel.h (revision 15015) @@ -100,5 +100,6 @@ void UpdateConstraintsx(void); int UpdateVertexPositionsx(void); - void BasisIntegralsx(void); + void BasisIntegralsx(void); + void HydrologyTransferx(void); }; Index: /issm/trunk-jpl/src/c/classes/Materials/Matpar.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Materials/Matpar.cpp (revision 15014) +++ /issm/trunk-jpl/src/c/classes/Materials/Matpar.cpp (revision 15015) @@ -370,4 +370,9 @@ } /*}}}*/ +/*FUNCTION Matpar::GetSedimentThickness {{{*/ +IssmDouble Matpar::GetSedimentThickness(){ + return sediment_thickness; +} +/*}}}*/ /*FUNCTION Matpar::GetEplTransitivity {{{*/ IssmDouble Matpar::GetEplTransmitivity(){ Index: /issm/trunk-jpl/src/c/classes/Materials/Matpar.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Materials/Matpar.h (revision 15014) +++ /issm/trunk-jpl/src/c/classes/Materials/Matpar.h (revision 15015) @@ -123,4 +123,5 @@ IssmDouble GetEplStoring(); IssmDouble GetSedimentTransmitivity(); + IssmDouble GetSedimentThickness(); IssmDouble GetEplTransmitivity(); IssmDouble TMeltingPoint(IssmDouble pressure); Index: /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h =================================================================== --- /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h (revision 15014) +++ /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h (revision 15015) @@ -104,4 +104,6 @@ HydrologydcSedimentlimitFlagEnum, HydrologydcSedimentlimitEnum, + HydrologydcTransferFlagEnum, + HydrologydcLeakageFactorEnum, HydrologydcPenaltyFactorEnum, HydrologyLayerEnum, @@ -110,4 +112,5 @@ HydrologySedimentKmaxEnum, BasisIntegralEnum, + WaterTransferEnum, IndependentObjectEnum, InversionControlParametersEnum, Index: /issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp =================================================================== --- /issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp (revision 15014) +++ /issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp (revision 15015) @@ -112,4 +112,6 @@ case HydrologydcSedimentlimitFlagEnum : return "HydrologydcSedimentlimitFlag"; case HydrologydcSedimentlimitEnum : return "HydrologydcSedimentlimit"; + case HydrologydcTransferFlagEnum : return "HydrologydcTransferFlag"; + case HydrologydcLeakageFactorEnum : return "HydrologydcLeakageFactor"; case HydrologydcPenaltyFactorEnum : return "HydrologydcPenaltyFactor"; case HydrologyLayerEnum : return "HydrologyLayer"; @@ -118,4 +120,5 @@ case HydrologySedimentKmaxEnum : return "HydrologySedimentKmax"; case BasisIntegralEnum : return "BasisIntegral"; + case WaterTransferEnum : return "WaterTransfer"; case IndependentObjectEnum : return "IndependentObject"; case InversionControlParametersEnum : return "InversionControlParameters"; Index: /issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp =================================================================== --- /issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp (revision 15014) +++ /issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp (revision 15015) @@ -112,4 +112,6 @@ else if (strcmp(name,"HydrologydcSedimentlimitFlag")==0) return HydrologydcSedimentlimitFlagEnum; else if (strcmp(name,"HydrologydcSedimentlimit")==0) return HydrologydcSedimentlimitEnum; + else if (strcmp(name,"HydrologydcTransferFlag")==0) return HydrologydcTransferFlagEnum; + else if (strcmp(name,"HydrologydcLeakageFactor")==0) return HydrologydcLeakageFactorEnum; else if (strcmp(name,"HydrologydcPenaltyFactor")==0) return HydrologydcPenaltyFactorEnum; else if (strcmp(name,"HydrologyLayer")==0) return HydrologyLayerEnum; @@ -118,4 +120,5 @@ else if (strcmp(name,"HydrologySedimentKmax")==0) return HydrologySedimentKmaxEnum; else if (strcmp(name,"BasisIntegral")==0) return BasisIntegralEnum; + else if (strcmp(name,"WaterTransfer")==0) return WaterTransferEnum; else if (strcmp(name,"IndependentObject")==0) return IndependentObjectEnum; else if (strcmp(name,"InversionControlParameters")==0) return InversionControlParametersEnum; @@ -134,11 +137,11 @@ else if (strcmp(name,"InversionNumControlParameters")==0) return InversionNumControlParametersEnum; else if (strcmp(name,"InversionNumCostFunctions")==0) return InversionNumCostFunctionsEnum; - else if (strcmp(name,"InversionStepThreshold")==0) return InversionStepThresholdEnum; - else if (strcmp(name,"InversionThicknessObs")==0) return InversionThicknessObsEnum; - else if (strcmp(name,"InversionVelObs")==0) return InversionVelObsEnum; else stage=2; } if(stage==2){ - if (strcmp(name,"InversionVxObs")==0) return InversionVxObsEnum; + if (strcmp(name,"InversionStepThreshold")==0) return InversionStepThresholdEnum; + else if (strcmp(name,"InversionThicknessObs")==0) return InversionThicknessObsEnum; + else if (strcmp(name,"InversionVelObs")==0) return InversionVelObsEnum; + else if (strcmp(name,"InversionVxObs")==0) return InversionVxObsEnum; else if (strcmp(name,"InversionVyObs")==0) return InversionVyObsEnum; else if (strcmp(name,"InversionVzObs")==0) return InversionVzObsEnum; @@ -257,11 +260,11 @@ else if (strcmp(name,"TimesteppingCflCoefficient")==0) return TimesteppingCflCoefficientEnum; else if (strcmp(name,"TimesteppingTimeAdapt")==0) return TimesteppingTimeAdaptEnum; - else if (strcmp(name,"TimesteppingTimeStep")==0) return TimesteppingTimeStepEnum; - else if (strcmp(name,"TransientIsdiagnostic")==0) return TransientIsdiagnosticEnum; - else if (strcmp(name,"TransientIsgroundingline")==0) return TransientIsgroundinglineEnum; else stage=3; } if(stage==3){ - if (strcmp(name,"TransientIsprognostic")==0) return TransientIsprognosticEnum; + if (strcmp(name,"TimesteppingTimeStep")==0) return TimesteppingTimeStepEnum; + else if (strcmp(name,"TransientIsdiagnostic")==0) return TransientIsdiagnosticEnum; + else if (strcmp(name,"TransientIsgroundingline")==0) return TransientIsgroundinglineEnum; + else if (strcmp(name,"TransientIsprognostic")==0) return TransientIsprognosticEnum; else if (strcmp(name,"TransientIsthermal")==0) return TransientIsthermalEnum; else if (strcmp(name,"TransientIsgia")==0) return TransientIsgiaEnum; @@ -380,11 +383,11 @@ else if (strcmp(name,"StokesIceFront")==0) return StokesIceFrontEnum; else if (strcmp(name,"StringArrayParam")==0) return StringArrayParamEnum; - else if (strcmp(name,"StringParam")==0) return StringParamEnum; - else if (strcmp(name,"Tria")==0) return TriaEnum; - else if (strcmp(name,"TriaP1Input")==0) return TriaP1InputEnum; else stage=4; } if(stage==4){ - if (strcmp(name,"Vertex")==0) return VertexEnum; + if (strcmp(name,"StringParam")==0) return StringParamEnum; + else if (strcmp(name,"Tria")==0) return TriaEnum; + else if (strcmp(name,"TriaP1Input")==0) return TriaP1InputEnum; + else if (strcmp(name,"Vertex")==0) return VertexEnum; else if (strcmp(name,"VertexPId")==0) return VertexPIdEnum; else if (strcmp(name,"VertexSId")==0) return VertexSIdEnum; @@ -503,11 +506,11 @@ else if (strcmp(name,"PentaP1ElementResult")==0) return PentaP1ElementResultEnum; else if (strcmp(name,"StringExternalResult")==0) return StringExternalResultEnum; - else if (strcmp(name,"Step")==0) return StepEnum; - else if (strcmp(name,"Time")==0) return TimeEnum; - else if (strcmp(name,"TriaP1ElementResult")==0) return TriaP1ElementResultEnum; else stage=5; } if(stage==5){ - if (strcmp(name,"WaterColumnOld")==0) return WaterColumnOldEnum; + if (strcmp(name,"Step")==0) return StepEnum; + else if (strcmp(name,"Time")==0) return TimeEnum; + else if (strcmp(name,"TriaP1ElementResult")==0) return TriaP1ElementResultEnum; + else if (strcmp(name,"WaterColumnOld")==0) return WaterColumnOldEnum; else if (strcmp(name,"MinVel")==0) return MinVelEnum; else if (strcmp(name,"MaxVel")==0) return MaxVelEnum; Index: /issm/trunk-jpl/src/m/classes/hydrologydc.m =================================================================== --- /issm/trunk-jpl/src/m/classes/hydrologydc.m (revision 15014) +++ /issm/trunk-jpl/src/m/classes/hydrologydc.m (revision 15015) @@ -16,4 +16,6 @@ sedimentlimit_flag = 0; sedimentlimit = 0; + transfer_flag = 0; + leakage_factor = 0; spcepl_head = NaN; @@ -45,4 +47,6 @@ obj.sedimentlimit_flag = 0; obj.sedimentlimit = 0; + obj.transfer_flag = 0; + obj.leakage_factor = 10.0; obj.epl_compressibility = 1.0e-08; @@ -68,9 +72,14 @@ md = checkfield(md,'hydrology.penalty_factor','>',0,'numel',1); md = checkfield(md,'hydrology.sedimentlimit_flag','numel',[1],'values',[0 1 2 3]); - + md = checkfield(md,'hydrology.transfer_flag','numel',[1],'values',[0 1]); + if obj.sedimentlimit_flag==1, md = checkfield(md,'hydrology.sedimentlimit','>',0,'numel',1); end - + + if obj.transfer_flag==1, + md = checkfield(md,'hydrology.leakage_factor','>',0,'numel',1); + end + if obj.isefficientlayer==1, @@ -102,4 +111,10 @@ fielddisplay(obj,'sedimentlimit','user defined upper limit for the inefficient layer [m]'); end + fielddisplay(obj,'transfer_flag',['what kind of transfer method is applied between the layers']); + disp(sprintf('%55s 0: no transfer',' ')); + disp(sprintf('%55s 1: constant leakage factor: %s',' ','leakage_factor')); + if obj.transfer_flag==1, + fielddisplay(obj,'leakage_factor','user defined leakage factor [m]'); + end if obj.isefficientlayer==1, @@ -125,7 +140,12 @@ WriteData(fid,'object',obj,'fieldname','penalty_factor','format','Double'); WriteData(fid,'object',obj,'fieldname','sedimentlimit_flag','format','Integer'); + WriteData(fid,'object',obj,'fieldname','transfer_flag','format','Integer'); if obj.sedimentlimit_flag==1, WriteData(fid,'object',obj,'fieldname','sedimentlimit','format','Double'); + end + + if obj.transfer_flag==1, + WriteData(fid,'object',obj,'fieldname','leakage_factor','format','Double'); end Index: /issm/trunk-jpl/src/m/enum/EnumDefinitions.py =================================================================== --- /issm/trunk-jpl/src/m/enum/EnumDefinitions.py (revision 15014) +++ /issm/trunk-jpl/src/m/enum/EnumDefinitions.py (revision 15015) @@ -1339,4 +1339,32 @@ return StringToEnum('HydrologydcSedimentlimit')[0] +def HydrologydcTransferFlagEnum(): + """ + HYDROLOGYDCTRANSFERFLAGENUM - Enum of HydrologydcTransferFlag + + WARNING: DO NOT MODIFY THIS FILE + this file has been automatically generated by src/c/shared/Enum/Synchronize.sh + Please read src/c/shared/Enum/README for more information + + Usage: + macro=HydrologydcTransferFlagEnum() + """ + + return StringToEnum('HydrologydcTransferFlag')[0] + +def HydrologydcLeakageFactorEnum(): + """ + HYDROLOGYDCLEAKAGEFACTORENUM - Enum of HydrologydcLeakageFactor + + WARNING: DO NOT MODIFY THIS FILE + this file has been automatically generated by src/c/shared/Enum/Synchronize.sh + Please read src/c/shared/Enum/README for more information + + Usage: + macro=HydrologydcLeakageFactorEnum() + """ + + return StringToEnum('HydrologydcLeakageFactor')[0] + def HydrologydcPenaltyFactorEnum(): """ @@ -1423,4 +1451,18 @@ return StringToEnum('BasisIntegral')[0] +def WaterTransferEnum(): + """ + WATERTRANSFERENUM - Enum of WaterTransfer + + WARNING: DO NOT MODIFY THIS FILE + this file has been automatically generated by src/c/shared/Enum/Synchronize.sh + Please read src/c/shared/Enum/README for more information + + Usage: + macro=WaterTransferEnum() + """ + + return StringToEnum('WaterTransfer')[0] + def IndependentObjectEnum(): """ @@ -7721,4 +7763,4 @@ """ - return 550 - + return 553 + Index: /issm/trunk-jpl/src/m/enum/HydrologydcLeakageFactorEnum.m =================================================================== --- /issm/trunk-jpl/src/m/enum/HydrologydcLeakageFactorEnum.m (revision 15015) +++ /issm/trunk-jpl/src/m/enum/HydrologydcLeakageFactorEnum.m (revision 15015) @@ -0,0 +1,11 @@ +function macro=HydrologydcLeakageFactorEnum() +%HYDROLOGYDCLEAKAGEFACTORENUM - Enum of HydrologydcLeakageFactor +% +% WARNING: DO NOT MODIFY THIS FILE +% this file has been automatically generated by src/c/shared/Enum/Synchronize.sh +% Please read src/c/shared/Enum/README for more information +% +% Usage: +% macro=HydrologydcLeakageFactorEnum() + +macro=StringToEnum('HydrologydcLeakageFactor'); Index: /issm/trunk-jpl/src/m/enum/HydrologydcTransferFlagEnum.m =================================================================== --- /issm/trunk-jpl/src/m/enum/HydrologydcTransferFlagEnum.m (revision 15015) +++ /issm/trunk-jpl/src/m/enum/HydrologydcTransferFlagEnum.m (revision 15015) @@ -0,0 +1,11 @@ +function macro=HydrologydcTransferFlagEnum() +%HYDROLOGYDCTRANSFERFLAGENUM - Enum of HydrologydcTransferFlag +% +% WARNING: DO NOT MODIFY THIS FILE +% this file has been automatically generated by src/c/shared/Enum/Synchronize.sh +% Please read src/c/shared/Enum/README for more information +% +% Usage: +% macro=HydrologydcTransferFlagEnum() + +macro=StringToEnum('HydrologydcTransferFlag'); Index: /issm/trunk-jpl/src/m/enum/MaximumNumberOfEnums.m =================================================================== --- /issm/trunk-jpl/src/m/enum/MaximumNumberOfEnums.m (revision 15014) +++ /issm/trunk-jpl/src/m/enum/MaximumNumberOfEnums.m (revision 15015) @@ -9,3 +9,3 @@ % macro=MaximumNumberOfEnums() -macro=550; +macro=553; Index: /issm/trunk-jpl/src/m/enum/WaterTransferEnum.m =================================================================== --- /issm/trunk-jpl/src/m/enum/WaterTransferEnum.m (revision 15015) +++ /issm/trunk-jpl/src/m/enum/WaterTransferEnum.m (revision 15015) @@ -0,0 +1,11 @@ +function macro=WaterTransferEnum() +%WATERTRANSFERENUM - Enum of WaterTransfer +% +% WARNING: DO NOT MODIFY THIS FILE +% this file has been automatically generated by src/c/shared/Enum/Synchronize.sh +% Please read src/c/shared/Enum/README for more information +% +% Usage: +% macro=WaterTransferEnum() + +macro=StringToEnum('WaterTransfer');