Index: ../trunk-jpl/src/c/analyses/LevelsetAnalysis.cpp =================================================================== --- ../trunk-jpl/src/c/analyses/LevelsetAnalysis.cpp (revision 18757) +++ ../trunk-jpl/src/c/analyses/LevelsetAnalysis.cpp (revision 18758) @@ -19,10 +19,9 @@ } /*}}}*/ void LevelsetAnalysis::UpdateElements(Elements* elements,IoModel* iomodel,int analysis_counter,int analysis_type){/*{{{*/ - int finiteelement; /*Finite element type*/ - finiteelement = P1Enum; + int finiteelement = P1Enum; /*Update elements: */ int counter=0; @@ -37,7 +36,20 @@ iomodel->FetchDataToInput(elements,MaskIceLevelsetEnum); iomodel->FetchDataToInput(elements,VxEnum); iomodel->FetchDataToInput(elements,VyEnum); - iomodel->FetchDataToInput(elements,MasstransportCalvingrateEnum); + + /*Get calving parameters*/ + int calvinglaw; + iomodel->Constant(&calvinglaw,CalvingLawEnum); + switch(calvinglaw){ + case DefaultCalvingEnum: + iomodel->FetchDataToInput(elements,CalvingCalvingrateEnum); + break; + case CalvingLevermannEnum: + iomodel->FetchDataToInput(elements,CalvinglevermannCoeffEnum); + break; + default: + _error_("Calving law "<SpawnBasalElement(); /*Intermediaries */ - int dim, domaintype; - bool iscalvingrate; + int stabilization=2; + int dim, domaintype, calvinglaw; + bool iscalving; int i, row, col; IssmDouble kappa; IssmDouble Jdet, dt, D_scalar; IssmDouble h,hx,hy,hz; IssmDouble vel; -// IssmDouble norm_dlsf; + IssmDouble norm_dlsf, calvingrate; IssmDouble* xyz_list = NULL; /*Get problem dimension and whether there is calving or not*/ - basalelement->FindParam(&iscalvingrate,MasstransportIscalvingrateEnum); + basalelement->FindParam(&iscalving,TransientIscalvingEnum); basalelement->FindParam(&domaintype,DomainTypeEnum); + basalelement->FindParam(&calvinglaw,CalvingLawEnum); switch(domaintype){ case Domain2DverticalEnum: dim = 1; break; case Domain2DhorizontalEnum: dim = 2; break; @@ -130,45 +144,63 @@ IssmDouble* D = xNew(dim*dim); IssmDouble* v = xNew(dim); IssmDouble* w = xNew(dim); - IssmDouble* c = xNew(dim); - //IssmDouble* dlsf = xNew(dim); + IssmDouble* c = xNewZeroInit(dim); + IssmDouble* dlsf = xNew(dim); /*Retrieve all inputs and parameters*/ basalelement->GetVerticesCoordinates(&xyz_list); basalelement->FindParam(&dt,TimesteppingTimeStepEnum); - Input* vx_input=NULL; - Input* vy_input=NULL; - Input* calvingratex_input=NULL; - Input* calvingratey_input=NULL; + Input* vx_input = NULL; + Input* vy_input = NULL; + Input* calvingratex_input = NULL; + Input* calvingratey_input = NULL; + Input* lsf_slopex_input = NULL; + Input* lsf_slopey_input = NULL; + Input* calvingrate_input = NULL; + + /*Load velocities*/ if(domaintype==Domain2DhorizontalEnum){ vx_input=basalelement->GetInput(VxEnum); _assert_(vx_input); vy_input=basalelement->GetInput(VyEnum); _assert_(vy_input); - if(iscalvingrate){ - calvingratex_input=basalelement->GetInput(CalvingratexEnum); _assert_(calvingratex_input); - calvingratey_input=basalelement->GetInput(CalvingrateyEnum); _assert_(calvingratey_input); - } } else{ if(dim==1){ vx_input=basalelement->GetInput(VxEnum); _assert_(vx_input); - if(iscalvingrate){ - calvingratex_input=basalelement->GetInput(CalvingratexEnum); _assert_(calvingratex_input); - } } if(dim==2){ vx_input=basalelement->GetInput(VxAverageEnum); _assert_(vx_input); vy_input=basalelement->GetInput(VyAverageEnum); _assert_(vy_input); - if(iscalvingrate){ - calvingratex_input=basalelement->GetInput(CalvingratexAverageEnum); _assert_(calvingratex_input); - calvingratey_input=basalelement->GetInput(CalvingrateyAverageEnum); _assert_(calvingratey_input); - } } } -// Input* lsf_slopex_input = basalelement->GetInput(LevelsetfunctionSlopeXEnum); _assert_(lsf_slopex_input); -// Input* lsf_slopey_input = basalelement->GetInput(LevelsetfunctionSlopeYEnum); _assert_(lsf_slopey_input); - //Input* calvingrate_input = basalelement->GetInput(MasstransportCalvingrateEnum); _assert_(calvingrate_input); - + /*Load calving inputs*/ + if(iscalving){ + switch(calvinglaw){ + case DefaultCalvingEnum: + lsf_slopex_input = basalelement->GetInput(LevelsetfunctionSlopeXEnum); _assert_(lsf_slopex_input); + if(dim==2) lsf_slopey_input = basalelement->GetInput(LevelsetfunctionSlopeYEnum); _assert_(lsf_slopey_input); + calvingrate_input = basalelement->GetInput(CalvingCalvingrateEnum); _assert_(calvingrate_input); + break; + case CalvingLevermannEnum: + if(domaintype==Domain2DhorizontalEnum){ + calvingratex_input=basalelement->GetInput(CalvingratexEnum); _assert_(calvingratex_input); + calvingratey_input=basalelement->GetInput(CalvingrateyEnum); _assert_(calvingratey_input); + } + else{ + if(dim==1){ + calvingratex_input=basalelement->GetInput(CalvingratexEnum); _assert_(calvingratex_input); + } + if(dim==2){ + calvingratex_input=basalelement->GetInput(CalvingratexAverageEnum); _assert_(calvingratex_input); + calvingratey_input=basalelement->GetInput(CalvingrateyAverageEnum); _assert_(calvingratey_input); + } + } + break; + default: + _error_("Calving law "<NewGauss(2); for(int ig=gauss->begin();igend();ig++){ @@ -192,34 +224,45 @@ GetBprime(Bprime,basalelement,xyz_list,gauss); vx_input->GetInputValue(&v[0],gauss); // in 3D case, add mesh velocity vy_input->GetInputValue(&v[1],gauss); - if(iscalvingrate){ - calvingratex_input->GetInputValue(&c[0],gauss); // in 3D case, add mesh velocity - calvingratey_input->GetInputValue(&c[1],gauss); - for(i=0;iGetInputValue(&dlsf[0],gauss); + if(dim==2) lsf_slopey_input->GetInputValue(&dlsf[1],gauss); + calvingrate_input->GetInputValue(&calvingrate,gauss); + + norm_dlsf=0.; + for(i=0;i1.e-10) + for(i=0;iGetInputValue(&c[0],gauss); // in 3D case, add mesh velocity + if(dim==2) calvingratey_input->GetInputValue(&c[1],gauss); + break; + default: + _error_("Calving law "<GetInputValue(&dlsf[0],gauss); - //lsf_slopey_input->GetInputValue(&dlsf[1],gauss); - //calvingrate_input->GetInputValue(&calvingrate,gauss); - //norm_dlsf=0.; - //for(i=0;i1.e-10) - // for(i=0;ivalues[0],1); /* Stabilization */ - int stabilization=2; vel=0.; - for(i=0;iElementSizes(&hx,&hy,&hz); - h=sqrt( pow(hx*w[0]/vel,2) + pow(hy*w[1]/vel,2) ); + h=sqrt( pow(hx*v[0]/vel,2) + pow(hy*v[1]/vel,2) ); for(row=0;row(v); xDelete(w); xDelete(c); - //xDelete(dlsf); + xDelete(dlsf); delete gauss; if(domaintype!=Domain2DhorizontalEnum){basalelement->DeleteMaterials(); delete basalelement;}; return Ke; Index: ../trunk-jpl/src/c/shared/Enum/EnumDefinitions.h =================================================================== --- ../trunk-jpl/src/c/shared/Enum/EnumDefinitions.h (revision 18757) +++ ../trunk-jpl/src/c/shared/Enum/EnumDefinitions.h (revision 18758) @@ -208,6 +208,12 @@ DamageEnum, NewDamageEnum, StressIntensityFactorEnum, + CalvingLawEnum, + CalvingCalvingrateEnum, + CalvingLevermannEnum, + DefaultCalvingEnum, + CalvingRequestedOutputsEnum, + CalvinglevermannCoeffEnum, CalvingratexEnum, CalvingrateyEnum, CalvingratexAverageEnum, @@ -249,13 +255,10 @@ Domain3DEnum, MiscellaneousNameEnum, //FIXME: only used by qmu, should not be marshalled (already in queueing script) MasstransportHydrostaticAdjustmentEnum, - MasstransportIscalvingrateEnum, - MasstransportLevermannCalvingCoeffEnum, MasstransportIsfreesurfaceEnum, MasstransportMinThicknessEnum, MasstransportPenaltyFactorEnum, MasstransportSpcthicknessEnum, - MasstransportCalvingrateEnum, MasstransportStabilizationEnum, MasstransportVertexPairingEnum, MasstransportNumRequestedOutputsEnum, @@ -313,6 +316,7 @@ TransientIsgiaEnum, TransientIsdamageevolutionEnum, TransientIshydrologyEnum, + TransientIscalvingEnum, TransientNumRequestedOutputsEnum, TransientRequestedOutputsEnum, PotentialEnum, Index: ../trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp =================================================================== --- ../trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp (revision 18757) +++ ../trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp (revision 18758) @@ -216,6 +216,12 @@ case DamageEnum : return "Damage"; case NewDamageEnum : return "NewDamage"; case StressIntensityFactorEnum : return "StressIntensityFactor"; + case CalvingLawEnum : return "CalvingLaw"; + case CalvingCalvingrateEnum : return "CalvingCalvingrate"; + case CalvingLevermannEnum : return "CalvingLevermann"; + case DefaultCalvingEnum : return "DefaultCalving"; + case CalvingRequestedOutputsEnum : return "CalvingRequestedOutputs"; + case CalvinglevermannCoeffEnum : return "CalvinglevermannCoeff"; case CalvingratexEnum : return "Calvingratex"; case CalvingrateyEnum : return "Calvingratey"; case CalvingratexAverageEnum : return "CalvingratexAverage"; @@ -257,13 +263,10 @@ case Domain3DEnum : return "Domain3D"; case MiscellaneousNameEnum : return "MiscellaneousName"; case MasstransportHydrostaticAdjustmentEnum : return "MasstransportHydrostaticAdjustment"; - case MasstransportIscalvingrateEnum : return "MasstransportIscalvingrate"; - case MasstransportLevermannCalvingCoeffEnum : return "MasstransportLevermannCalvingCoeff"; case MasstransportIsfreesurfaceEnum : return "MasstransportIsfreesurface"; case MasstransportMinThicknessEnum : return "MasstransportMinThickness"; case MasstransportPenaltyFactorEnum : return "MasstransportPenaltyFactor"; case MasstransportSpcthicknessEnum : return "MasstransportSpcthickness"; - case MasstransportCalvingrateEnum : return "MasstransportCalvingrate"; case MasstransportStabilizationEnum : return "MasstransportStabilization"; case MasstransportVertexPairingEnum : return "MasstransportVertexPairing"; case MasstransportNumRequestedOutputsEnum : return "MasstransportNumRequestedOutputs"; @@ -321,6 +324,7 @@ case TransientIsgiaEnum : return "TransientIsgia"; case TransientIsdamageevolutionEnum : return "TransientIsdamageevolution"; case TransientIshydrologyEnum : return "TransientIshydrology"; + case TransientIscalvingEnum : return "TransientIscalving"; case TransientNumRequestedOutputsEnum : return "TransientNumRequestedOutputs"; case TransientRequestedOutputsEnum : return "TransientRequestedOutputs"; case PotentialEnum : return "Potential"; Index: ../trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp =================================================================== --- ../trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp (revision 18757) +++ ../trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp (revision 18758) @@ -219,6 +219,12 @@ else if (strcmp(name,"Damage")==0) return DamageEnum; else if (strcmp(name,"NewDamage")==0) return NewDamageEnum; else if (strcmp(name,"StressIntensityFactor")==0) return StressIntensityFactorEnum; + else if (strcmp(name,"CalvingLaw")==0) return CalvingLawEnum; + else if (strcmp(name,"CalvingCalvingrate")==0) return CalvingCalvingrateEnum; + else if (strcmp(name,"CalvingLevermann")==0) return CalvingLevermannEnum; + else if (strcmp(name,"DefaultCalving")==0) return DefaultCalvingEnum; + else if (strcmp(name,"CalvingRequestedOutputs")==0) return CalvingRequestedOutputsEnum; + else if (strcmp(name,"CalvinglevermannCoeff")==0) return CalvinglevermannCoeffEnum; else if (strcmp(name,"Calvingratex")==0) return CalvingratexEnum; else if (strcmp(name,"Calvingratey")==0) return CalvingrateyEnum; else if (strcmp(name,"CalvingratexAverage")==0) return CalvingratexAverageEnum; @@ -253,23 +259,20 @@ else if (strcmp(name,"MeshY")==0) return MeshYEnum; else if (strcmp(name,"MeshZ")==0) return MeshZEnum; else if (strcmp(name,"MeshElementtype")==0) return MeshElementtypeEnum; - else if (strcmp(name,"DomainType")==0) return DomainTypeEnum; + else stage=3; + } + if(stage==3){ + if (strcmp(name,"DomainType")==0) return DomainTypeEnum; else if (strcmp(name,"DomainDimension")==0) return DomainDimensionEnum; else if (strcmp(name,"Domain2Dhorizontal")==0) return Domain2DhorizontalEnum; else if (strcmp(name,"Domain2Dvertical")==0) return Domain2DverticalEnum; else if (strcmp(name,"Domain3D")==0) return Domain3DEnum; else if (strcmp(name,"MiscellaneousName")==0) return MiscellaneousNameEnum; - else stage=3; - } - if(stage==3){ - if (strcmp(name,"MasstransportHydrostaticAdjustment")==0) return MasstransportHydrostaticAdjustmentEnum; - else if (strcmp(name,"MasstransportIscalvingrate")==0) return MasstransportIscalvingrateEnum; - else if (strcmp(name,"MasstransportLevermannCalvingCoeff")==0) return MasstransportLevermannCalvingCoeffEnum; + else if (strcmp(name,"MasstransportHydrostaticAdjustment")==0) return MasstransportHydrostaticAdjustmentEnum; else if (strcmp(name,"MasstransportIsfreesurface")==0) return MasstransportIsfreesurfaceEnum; else if (strcmp(name,"MasstransportMinThickness")==0) return MasstransportMinThicknessEnum; else if (strcmp(name,"MasstransportPenaltyFactor")==0) return MasstransportPenaltyFactorEnum; else if (strcmp(name,"MasstransportSpcthickness")==0) return MasstransportSpcthicknessEnum; - else if (strcmp(name,"MasstransportCalvingrate")==0) return MasstransportCalvingrateEnum; else if (strcmp(name,"MasstransportStabilization")==0) return MasstransportStabilizationEnum; else if (strcmp(name,"MasstransportVertexPairing")==0) return MasstransportVertexPairingEnum; else if (strcmp(name,"MasstransportNumRequestedOutputs")==0) return MasstransportNumRequestedOutputsEnum; @@ -327,6 +330,7 @@ else if (strcmp(name,"TransientIsgia")==0) return TransientIsgiaEnum; else if (strcmp(name,"TransientIsdamageevolution")==0) return TransientIsdamageevolutionEnum; else if (strcmp(name,"TransientIshydrology")==0) return TransientIshydrologyEnum; + else if (strcmp(name,"TransientIscalving")==0) return TransientIscalvingEnum; else if (strcmp(name,"TransientNumRequestedOutputs")==0) return TransientNumRequestedOutputsEnum; else if (strcmp(name,"TransientRequestedOutputs")==0) return TransientRequestedOutputsEnum; else if (strcmp(name,"Potential")==0) return PotentialEnum; @@ -378,14 +382,14 @@ else if (strcmp(name,"AdjointBalancethicknessAnalysis")==0) return AdjointBalancethicknessAnalysisEnum; else if (strcmp(name,"AdjointBalancethickness2Analysis")==0) return AdjointBalancethickness2AnalysisEnum; else if (strcmp(name,"AdjointHorizAnalysis")==0) return AdjointHorizAnalysisEnum; - else if (strcmp(name,"AnalysisCounter")==0) return AnalysisCounterEnum; + else stage=4; + } + if(stage==4){ + if (strcmp(name,"AnalysisCounter")==0) return AnalysisCounterEnum; else if (strcmp(name,"DefaultAnalysis")==0) return DefaultAnalysisEnum; else if (strcmp(name,"BalancethicknessAnalysis")==0) return BalancethicknessAnalysisEnum; else if (strcmp(name,"BalancethicknessSolution")==0) return BalancethicknessSolutionEnum; - else stage=4; - } - if(stage==4){ - if (strcmp(name,"Balancethickness2Analysis")==0) return Balancethickness2AnalysisEnum; + else if (strcmp(name,"Balancethickness2Analysis")==0) return Balancethickness2AnalysisEnum; else if (strcmp(name,"Balancethickness2Solution")==0) return Balancethickness2SolutionEnum; else if (strcmp(name,"BalancethicknessSoftAnalysis")==0) return BalancethicknessSoftAnalysisEnum; else if (strcmp(name,"BalancethicknessSoftSolution")==0) return BalancethicknessSoftSolutionEnum; @@ -501,14 +505,14 @@ else if (strcmp(name,"SpcTransient")==0) return SpcTransientEnum; else if (strcmp(name,"StringArrayParam")==0) return StringArrayParamEnum; else if (strcmp(name,"StringParam")==0) return StringParamEnum; - else if (strcmp(name,"Seg")==0) return SegEnum; + else stage=5; + } + if(stage==5){ + if (strcmp(name,"Seg")==0) return SegEnum; else if (strcmp(name,"SegInput")==0) return SegInputEnum; else if (strcmp(name,"Tria")==0) return TriaEnum; else if (strcmp(name,"TriaInput")==0) return TriaInputEnum; - else stage=5; - } - if(stage==5){ - if (strcmp(name,"Tetra")==0) return TetraEnum; + else if (strcmp(name,"Tetra")==0) return TetraEnum; else if (strcmp(name,"TetraInput")==0) return TetraInputEnum; else if (strcmp(name,"Penta")==0) return PentaEnum; else if (strcmp(name,"PentaInput")==0) return PentaInputEnum; @@ -624,14 +628,14 @@ else if (strcmp(name,"P2")==0) return P2Enum; else if (strcmp(name,"P2bubble")==0) return P2bubbleEnum; else if (strcmp(name,"P2bubblecondensed")==0) return P2bubblecondensedEnum; - else if (strcmp(name,"P2xP1")==0) return P2xP1Enum; + else stage=6; + } + if(stage==6){ + if (strcmp(name,"P2xP1")==0) return P2xP1Enum; else if (strcmp(name,"P1xP2")==0) return P1xP2Enum; else if (strcmp(name,"P1xP3")==0) return P1xP3Enum; else if (strcmp(name,"P2xP4")==0) return P2xP4Enum; - else stage=6; - } - if(stage==6){ - if (strcmp(name,"P1P1")==0) return P1P1Enum; + else if (strcmp(name,"P1P1")==0) return P1P1Enum; else if (strcmp(name,"P1P1GLS")==0) return P1P1GLSEnum; else if (strcmp(name,"MINI")==0) return MINIEnum; else if (strcmp(name,"MINIcondensed")==0) return MINIcondensedEnum; @@ -747,14 +751,14 @@ else if (strcmp(name,"Arrhenius")==0) return ArrheniusEnum; else if (strcmp(name,"LliboutryDuval")==0) return LliboutryDuvalEnum; else if (strcmp(name,"TransientIslevelset")==0) return TransientIslevelsetEnum; - else if (strcmp(name,"ExtrapolationVariable")==0) return ExtrapolationVariableEnum; + else stage=7; + } + if(stage==7){ + if (strcmp(name,"ExtrapolationVariable")==0) return ExtrapolationVariableEnum; else if (strcmp(name,"IceMaskNodeActivation")==0) return IceMaskNodeActivationEnum; else if (strcmp(name,"LevelsetfunctionSlopeX")==0) return LevelsetfunctionSlopeXEnum; else if (strcmp(name,"LevelsetfunctionSlopeY")==0) return LevelsetfunctionSlopeYEnum; - else stage=7; - } - if(stage==7){ - if (strcmp(name,"LevelsetfunctionPicard")==0) return LevelsetfunctionPicardEnum; + else if (strcmp(name,"LevelsetfunctionPicard")==0) return LevelsetfunctionPicardEnum; else if (strcmp(name,"Seaiceatm")==0) return SeaiceatmEnum; else if (strcmp(name,"Seaiceocean")==0) return SeaiceoceanEnum; else if (strcmp(name,"SeaiceThickness")==0) return SeaiceThicknessEnum; Index: ../trunk-jpl/src/c/modules/ModelProcessorx/CreateParameters.cpp =================================================================== --- ../trunk-jpl/src/c/modules/ModelProcessorx/CreateParameters.cpp (revision 18757) +++ ../trunk-jpl/src/c/modules/ModelProcessorx/CreateParameters.cpp (revision 18758) @@ -66,6 +66,7 @@ parameters->AddObject(iomodel->CopyConstantObject(QmuIsdakotaEnum)); parameters->AddObject(iomodel->CopyConstantObject(InversionIscontrolEnum)); parameters->AddObject(iomodel->CopyConstantObject(InversionTypeEnum)); + parameters->AddObject(iomodel->CopyConstantObject(CalvingLawEnum)); if(solution_type==SeaiceSolutionEnum){ } else{ @@ -82,8 +83,7 @@ parameters->AddObject(iomodel->CopyConstantObject(TransientIslevelsetEnum)); parameters->AddObject(iomodel->CopyConstantObject(TransientIsdamageevolutionEnum)); parameters->AddObject(iomodel->CopyConstantObject(TransientIshydrologyEnum)); - parameters->AddObject(iomodel->CopyConstantObject(MasstransportIscalvingrateEnum)); - parameters->AddObject(iomodel->CopyConstantObject(MasstransportLevermannCalvingCoeffEnum)); + parameters->AddObject(iomodel->CopyConstantObject(TransientIscalvingEnum)); parameters->AddObject(iomodel->CopyConstantObject(MaterialsRheologyLawEnum)); parameters->AddObject(iomodel->CopyConstantObject(GiaCrossSectionShapeEnum)); Index: ../trunk-jpl/src/c/cores/transient_core.cpp =================================================================== --- ../trunk-jpl/src/c/cores/transient_core.cpp (revision 18757) +++ ../trunk-jpl/src/c/cores/transient_core.cpp (revision 18758) @@ -20,11 +20,11 @@ /*parameters: */ IssmDouble starttime,finaltime,dt,yts; - bool isstressbalance,ismasstransport,isFS,isthermal,isgroundingline,isgia,islevelset,isdamageevolution,ishydrology,iscalvingrate; + bool isstressbalance,ismasstransport,isFS,isthermal,isgroundingline,isgia,islevelset,isdamageevolution,ishydrology,iscalving; bool save_results,dakota_analysis; bool time_adapt=false; int output_frequency; - int domaintype,groundingline_migration; + int domaintype,groundingline_migration,calvinglaw; int numoutputs = 0; Analysis *analysis = NULL; char** requested_outputs = NULL; @@ -52,7 +52,8 @@ femmodel->parameters->FindParam(&isdamageevolution,TransientIsdamageevolutionEnum); femmodel->parameters->FindParam(&ishydrology,TransientIshydrologyEnum); femmodel->parameters->FindParam(&isFS,FlowequationIsFSEnum); - femmodel->parameters->FindParam(&iscalvingrate,MasstransportIscalvingrateEnum); + femmodel->parameters->FindParam(&iscalving,TransientIscalvingEnum); + femmodel->parameters->FindParam(&calvinglaw,CalvingLawEnum); if(isgroundingline) femmodel->parameters->FindParam(&groundingline_migration,GroundinglineMigrationEnum); femmodel->parameters->FindParam(&numoutputs,TransientNumRequestedOutputsEnum); if(numoutputs) femmodel->parameters->FindParam(&requested_outputs,&numoutputs,TransientRequestedOutputsEnum); @@ -101,7 +102,7 @@ } if(islevelset){ - if(iscalvingrate){ + if(iscalving && calvinglaw==CalvingLevermannEnum){ if(VerboseSolution()) _printf0_(" computing calving rate\n"); femmodel->StrainRateparallelx(); femmodel->StrainRateperpendicularx(); Index: ../trunk-jpl/src/c/classes/Elements/Tria.cpp =================================================================== --- ../trunk-jpl/src/c/classes/Elements/Tria.cpp (revision 18757) +++ ../trunk-jpl/src/c/classes/Elements/Tria.cpp (revision 18758) @@ -438,7 +438,7 @@ Input* vy_input=inputs->GetInput(VyEnum); _assert_(vy_input); Input* strainparallel_input=inputs->GetInput(StrainRateparallelEnum); _assert_(strainparallel_input); Input* strainperpendicular_input=inputs->GetInput(StrainRateperpendicularEnum); _assert_(strainperpendicular_input); - this->parameters->FindParam(&propcoeff,MasstransportLevermannCalvingCoeffEnum); + Input* levermanncoeff_input=inputs->GetInput(CalvinglevermannCoeffEnum); _assert_(levermanncoeff_input); /* Start looping on the number of vertices: */ gauss=new GaussTria(); @@ -451,6 +451,7 @@ vel=vx*vx+vy*vy; strainparallel_input->GetInputValue(&strainparallel,gauss); strainperpendicular_input->GetInputValue(&strainperpendicular,gauss); + levermanncoeff_input->GetInputValue(&propcoeff,gauss); /*Calving rate proportionnal to the positive product of the strain rate along the ice flow direction and the strain rate perpendicular to the ice flow */ calvingrate[iv]=propcoeff*strainparallel*strainperpendicular; @@ -464,7 +465,7 @@ /*Add input*/ this->inputs->AddInput(new TriaInput(CalvingratexEnum,&calvingratex[0],P1Enum)); this->inputs->AddInput(new TriaInput(CalvingrateyEnum,&calvingratey[0],P1Enum)); - this->inputs->AddInput(new TriaInput(MasstransportCalvingrateEnum,&calvingrate[0],P1Enum)); + this->inputs->AddInput(new TriaInput(CalvingCalvingrateEnum,&calvingrate[0],P1Enum)); /*Clean up and return*/ delete gauss; Index: ../trunk-jpl/src/c/classes/Elements/Penta.cpp =================================================================== --- ../trunk-jpl/src/c/classes/Elements/Penta.cpp (revision 18757) +++ ../trunk-jpl/src/c/classes/Elements/Penta.cpp (revision 18758) @@ -366,8 +366,8 @@ Input* vx_input=inputs->GetInput(VxEnum); _assert_(vx_input); Input* vy_input=inputs->GetInput(VyEnum); _assert_(vy_input); Input* strainparallel_input=inputs->GetInput(StrainRateparallelEnum); _assert_(strainparallel_input); - Input* strainperpendicular_input=inputs->GetInput(StrainRateperpendicularEnum); _assert_(strainperpendicular_input); - this->parameters->FindParam(&propcoeff,MasstransportLevermannCalvingCoeffEnum); + Input* strainperpendicular_input=inputs->GetInput(StrainRateperpendicularEnum); _assert_(strainperpendicular_input); + Input* levermanncoeff_input=inputs->GetInput(CalvinglevermannCoeffEnum); _assert_(levermanncoeff_input); /* Start looping on the number of vertices: */ gauss=new GaussPenta(); @@ -380,6 +380,7 @@ vel=vx*vx+vy*vy; strainparallel_input->GetInputValue(&strainparallel,gauss); strainperpendicular_input->GetInputValue(&strainperpendicular,gauss); + levermanncoeff_input->GetInputValue(&propcoeff,gauss); /*Calving rate proportionnal to the positive product of the strain rate along the ice flow direction and the strain rate perpendicular to the ice flow */ calvingrate[iv]=propcoeff*strainparallel*strainperpendicular; @@ -393,7 +394,7 @@ /*Add input*/ this->inputs->AddInput(new PentaInput(CalvingratexEnum,&calvingratex[0],P1Enum)); this->inputs->AddInput(new PentaInput(CalvingrateyEnum,&calvingratey[0],P1Enum)); - this->inputs->AddInput(new PentaInput(MasstransportCalvingrateEnum,&calvingrate[0],P1Enum)); + this->inputs->AddInput(new PentaInput(CalvingCalvingrateEnum,&calvingrate[0],P1Enum)); /*Clean up and return*/ delete gauss; Index: ../trunk-jpl/src/c/classes/Elements/Element.cpp =================================================================== --- ../trunk-jpl/src/c/classes/Elements/Element.cpp (revision 18757) +++ ../trunk-jpl/src/c/classes/Elements/Element.cpp (revision 18758) @@ -1091,7 +1091,7 @@ name==LevelsetfunctionSlopeXEnum || name==LevelsetfunctionSlopeYEnum || name==LevelsetfunctionPicardEnum || - name==MasstransportCalvingrateEnum || + //name==CalvingCalvingrateEnum || name==GradientEnum || name==OldGradientEnum || name==ConvergedEnum || @@ -1211,7 +1211,7 @@ break; case CalvingratexEnum: case CalvingrateyEnum: - case MasstransportCalvingrateEnum: + case CalvingCalvingrateEnum: this->StrainRateparallel(); this->StrainRateperpendicular(); this->CalvingRateLevermann();