Index: /issm/trunk-jpl/src/c/analyses/DebrisAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/DebrisAnalysis.cpp (revision 27851) +++ /issm/trunk-jpl/src/c/analyses/DebrisAnalysis.cpp (revision 27852) @@ -10,7 +10,5 @@ #define FINITEELEMENT P1Enum -#define EPS 1e-8 -#define RHO 2 -//#define KAPPA 5 +#define EPS 1e-14 /*Model processing*/ @@ -117,5 +115,4 @@ parameters->AddObject(iomodel->CopyConstantObject("md.debris.min_thickness",DebrisMinThicknessEnum)); parameters->AddObject(iomodel->CopyConstantObject("md.debris.packingfraction",DebrisPackingFractionEnum)); - parameters->AddObject(iomodel->CopyConstantObject("md.debris.max_displacementvelocity",DebrisMaxdisplacementvelocityEnum)); parameters->AddObject(iomodel->CopyConstantObject("md.debris.removal_slope_threshold",DebrisRemovalSlopeThresholdEnum)); parameters->AddObject(iomodel->CopyConstantObject("md.debris.removal_stress_threshold",DebrisRemovalStressThresholdEnum)); @@ -131,10 +128,11 @@ void DebrisAnalysis::Core(FemModel* femmodel){/*{{{*/ - PreProcessing(femmodel); + //PreProcessing(femmodel); + //femmodel->parameters->SetParam(VxDebrisEnum,InputToExtrudeEnum); + //extrudefromtop_core(femmodel); femmodel->SetCurrentConfiguration(DebrisAnalysisEnum); solutionsequence_linear(femmodel); - femmodel->inputs->DuplicateInput(DebrisThicknessEnum,DebrisThicknessOldEnum); PostProcessing(femmodel); - + }/*}}}*/ void DebrisAnalysis::PreCore(FemModel* femmodel){/*{{{*/ @@ -199,11 +197,5 @@ Input* vy_input=NULL; if(dim>1){vy_input = topelement->GetInput(VyDebrisEnum); _assert_(vy_input);} - Input* slopex_input = topelement->GetInput(SurfaceSlopeXEnum); _assert_(slopex_input); - Input* slopey_input = topelement->GetInput(SurfaceSlopeYEnum); _assert_(slopey_input); - Input* thickness_input = topelement->GetInput(DebrisThicknessEnum); _assert_(thickness_input); h=topelement->CharacteristicLength(); - //IssmDouble hx,hy,hz; - //topelement->ElementSizes(&hx,&hy,&hz); - //h=sqrt(hx*hx+hy*hy); /* Start looping on the number of gaussian points: */ @@ -217,5 +209,5 @@ vx_input->GetInputValue(&vx,gauss); if(dim==2) vy_input->GetInputValue(&vy,gauss); - + D_scalar=gauss->weight*Jdet; for(int i=0;ivalues[i*numnodes+j] += D_scalar*basis[i]*basis[j]; @@ -223,8 +215,8 @@ /*Advection terms*/ vx_input->GetInputDerivativeValue(&dvx[0],xyz_list,gauss); - dvxdx=dvx[0]; D_scalar=dt*gauss->weight*Jdet; if(dim==2){ vy_input->GetInputDerivativeValue(&dvy[0],xyz_list,gauss); + dvxdx=dvx[0]; dvydy=dvy[1]; for(int i=0;iFindParam(&step,StepEnum); - IssmDouble slope_threshold; - topelement->FindParam(&slope_threshold,DebrisRemovalSlopeThresholdEnum); - IssmDouble rad2deg=180./M_PI; - IssmDouble kappa,f,smb,debristhickness;//,slopex; - IssmDouble M=1,C=-0.2,m,alpha_t_deg=25.0; - IssmDouble alpha_t_rad=tan(alpha_t_deg/rad2deg); - IssmDouble Diff,fraction; - Diff=3.2/3e7; - //Input* slopex_input=topelement->GetInput(SurfaceSlopeXEnum); _assert_(slopex_input); - Input* smb_input=topelement->GetInput(SmbMassBalanceEnum); _assert_(smb_input); - Input* debristhickness_input=topelement->GetInput(DebrisThicknessEnum); _assert_(debristhickness_input); - - /*if(isdisplacement){ - - slopex_input->GetInputValue(&slopex, gauss); - smb_input->GetInputValue(&smb, gauss); - debristhickness_input->GetInputValue(&debristhickness, gauss); - if((atan(fabs(slopex))*rad2deg)=alpha_t_deg & fabs(slopex)<=((1-C)/M)){ - m=1/(((1-C)/M)-alpha_t_rad); - f=m*fabs(slopex)-m*alpha_t_rad; - }else{ - f=1; - } - //f=1; - //kappa=-5.6e16*smb*debristhickness*f; - //kappa=debristhickness/h*4e9*f; - //kappa=14.2809e8*f; // 25° - kappa=247*f; // 35° - if(dim==2){ - for(int i=0;ivalues[i*numnodes+j] += D_scalar*kappa*( - dbasis[0*numnodes+j]*dbasis[0*numnodes+i] + dbasis[1*numnodes+j]*dbasis[1*numnodes+i] + dbasis[2*numnodes+j]*dbasis[2*numnodes+i] - ); - } - } - }else{ - for(int i=0;ivalues[i*numnodes+j] += D_scalar*kappa*(dbasis[0*numnodes+j]*dbasis[0*numnodes+i]); - } - } - } - } - /*******************************************************************/ - - IssmDouble rho,KAPPA,slopex,slopey, thickness; + IssmDouble rho; if(FINITEELEMENT==P1Enum){ - rho=RHO; + rho=2.; }else if(FINITEELEMENT==P2Enum){ - rho=RHO*2.; + rho=4.; } @@ -309,6 +249,5 @@ /*SSA*/ if(dim==1){ - //vx_input->GetInputValue(&vx,gauss); - vx_input->GetInputAverage(&vx); + vx_input->GetInputValue(&vx,gauss); D[0]=h/rho*fabs(vx); }else{ @@ -323,12 +262,9 @@ if(dim==1){ vx_input->GetInputValue(&vx,gauss); - //vx_input->GetInputAverage(&vx); vel=fabs(vx)+EPS; D[0] = h/(rho*vel)*vx*vx; }else{ vx_input->GetInputValue(&vx,gauss); - //vx_input->GetInputAverage(&vx); vy_input->GetInputValue(&vy,gauss); - //vy_input->GetInputAverage(&vy); vel=sqrt(vx*vx+vy*vy)+EPS; D[0*dim+0]=h/(rho*vel)*vx*vx; @@ -341,28 +277,9 @@ if(dim==1){ vx_input->GetInputValue(&vx,gauss); - //vx_input->GetInputAverage(&vx); - //tau=h/(rho*max(fabs(vx),EPS)); - tau=h/(rho*(fabs(vx)+EPS)); + tau=h/(rho*max(fabs(vx),EPS)); }else{ vx_input->GetInputValue(&vx,gauss); vy_input->GetInputValue(&vy,gauss); - vx_input->GetInputAverage(&vx); - vy_input->GetInputAverage(&vy); - vel=sqrt(vx*vx+vy*vy)+EPS; - - KAPPA=vel; - // adapt diffusion coefficient from here: http://websrv.cs.umt.edu/isis/index.php/Solving_the_equation_for_thickness_evolution - slopex_input->GetInputValue(&slopex,gauss); - slopey_input->GetInputValue(&slopey,gauss); - thickness_input->GetInputValue(&thickness,gauss); - //KAPPA=vel*thickness/(pow(pow(slopex,2)+pow(slopey,2),0.5)+EPS); - tau=h/(rho*vel)*(cosh(vel*h/2/KAPPA)/sinh(vel*h/2/KAPPA)-1./(vel*h/2/KAPPA)); - - /*kappa=KAPPA/yts; - if(vel*h/(3*2*kappa)<1){ - tau=pow(h,2)/(3*2*2*kappa); - } - else tau=h/(2*vel);*/ - //tau=dt*3; + tau=h/(rho*sqrt(vx*vx+vy*vy)+EPS); } } @@ -410,4 +327,5 @@ } + /*Advection matrix - part 2, B*/ for(int i=0;ivalues[i*numnodes+j]+=dt*gauss->weight*Jdet*tau*(basis[j]*dvxdx)*(basis[i]*dvxdx); + Ke->values[i*numnodes+j]+=dt*gauss->weight*Jdet*tau*(basis[j]*dvxdx)*(basis[i]*dvxdx);; } } @@ -524,6 +442,4 @@ Input* smb_input = topelement->GetInput(SmbMassBalanceEnum); _assert_(smb_input); Input* thickness_input = topelement->GetInput(DebrisThicknessEnum); _assert_(thickness_input); - Input* slopex_input = topelement->GetInput(SurfaceSlopeXEnum); _assert_(slopex_input); - Input* slopey_input = topelement->GetInput(SurfaceSlopeYEnum); _assert_(slopey_input); Input* vx_input = topelement->GetInput(VxDebrisEnum); _assert_(vx_input); Input* vy_input=NULL; @@ -532,13 +448,10 @@ } h=topelement->CharacteristicLength(); - //IssmDouble hx,hy,hz; - //topelement->ElementSizes(&hx,&hy,&hz); - //h=sqrt(hx*hx+hy*hy); - - IssmDouble rho,KAPPA,slopex,slopey; + + IssmDouble rho; if(FINITEELEMENT==P1Enum){ - rho=RHO; + rho=2.; }else if(FINITEELEMENT==P2Enum){ - rho=RHO*2.; + rho=4.; } @@ -553,9 +466,9 @@ thickness_input->GetInputValue(&thickness,gauss); if(smb>0.){ - psi=thickness-0*dt*smb*pf; + psi=thickness-0.*dt*smb*pf; }else{ psi=thickness-dt*smb*pf; } - //psi=thickness-dt*smb*pf; + for(int i=0;ivalues[i]+=Jdet*gauss->weight*psi*basis[i]; @@ -567,7 +480,5 @@ if(dim==1){ vx_input->GetInputValue(&vx,gauss); - //vx_input->GetInputAverage(&vx); - //tau=h/(rho*max(fabs(vx),EPS)); - tau=h/(rho*(fabs(vx)+EPS)); + tau=h/(rho*max(fabs(vx),EPS)); /*Force vector - part 2*/ @@ -583,24 +494,8 @@ vx_input->GetInputValue(&vx,gauss); vy_input->GetInputValue(&vy,gauss); - vx_input->GetInputAverage(&vx); - vy_input->GetInputAverage(&vy); vy_input->GetInputDerivativeValue(&dvy[0],xyz_list,gauss); - dvydy=dvy[1]; - - vel=sqrt(vx*vx+vy*vy)+EPS; - // tuned value - KAPPA=vel; - // adapt diffusion coefficient from here: http://websrv.cs.umt.edu/isis/index.php/Solving_the_equation_for_thickness_evolution - slopex_input->GetInputValue(&slopex,gauss); - slopey_input->GetInputValue(&slopey,gauss); - //KAPPA=vel*thickness/(pow(pow(slopex,2)+pow(slopey,2),0.5)+EPS); - tau=h/(rho*vel)*(cosh(vel*h/2/KAPPA)/sinh(vel*h/2/KAPPA)-1./(vel*h/2/KAPPA)); - - /*kappa=KAPPA/yts; - if(vel*h/(3*2*kappa)<1){ - tau=pow(h,2)/(3*2*2*kappa); - } - else tau=h/(2*vel);*/ - //tau=dt*3; + vel=sqrt(vx*vx+vy*vy); + dvydy=dvy[1]; + tau=h/(rho*vel+EPS); /*Force vector - part 2*/ @@ -636,11 +531,8 @@ int numnodes = element->GetNumberOfNodes(); IssmDouble* newthickness = xNew(numnodes); - IssmDouble* icethickness = xNew(numnodes); - /*Use the dof list to index into the solution vector: */ IssmDouble minthickness = element->FindParam(DebrisMinThicknessEnum); element->GetDofListLocal(&ddoflist,NoneApproximationEnum,GsetEnum); - element->GetInputListOnVertices(&icethickness[0],ThicknessEnum); for(int i=0;i(newthickness); - xDelete(icethickness); xDelete(ddoflist); //*/ @@ -708,4 +599,5 @@ femmodel->parameters->FindParam(&removalmodel,DebrisRemovalmodelEnum); Element* element= NULL; + Element* topelement = NULL; if(removalmodel==0){ @@ -724,9 +616,9 @@ IssmDouble* slopey = xNew(numnodes); IssmDouble* onsurface = xNew(numnodes); - //IssmDouble* ls_active = xNew(numnodes); - element->GetInputListOnNodes(debristhickness,DebrisThicknessOldEnum); + IssmDouble* ls_active = xNew(numnodes); + element->GetInputListOnNodes(debristhickness,DebrisThicknessEnum); element->GetInputListOnNodes(icethickness,ThicknessEnum); element->GetInputListOnNodes(onsurface,MeshVertexonsurfaceEnum); - //element->GetInputListOnNodes(ls_active,DebrisMaskNodeActivationEnum); + element->GetInputListOnNodes(ls_active,DebrisMaskNodeActivationEnum); dim=1; @@ -751,7 +643,6 @@ if(icethickness[k]<=(iceminthickness+0.00001)) kk++; } - //isminthicknessinelement=true; - //if(kkAddInput(DebrisThicknessEnum,debristhickness,P1DGEnum); + element->AddInput(DebrisThicknessEnum,debristhickness,P1Enum); xDelete(debristhickness); - xDelete(onsurface); xDelete(icethickness); xDelete(slopex); xDelete(slopey); - //xDelete(ls_active); break; } case 2:{ - IssmDouble stress_threshold=element->FindParam(DebrisRemovalStressThresholdEnum); - IssmDouble gravity=element->FindParam(ConstantsGEnum); + IssmDouble stress_threshold = element->FindParam(DebrisRemovalStressThresholdEnum); + IssmDouble gravity = element->FindParam(ConstantsGEnum); IssmDouble stress,rhod=1900.; int kk=0; @@ -787,21 +676,24 @@ if(icethickness[k]<=(iceminthickness+0.00001)) kk++; } - //isminthicknessinelement=true; - //if(isminthicknessinelement){ + isminthicknessinelement=true; if(kkstress_threshold) debristhickness[k]=0; + stress=rhod*gravity*debristhickness[k]*slope;//pow(slope*slope/(slope*slope+1),0.5);//sin(slope/rad2deg); + //stress_sum=stress_sum+stress; if(stress>stress_threshold) remove_debris=true; - } - if(remove_debris){ - for(k=0; kstress_threshold) debristhickness[k]=0.; + } + //if((stress_sum/double(kk))>stress_threshold) remove_debris=true; + if(remove_debris){ + for(k=0; kAddInput(DebrisThicknessEnum,debristhickness,P1DGEnum); + element->AddInput(DebrisThicknessEnum,debristhickness,P1Enum); xDelete(debristhickness); @@ -809,4 +701,5 @@ xDelete(slopex); xDelete(slopey); + xDelete(ls_active); break; } @@ -821,129 +714,70 @@ if(VerboseSolution()) _printf0_(" Debris preprocessing\n"); - int removalmodel, displacementmodel; - femmodel->parameters->FindParam(&removalmodel,DebrisRemovalmodelEnum); - femmodel->parameters->FindParam(&displacementmodel,DebrisDisplacementmodelEnum); - Element* element= NULL; for(Object* & object : femmodel->elements->objects){ element=xDynamicCast(object); - int numvertices=element->GetNumberOfVertices(); - int numnodes=element->GetNumberOfNodes(); - numvertices=numnodes; + int numvertices = element->GetNumberOfVertices(); IssmDouble* vx = xNew(numvertices); - IssmDouble* vy = xNew(numvertices); - IssmDouble* vx_new = xNew(numvertices); - IssmDouble* vy_new = xNew(numvertices); IssmDouble* debristhickness = xNew(numvertices); IssmDouble* slopex = xNew(numvertices); - IssmDouble* slopey = xNew(numvertices); - IssmDouble* icethickness = xNew(numvertices); - int dim=1,domaintype; - femmodel->parameters->FindParam(&domaintype,DomainTypeEnum); - if(domaintype!=Domain2DverticalEnum){ - dim=2; - } + IssmDouble* onsurface = xNew(numvertices); + IssmDouble* icethickness = xNew(numvertices); + element->GetInputListOnVertices(&debristhickness[0],DebrisThicknessEnum); - element->GetInputListOnVertices(&vx[0],VxEnum); + element->GetInputListOnVertices(&vx[0],VxDebrisEnum); element->GetInputListOnVertices(&slopex[0],SurfaceSlopeXEnum); + element->GetInputListOnVertices(&onsurface[0],MeshVertexonsurfaceEnum); element->GetInputListOnVertices(&icethickness[0],ThicknessEnum); - if(dim==2){ - element->GetInputListOnVertices(&vy[0],VyEnum); - element->GetInputListOnVertices(&slopey[0],SurfaceSlopeYEnum); - } - - IssmDouble rad2deg=180./M_PI; //=57.2958 - IssmDouble vslipx,vslipy,stress,rhod=1900.; - IssmDouble taud_plus, taud_minus, taudx, taudy, taud; + + IssmDouble slope,rad2deg=180./M_PI; //=57.2958 + IssmDouble vslipx,rhod=1900.; IssmDouble gravity=element->FindParam(ConstantsGEnum); IssmDouble slope_threshold=element->FindParam(DebrisRemovalSlopeThresholdEnum); - IssmDouble stress_threshold = element->FindParam(DebrisRemovalStressThresholdEnum); - IssmDouble maxv = element->FindParam(DebrisMaxdisplacementvelocityEnum); IssmDouble iceminthickness=element->FindParam(MasstransportMinThicknessEnum); - IssmDouble hx,hy,hz; - element->ElementSizes(&hx,&hy,&hz); int step; - IssmDouble dt; + IssmDouble dt, maxv; IssmDouble yts=31536000.; femmodel->parameters->FindParam(&step,StepEnum); femmodel->parameters->FindParam(&dt,TimesteppingTimeStepEnum); + bool isminthicknessinelement=false; - int kk=0; - - if(displacementmodel==1){ - + for(int i=0;i=taud_plus){ - vx[i]=vx[i]+vslipx; - if(dim==2) vy[i]=vy[i]+vslipy; - }else if(taud>=taud_minus & taudmaxv/yts) vx[i]=maxv/yts*fabs(vx[i])/(vx[i]+1e-10); - if(fabs(vy[i])>maxv/yts) vy[i]=maxv/yts*fabs(vy[i])/(vy[i]+1e-10); - } - } - element->AddInput(VxDebrisEnum,vx,P1Enum); - if(dim==2) element->AddInput(VyDebrisEnum,vy,P1Enum); - } + //if(onsurface[i]>.5){ + slope=fabs(slopex[i]); + //if((atan(slope)*rad2deg)>25.){ + //if(debristhickness[i]>0.01){ + vslipx=1.0/yts; + //maxv=10.0/2./dt; + //vslipx=-slope_threshold*rhod*gravity*debristhickness[i]*slopex[i]/yts; + vx[i]=vx[i]+vslipx; + //debristhickness[i]=debristhickness[i]; + //if(vx[i]>maxv) vx[i]=maxv; + //} + //} + //} + } + } + //if(step%100==0) + element->AddInput(VxDebrisEnum,vx,P1Enum); + /* Free resources */ xDelete(debristhickness); xDelete(icethickness); xDelete(vx); - xDelete(vx_new); xDelete(slopex); - xDelete(vy); - xDelete(vy_new); - xDelete(slopey); - } -}/*}}}*/ + xDelete(onsurface); + } + +}/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/SmbAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/SmbAnalysis.cpp (revision 27851) +++ /issm/trunk-jpl/src/c/analyses/SmbAnalysis.cpp (revision 27852) @@ -242,21 +242,7 @@ } break; - case SMBdebrisEvattEnum: - iomodel->FetchDataToInput(inputs,elements,"md.initialization.debris",DebrisThicknessEnum); - iomodel->FetchDataToInput(inputs,elements,"md.smb.s0t",SmbS0tEnum); - iomodel->FetchDataToInput(inputs,elements,"md.smb.snowheight",SmbSnowheightEnum); - iomodel->FetchDataToDatasetInput(inputs,elements,"md.smb.monthlytemperatures",SmbMonthlytemperaturesEnum); - iomodel->FetchDataToDatasetInput(inputs,elements,"md.smb.monthlyprecipitation",SmbPrecipitationEnum); - iomodel->FetchDataToDatasetInput(inputs,elements,"md.smb.monthlydsradiation",SmbMonthlydsradiationEnum); - iomodel->FetchDataToDatasetInput(inputs,elements,"md.smb.monthlydlradiation",SmbMonthlydlradiationEnum); - iomodel->FetchDataToDatasetInput(inputs,elements,"md.smb.monthlywindspeed",SmbMonthlywindspeedEnum); - iomodel->FetchDataToDatasetInput(inputs,elements,"md.smb.monthlyairhumidity",SmbMonthlyairhumidityEnum); - iomodel->FetchDataToInput(inputs,elements,"md.smb.precipitation_anomaly",SmbPrecipitationsAnomalyEnum); - iomodel->FetchDataToInput(inputs,elements,"md.smb.temperature_anomaly",SmbTemperaturesAnomalyEnum); - iomodel->FetchDataToInput(inputs,elements,"md.smb.dsradiation_anomaly",SmbDsradiationAnomalyEnum); - iomodel->FetchDataToInput(inputs,elements,"md.smb.dlradiation_anomaly",SmbDlradiationAnomalyEnum); - iomodel->FetchDataToInput(inputs,elements,"md.smb.windspeed_anomaly",SmbWindspeedAnomalyEnum); - iomodel->FetchDataToInput(inputs,elements,"md.smb.airhumidity_anomaly",SmbAirhumidityAnomalyEnum); - break; + case SMBdebrisMLEnum: + iomodel->FetchDataToInput(inputs,elements,"md.initialization.debris",DebrisThicknessEnum); + break; default: _error_("Surface mass balance model "< /*}}}*/ #define MAXVERTICES 6 /*Maximum number of vertices per element, currently Penta, to avoid dynamic mem allocation*/ @@ -2351,7 +2350,7 @@ IssmDouble minthickness = this->FindParam(MasstransportMinThicknessEnum); - + Input* input=this->GetInput(ThicknessEnum); _assert_(input); - if(input->GetInputMin()<=(minthickness+0.00000001)){ + if(input->GetInputMax()<=(minthickness+0.00000001)){ return true; } @@ -3401,5 +3400,5 @@ monthlytemperatures[iv*12+month]=monthlytemperatures[iv*12+month]-273.15; // conversion from Kelvin to celcius for PDD module dinput2->GetInputValue(&monthlyprec[iv*12+month],gauss,month); - monthlyprec[iv*12+month]=monthlyprec[iv*12+month];//*yts; + monthlyprec[iv*12+month]=monthlyprec[iv*12+month]*yts; } } @@ -3595,4 +3594,5 @@ Gauss* gauss=this->NewGauss(); for(int month=0;month<12;month++){ + for(int iv=0;ivGaussVertex(iv); @@ -3727,334 +3727,4 @@ xDelete(smbcorr); xDelete(melt_star); -} -/*}}}*/ -void Element::SmbDebrisEvatt(){/*{{{*/ - - const int NUM_VERTICES = this->GetNumberOfVertices(); - const int NUM_VERTICES_MONTHS_PER_YEAR = NUM_VERTICES * 365; - - int i,vertexlids[MAXVERTICES];; - IssmDouble* smb=xNew(NUM_VERTICES); // surface mass balance - IssmDouble* melt=xNew(NUM_VERTICES); // melting comp. of surface mass balance - IssmDouble* summermelt=xNew(NUM_VERTICES); - IssmDouble* albedo=xNew(NUM_VERTICES); // melting comp. of surface mass balance - IssmDouble* summeralbedo=xNew(NUM_VERTICES); - IssmDouble* accu=xNew(NUM_VERTICES); // accuumulation comp. of surface mass balance - IssmDouble* monthlytemperatures=xNew(NUM_VERTICES_MONTHS_PER_YEAR); - IssmDouble* monthlyprecip=xNew(NUM_VERTICES_MONTHS_PER_YEAR); - IssmDouble* monthlylw=xNew(NUM_VERTICES_MONTHS_PER_YEAR); - IssmDouble* monthlysw=xNew(NUM_VERTICES_MONTHS_PER_YEAR); - IssmDouble* monthlywind=xNew(NUM_VERTICES_MONTHS_PER_YEAR); - IssmDouble* monthlyhumidity=xNew(NUM_VERTICES_MONTHS_PER_YEAR); - IssmDouble* yearlytemperatures=xNew(NUM_VERTICES); memset(yearlytemperatures, 0., NUM_VERTICES*sizeof(IssmDouble)); - IssmDouble* p_ampl=xNew(NUM_VERTICES); // precip anomaly - IssmDouble* t_ampl=xNew(NUM_VERTICES); // temperature anomaly - IssmDouble* lw_ampl=xNew(NUM_VERTICES); - IssmDouble* sw_ampl=xNew(NUM_VERTICES); - IssmDouble* wind_ampl=xNew(NUM_VERTICES); - IssmDouble* humidity_ampl=xNew(NUM_VERTICES); - - IssmDouble* surface=xNew(NUM_VERTICES); // actual surface height - IssmDouble* s0t=xNew(NUM_VERTICES); // reference elevation - IssmDouble* snowheight=xNew(NUM_VERTICES); - IssmDouble* debriscover=xNew(NUM_VERTICES); - IssmDouble rho_water,rho_ice,Tf,debris,debris_here; - IssmDouble qlaps,rlaps,dslaps,dllaps,windspeedlaps,humiditylaps,Tm; - IssmDouble inv_twelve=1./365.; //factor for monthly average - IssmDouble time,yts,time_yr,lambda; - IssmDouble MonthlyMelt,CleanIceMonthlyMelt, CumMonthlyMelt=0,CleanIceMelt,CumMonthlySummerMelt=0; - IssmDouble MeanAlbedo=0, MeanSummerAlbedo=0; - bool isdebris,ismonthly,isAnderson,iscryokarst; - this->parameters->FindParam(&isdebris,TransientIsdebrisEnum); - this->parameters->FindParam(&ismonthly,SmbDebrisIsMonthlyEnum); - this->parameters->FindParam(&isAnderson,SmbDebrisIsAndersonEnum); - this->parameters->FindParam(&iscryokarst,SmbDebrisIsCryokarstEnum); - IssmDouble PhiD=0.,p; - IssmDouble icealbedo=this->FindParam(SmbIcealbedoEnum); - IssmDouble snowalbedo=this->FindParam(SmbSnowalbedoEnum); - IssmDouble debrisalbedo=this->FindParam(SmbDebrisalbedoEnum); - IssmDouble Lm=this->FindParam(MaterialsLatentheatEnum); - IssmDouble D0=this->FindParam(SmbDebrisAndersonD0Enum); - int step; - this->FindParam(&step,StepEnum); - - // cryokarst - int dim=1,domaintype; - this->parameters->FindParam(&domaintype,DomainTypeEnum); - if(domaintype!=Domain2DverticalEnum) dim=2; - IssmDouble taud_plus=110e3, taud_minus=60e3; - IssmDouble taud, slope, gravity, taudx, taudy; - this->parameters->FindParam(&gravity,ConstantsGEnum); - IssmDouble* slopex = xNew(NUM_VERTICES); - IssmDouble* slopey = xNew(NUM_VERTICES); - IssmDouble* icethickness = xNew(NUM_VERTICES); - - /*Get material parameters :*/ - rho_water=this->FindParam(MaterialsRhoSeawaterEnum); - rho_ice=this->FindParam(MaterialsRhoIceEnum); - IssmDouble sconv=(rho_water/rho_ice); - Tf=this->FindParam(MaterialsMeltingpointEnum); - - /*Get parameters for height corrections*/ - qlaps=this->FindParam(SmbDesfacEnum); // comment MR; on alpine galciers we dont have the dertification effect - rlaps=this->FindParam(SmbRlapsEnum); - dslaps=this->FindParam(SmbSWlapsEnum); - dllaps=this->FindParam(SmbLWlapsEnum); - windspeedlaps=this->FindParam(SmbWindspeedlapsEnum); - humiditylaps=this->FindParam(SmbHumiditylapsEnum); - - /* Get time */ - this->parameters->FindParam(&time,TimeEnum); - this->parameters->FindParam(&yts,ConstantsYtsEnum); - time_yr=floor(time/yts)*yts; - - /*Get inputs*/ - DatasetInput* tempmonthly =this->GetDatasetInput(SmbMonthlytemperaturesEnum); _assert_(tempmonthly); - DatasetInput* precipmonthly =this->GetDatasetInput(SmbPrecipitationEnum); _assert_(precipmonthly); - DatasetInput* lwmonthly =this->GetDatasetInput(SmbMonthlydlradiationEnum); _assert_(lwmonthly); - DatasetInput* swmonthly =this->GetDatasetInput(SmbMonthlydsradiationEnum); _assert_(swmonthly); - DatasetInput* windmonthly =this->GetDatasetInput(SmbMonthlywindspeedEnum); _assert_(windmonthly); - DatasetInput* humiditymonthly =this->GetDatasetInput(SmbMonthlyairhumidityEnum); _assert_(humiditymonthly); - - /*loop over vertices: */ - Gauss* gauss=this->NewGauss(); - for(int month=0;month<365;month++){ - for(int iv=0;ivGaussVertex(iv); - tempmonthly->GetInputValue(&monthlytemperatures[iv*365+month],gauss,month); - monthlytemperatures[iv*365+month]=monthlytemperatures[iv*365+month]-Tf; // conversion from Kelvin to celcius for PDD module - precipmonthly->GetInputValue(&monthlyprecip[iv*365+month],gauss,month); - monthlyprecip[iv*365+month]=monthlyprecip[iv*365+month]*yts; // from m/s to m/a - lwmonthly->GetInputValue(&monthlylw[iv*365+month],gauss,month); - swmonthly->GetInputValue(&monthlysw[iv*365+month],gauss,month); - windmonthly->GetInputValue(&monthlywind[iv*365+month],gauss,month); - humiditymonthly->GetInputValue(&monthlyhumidity[iv*365+month],gauss,month); - } - } - - /*Recover info at the vertices: */ - GetInputListOnVertices(&surface[0],SurfaceEnum); - GetInputListOnVertices(&s0t[0],SmbS0tEnum); - GetInputListOnVertices(&snowheight[0],SmbSnowheightEnum); - GetInputListOnVertices(&debriscover[0],DebrisThicknessEnum); - GetInputListOnVertices(&t_ampl[0],SmbTemperaturesAnomalyEnum); - GetInputListOnVertices(&p_ampl[0],SmbPrecipitationsAnomalyEnum); - GetInputListOnVertices(&lw_ampl[0],SmbDsradiationAnomalyEnum); - GetInputListOnVertices(&sw_ampl[0],SmbDlradiationAnomalyEnum); - GetInputListOnVertices(&wind_ampl[0],SmbWindspeedAnomalyEnum); - GetInputListOnVertices(&humidity_ampl[0],SmbAirhumidityAnomalyEnum); - if(iscryokarst){ - GetInputListOnVertices(&slopex[0],SurfaceSlopeXEnum); - GetInputListOnVertices(&icethickness[0],ThicknessEnum); - if(dim==2){ - GetInputListOnVertices(&slopey[0],SurfaceSlopeYEnum); - } - taudx=rho_ice*gravity*icethickness[i]*slopex[i]; - if(dim==2) taudy=rho_ice*gravity*icethickness[i]*slopey[i]; - taud=sqrt(taudx*taudx+taudy*taudy); - } - IssmDouble Alphaeff,Alphaeff_cleanice; - - /*measure the surface mass balance*/ - for (int iv = 0; iv i.e - /*Partition of precip in solid and liquid parts */ - IssmDouble temp_plus=1; - IssmDouble temp_minus=-1.; - IssmDouble frac_solid; - if(Tm>=temp_plus){ - frac_solid=0; - }else if(Tm<=temp_minus){ - frac_solid=1; - }else{ - frac_solid=1*(1-cos(PI*(temp_plus-Tm)/(temp_plus-temp_minus)))/2; - } - - /*Get yearly temperatures and accumulation */ - yearlytemperatures[iv]=yearlytemperatures[iv]+((monthlytemperatures[iv*365+imonth]-rlaps*st+Tf+t_ampl[iv]))/365; // Has to be in Kelvin - accu[iv]=accu[iv]+P*frac_solid; - if(yearlytemperatures[iv]>Tf) yearlytemperatures[iv]=Tf; - - CleanIceMonthlyMelt=((In-(Eps*Sigma*(Tf*Tf*Tf*Tf))+ - Q*(1.-Alphaeff)+ - (Rhoaa*Ca*Ustar*Ustar)/(Um-Ustar*(2.-(exp(Gamma*Xr))))*Tm)/((1-PhiD)*rho_ice*Lm)/(1.+ - ((Rhoaa*Ca*Ustar*Ustar)/(Um-Ustar*(2.-(exp(Gamma*Xr))))+4.*Eps*Sigma*(Tf*Tf*Tf))/ - K*debris)-(Lv*Ustar*Ustar*((Humidity))*(exp(-Gamma*Xr)))/((1.-PhiD)* - rho_ice*Lm*Ustar)/(((Um - -2.*Ustar)*exp(-Gamma*Xr))/Ustar+exp(Gamma*debris))); - if(CleanIceMonthlyMelt<0) CleanIceMonthlyMelt=0.; - MonthlyMelt=CleanIceMonthlyMelt; - - if(ismb==1){ - snowheight[iv]=snowheight[iv]+(P-CleanIceMonthlyMelt*yts/365); - if(snowheight[iv]<=0) snowheight[iv]=0.; - if(snowheight[iv]<=0.0001){ - p=debris_here*PhiD/(2*0.2*0.01); //Eq. 51 from Evatt et al 2015 without source term g*t - if(p>1.) p=1.; - if(p>=0.999){ - Alphaeff=debrisalbedo; - } else { - Alphaeff=Alphaeff_cleanice+p*(debrisalbedo-Alphaeff_cleanice); - } - debris=debris_here; - MonthlyMelt=((In-(Eps*Sigma*(Tf*Tf*Tf*Tf))+ - Q*(1.-Alphaeff)+ - (Rhoaa*Ca*Ustar*Ustar)/(Um-Ustar*(2.-(exp(Gamma*Xr))))*Tm)/((1-PhiD)*rho_ice*Lm)/(1.+ - ((Rhoaa*Ca*Ustar*Ustar)/(Um-Ustar*(2.-(exp(Gamma*Xr))))+4.*Eps*Sigma*(Tf*Tf*Tf))/ - K*debris)-(Lv*Ustar*Ustar*((Humidity))*(exp(-Gamma*Xr)))/((1.-PhiD)* - rho_ice*Lm*Ustar)/(((Um-2.*Ustar)*exp(-Gamma*Xr))/Ustar+exp(Gamma*debris))); - if(MonthlyMelt<0) MonthlyMelt=0.; - MeanSummerAlbedo=MeanSummerAlbedo+Alphaeff; - CumMonthlySummerMelt=CumMonthlySummerMelt+MonthlyMelt/365; - } - } - CumMonthlyMelt=CumMonthlyMelt+MonthlyMelt/365; - } - MeanAlbedo=MeanAlbedo+Alphaeff; - if(ismb==0) CleanIceMelt=CumMonthlyMelt; - } - - if(iscryokarst){ - if(taud>=taud_plus){ - lambda=0; - }else if(taud>=taud_minus & taudAddInput(SmbMassBalanceEnum,smb,P1Enum); - this->AddInput(SmbAccumulationEnum,accu,P1Enum); - this->AddInput(SmbMeltEnum,melt,P1Enum); - this->AddInput(SmbSummerMeltEnum,summermelt,P1Enum); - this->AddInput(SmbSnowheightEnum,snowheight,P1Enum); - this->AddInput(SmbAlbedoEnum,albedo,P1Enum); - this->AddInput(SmbSummerAlbedoEnum,summeralbedo,P1Enum); - this->AddInput(TemperaturePDDEnum,yearlytemperatures,P1Enum); // TemperaturePDD is wrong here, but don't want to create new Enum ... - - /*clean-up*/ - xDelete(monthlytemperatures); - xDelete(monthlyprecip); - xDelete(monthlylw); - xDelete(monthlysw); - xDelete(monthlywind); - xDelete(monthlyhumidity); - xDelete(smb); - xDelete(surface); - xDelete(melt); - xDelete(summermelt); - xDelete(albedo); - xDelete(summeralbedo); - xDelete(accu); - xDelete(yearlytemperatures); - xDelete(s0t); - xDelete(snowheight); - xDelete(debriscover); - xDelete(t_ampl); - xDelete(p_ampl); - xDelete(lw_ampl); - xDelete(sw_ampl); - xDelete(humidity_ampl); - xDelete(wind_ampl); - xDelete(slopex); - xDelete(slopey); - xDelete(icethickness); } /*}}}*/ Index: /issm/trunk-jpl/src/c/classes/Elements/Element.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Element.h (revision 27851) +++ /issm/trunk-jpl/src/c/classes/Elements/Element.h (revision 27852) @@ -177,5 +177,4 @@ void PositiveDegreeDay(IssmDouble* pdds,IssmDouble* pds,IssmDouble signorm,bool ismungsm,bool issetpddfac); void PositiveDegreeDaySicopolis(bool isfirnwarming); - void SmbDebrisEvatt(); void RignotMeltParameterization(); void ResultInterpolation(int* pinterpolation,int*nodesperelement,int* parray_size, int output_enum); @@ -257,5 +256,4 @@ virtual void ComputeStressTensor(void)=0; virtual void ComputeEsaStrainAndVorticity(void)=0; - //virtual void ComputeMeanEla(IssmDouble* paltitude,int* pcounter)=0; virtual void Configure(Elements* elements,Loads* loads,Nodes* nodes,Vertices* vertices,Materials* materials,Parameters* parameters,Inputs* inputsin)=0; virtual void ControlInputSetGradient(IssmDouble* gradient,int enum_type,int control_index,int offset,int M,int N,int interp)=0; Index: /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp (revision 27851) +++ /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp (revision 27852) @@ -5283,26 +5283,2 @@ /*}}}*/ #endif -//void Penta::ComputeMeanEla(IssmDouble* paltitude, int* pcounter){/*{{{*/ -// -// if(!this->IsOnSurface() & !this->IsIceInElement()) return; -// int ipos=0,ineg=0, counter=0; -// IssmDouble mean_element_surface=0; - - /* input */ -// IssmDouble smb[NUMVERTICES],surface[NUMVERTICES]; -// Element::GetInputListOnVertices(&smb[0],SmbMassBalanceEnum); -// Element::GetInputListOnVertices(&surface[0],SurfaceEnum); - -// for(int iv=0;iv0) ipos=ipos+1; -// if(smb[iv]<0) ineg=ineg+1; -// } - - /* we define ELA if an element has pos and neg smb */ -// if(ipos>0 & ineg>0){ -// for(int iv=0;ivparameters->FindParam(&domaintype,DomainTypeEnum); femmodel->parameters->FindParam(&save_results,SaveResultsEnum); femmodel->parameters->FindParam(&solution_type,SolutionTypeEnum); @@ -44,4 +43,12 @@ femmodel->inputs->DuplicateInput(VyEnum,VyDebrisEnum); } + femmodel->parameters->SetParam(VxEnum,InputToDepthaverageInEnum); + femmodel->parameters->SetParam(VxAverageEnum,InputToDepthaverageOutEnum); + depthaverage_core(femmodel); + if(domaintype!=Domain2DverticalEnum){ + femmodel->parameters->SetParam(VyEnum,InputToDepthaverageInEnum); + femmodel->parameters->SetParam(VyAverageEnum,InputToDepthaverageOutEnum); + depthaverage_core(femmodel); + } debris_analysis = new DebrisAnalysis(); @@ -52,9 +59,4 @@ extrudefromtop_core(femmodel); - IssmDouble mean_ela; - //ComputeMeanElax(&mean_ela,femmodel); - //femmodel->parameters->SetParam(mean_ela,SmbMeanElaEnum); - _printf_("core ELA: " << mean_ela << "\n"); - if(save_results) femmodel->RequestedOutputsx(&femmodel->results,requested_outputs,numoutputs); if(solution_type==DebrisSolutionEnum)femmodel->RequestedDependentsx(); Index: /issm/trunk-jpl/src/c/modules/ModelProcessorx/CreateParameters.cpp =================================================================== --- /issm/trunk-jpl/src/c/modules/ModelProcessorx/CreateParameters.cpp (revision 27851) +++ /issm/trunk-jpl/src/c/modules/ModelProcessorx/CreateParameters.cpp (revision 27852) @@ -489,18 +489,5 @@ parameters->AddObject(iomodel->CopyConstantObject("md.smb.ismethod",SmbSemicMethodEnum)); break; - case SMBdebrisEvattEnum: - parameters->AddObject(iomodel->CopyConstantObject("md.smb.qlaps",SmbDesfacEnum)); - parameters->AddObject(iomodel->CopyConstantObject("md.smb.rlaps",SmbRlapsEnum)); - parameters->AddObject(iomodel->CopyConstantObject("md.smb.dslaps",SmbSWlapsEnum)); - parameters->AddObject(iomodel->CopyConstantObject("md.smb.dllaps",SmbLWlapsEnum)); - parameters->AddObject(iomodel->CopyConstantObject("md.smb.windspeedlaps",SmbWindspeedlapsEnum)); - parameters->AddObject(iomodel->CopyConstantObject("md.smb.humiditylaps",SmbHumiditylapsEnum)); - parameters->AddObject(iomodel->CopyConstantObject("md.smb.icealbedo",SmbIcealbedoEnum)); - parameters->AddObject(iomodel->CopyConstantObject("md.smb.snowalbedo",SmbSnowalbedoEnum)); - parameters->AddObject(iomodel->CopyConstantObject("md.smb.debrisalbedo",SmbDebrisalbedoEnum)); - parameters->AddObject(iomodel->CopyConstantObject("md.smb.ismonthly",SmbDebrisIsMonthlyEnum)); - parameters->AddObject(iomodel->CopyConstantObject("md.smb.isAnderson",SmbDebrisIsAndersonEnum)); - parameters->AddObject(iomodel->CopyConstantObject("md.smb.iscryokarst",SmbDebrisIsCryokarstEnum)); - parameters->AddObject(iomodel->CopyConstantObject("md.smb.AndersonD0",SmbDebrisAndersonD0Enum)); + case SMBdebrisMLEnum: break; default: Index: /issm/trunk-jpl/src/c/modules/SurfaceMassBalancex/SurfaceMassBalancex.cpp =================================================================== --- /issm/trunk-jpl/src/c/modules/SurfaceMassBalancex/SurfaceMassBalancex.cpp (revision 27851) +++ /issm/trunk-jpl/src/c/modules/SurfaceMassBalancex/SurfaceMassBalancex.cpp (revision 27852) @@ -504,9 +504,145 @@ }/*}}}*/ -void SmbDebrisEvattx(FemModel* femmodel){/*{{{*/ - for(Object* & object : femmodel->elements->objects){ - Element* element=xDynamicCast(object); - element->SmbDebrisEvatt(); - } +void SmbDebrisMLx(FemModel* femmodel){/*{{{*/ + + // The function is based on: + // Evatt GW, Abrahams ID, Heil M, Mayer C, Kingslake J, Mitchell SL, et al. Glacial melt under a porous debris layer. Journal of Glaciology 61 (2015) 825–836, doi:10.3189/2 + // Constants/Values are taken from Mayer, Licciulli (2021): https://www.frontiersin.org/articles/10.3389/feart.2021.710276/full#B7 + // function taken from https://github.com/carlolic/DebrisExp/blob/main/USFs/USF_DebrisCoverage.f90 + + /*Intermediaries*/ + // altitude gradients of the crucial parameters (radiation from Marty et al., TaAClimat; 2002) + IssmDouble LW=2.9; // W/m^2 /100m 2.9 + IssmDouble SW=1.3; // W/m^2 /100m 1.3 + IssmDouble HumidityG=0; // % /100m rough estimate + IssmDouble AirTemp=0.7; // C /100m + IssmDouble WindSpeed=0.02; // m/s /100m rough estimate 0.2 + + // accumulation follows a linear increase above the ELA up to a plateau + IssmDouble AccG=0.1; // m w.e. /100m + IssmDouble AccMax=1.; // m w.e. + IssmDouble ReferenceElevation; + IssmDouble AblationDays=120.; // + + IssmDouble In=100.; // Wm^-2 incoming long wave + IssmDouble Q=500.; // Wm^-2 incoming short wave + IssmDouble K=0.585; // Wm^-1K^-1 thermal conductivity 0.585 + IssmDouble Qm=0.0012; // kg m^-3 measured humiditiy level + IssmDouble Qh=0.006 ; // kg m^-3 saturated humidity level + IssmDouble Tm=2.; // C air temperature + IssmDouble Rhoaa=1.22; // kgm^-3 air densitiy + IssmDouble Um=1.5; // ms^-1 measured wind speed + IssmDouble Xm=1.5; // ms^-1 measurement height + IssmDouble Xr=0.01; // ms^-1 surface roughness 0.01 + IssmDouble Alphad=0.07; // debris albedo 0.07 + IssmDouble Alphai=0.4; // ice ablbedo + IssmDouble Alphaeff; + IssmDouble Ustar=0.16; // ms^-1 friction velocity 0.16 + IssmDouble Ca=1000.; // jkg^-1K^-1 specific heat capacity of air + IssmDouble Lm;//=3.34E+05; // jkg^-1K^-1 latent heat of ice melt + IssmDouble Lv=2.50E+06; // jkg^-1K^-1 latent heat of evaporation + IssmDouble Tf=273.; // K water freeezing temperature + IssmDouble Eps=0.95; // thermal emissivity + IssmDouble Rhoi=900.; // kgm^-3 ice density + IssmDouble Sigma=5.67E-08; // Wm^-2K^-4 Stefan Boltzmann constant + IssmDouble Kstar=0.4; // von kármán constant + IssmDouble Gamma=180.; // m^-1 wind speed attenuation 234 + IssmDouble PhiD;//=0.005; // debris packing fraction 0.01 + IssmDouble Humidity=0.2; // relative humidity + + IssmDouble smb,yts,z,debris; + IssmDouble MassBalanceCmDayDebris,MassBalanceMYearDebris; + bool isdebris; + int domaintype; + femmodel->parameters->FindParam(&isdebris,TransientIsdebrisEnum); + + /*Get material parameters and constants */ + //femmodel->parameters->FindParam(&Rhoi,MaterialsRhoIceEnum); // Note Carlo's model used as benchmark was run with different densities for debris and FS + femmodel->parameters->FindParam(&Lm,MaterialsLatentheatEnum); + femmodel->parameters->FindParam(&yts,ConstantsYtsEnum); + PhiD=0.; + if(isdebris) femmodel->parameters->FindParam(&PhiD,DebrisPackingFractionEnum); + + /* Loop over all the elements of this partition */ + for(Object* & object : femmodel->elements->objects){ + Element* element=xDynamicCast(object); + + /* Allocate all arrays */ + int numvertices=element->GetNumberOfVertices(); + IssmDouble* surfacelist=xNew(numvertices); + IssmDouble* smb=xNew(numvertices); + IssmDouble* debriscover=xNew(numvertices); + element->GetInputListOnVertices(surfacelist,SurfaceEnum); + + /* Get inputs */ + element->GetInputListOnVertices(debriscover,DebrisThicknessEnum); + element->FindParam(&domaintype,DomainTypeEnum); + + /*Loop over all vertices of element and calculate SMB as function of Debris Cover and z */ + for(int v=0;vnmax){ + Alphaeff=Alphad; + } else { + Alphaeff=Alphai+n*(Alphad-Alphai)/nmax; + } + ReferenceElevation=3200.; // m HEF + + + Alphaeff=Alphad; + ReferenceElevation=2200.; // m M&L + + MassBalanceCmDayDebris=(((In-(z-ReferenceElevation)*LW/100.)-(Eps*Sigma*(Tf*Tf*Tf*Tf))+ + (Q+(z-ReferenceElevation)*SW/100.)*(1.-Alphaeff)+ + (Rhoaa*Ca*Ustar*Ustar)/((Um-(z-ReferenceElevation)* + WindSpeed/100.)-Ustar*(2.-(exp(Gamma*Xr))))*(Tm-(z- + ReferenceElevation)*AirTemp/100.))/((1-PhiD)*Rhoi*Lm)/(1.+ + ((Rhoaa*Ca*Ustar*Ustar)/((Um-(z-ReferenceElevation)* + WindSpeed/100.)-Ustar*(2.-(exp(Gamma*Xr))))+4.*Eps*Sigma*(Tf*Tf*Tf))/ + K*debris)-(Lv*Ustar*Ustar*(Qh-(Qh*(Humidity-(z- + ReferenceElevation)*HumidityG/100.)))*(exp(-Gamma*Xr)))/((1.-PhiD)* + Rhoi*Lm*Ustar)/((((Um-(z-ReferenceElevation)*WindSpeed/100.) + -2.*Ustar)*exp(-Gamma*Xr))/Ustar+exp(Gamma*debris)))*100.*24.*60.*60.; + } + + /* account form ablation days, and convert to m/s */ + MassBalanceMYearDebris=-MassBalanceCmDayDebris/100.*AblationDays/yts; + + /*Update array accordingly*/ + smb[v]=MassBalanceMYearDebris; + } + + /*Add input to element and Free memory*/ + element->AddInput(SmbMassBalanceEnum,smb,P1Enum); + xDelete(surfacelist); + xDelete(smb); + xDelete(debriscover); + } }/*}}}*/ void SmbGradientsComponentsx(FemModel* femmodel){/*{{{*/ Index: /issm/trunk-jpl/src/c/modules/SurfaceMassBalancex/SurfaceMassBalancex.h =================================================================== --- /issm/trunk-jpl/src/c/modules/SurfaceMassBalancex/SurfaceMassBalancex.h (revision 27851) +++ /issm/trunk-jpl/src/c/modules/SurfaceMassBalancex/SurfaceMassBalancex.h (revision 27852) @@ -23,5 +23,5 @@ void SmbMeltComponentsx(FemModel* femmodel); void SmbGradientsComponentsx(FemModel* femmodel); -void SmbDebrisEvattx(FemModel* femmodel); +void SmbDebrisMLx(FemModel* femmodel); /* SEMIC: */ void SmbSemicx(FemModel* femmodel, int ismethod); Index: /issm/trunk-jpl/src/c/modules/modules.h =================================================================== --- /issm/trunk-jpl/src/c/modules/modules.h (revision 27851) +++ /issm/trunk-jpl/src/c/modules/modules.h (revision 27852) @@ -86,5 +86,4 @@ #include "./ResetConstraintsx/ResetConstraintsx.h" #include "./ResetFSBasalBoundaryConditionx/ResetFSBasalBoundaryConditionx.h" -//#include "./ComputeMeanElax/ComputeMeanElax.h" #include "./RheologyBbarAbsGradientx/RheologyBbarAbsGradientx.h" #include "./RheologyBAbsGradientx/RheologyBAbsGradientx.h" Index: /issm/trunk-jpl/src/c/shared/Enum/Enum.vim =================================================================== --- /issm/trunk-jpl/src/c/shared/Enum/Enum.vim (revision 27851) +++ /issm/trunk-jpl/src/c/shared/Enum/Enum.vim (revision 27852) @@ -176,5 +176,4 @@ syn keyword cConstant DebrisRemovalStressThresholdEnum syn keyword cConstant DebrisPackingFractionEnum -syn keyword cConstant DebrisMaxdisplacementvelocityEnum syn keyword cConstant DebugProfilingEnum syn keyword cConstant DomainDimensionEnum @@ -540,9 +539,4 @@ syn keyword cConstant SmbARMAmaOrderEnum syn keyword cConstant SmbAveragingEnum -syn keyword cConstant SmbDebrisalbedoEnum -syn keyword cConstant SmbDebrisIsMonthlyEnum -syn keyword cConstant SmbDebrisIsAndersonEnum -syn keyword cConstant SmbDebrisIsCryokarstEnum -syn keyword cConstant SmbDebrisAndersonD0Enum syn keyword cConstant SmbDesfacEnum syn keyword cConstant SmbDesfacElevEnum @@ -558,7 +552,5 @@ syn keyword cConstant SmbEIdxEnum syn keyword cConstant SmbFEnum -syn keyword cConstant SmbHumiditylapsEnum syn keyword cConstant SmbInitDensityScalingEnum -syn keyword cConstant SmbIcealbedoEnum syn keyword cConstant SmbIsaccumulationEnum syn keyword cConstant SmbIsalbedoEnum @@ -580,8 +572,4 @@ syn keyword cConstant SmbKEnum syn keyword cConstant SmbLapseRatesEnum -syn keyword cConstant SmbLWlapsEnum -syn keyword cConstant SmbSWlapsEnum -syn keyword cConstant SmbSnowalbedoEnum -syn keyword cConstant SmbMeanElaEnum syn keyword cConstant SmbNumBasinsEnum syn keyword cConstant SmbNumBreaksEnum @@ -624,5 +612,4 @@ syn keyword cConstant SmbTemperaturesReconstructedYearsEnum syn keyword cConstant SmbPrecipitationsReconstructedYearsEnum -syn keyword cConstant SmbWindspeedlapsEnum syn keyword cConstant SmoothThicknessMultiplierEnum syn keyword cConstant SolutionTypeEnum @@ -792,5 +779,4 @@ syn keyword cConstant DamageFEnum syn keyword cConstant DebrisThicknessEnum -syn keyword cConstant DebrisThicknessOldEnum syn keyword cConstant DegreeOfChannelizationEnum syn keyword cConstant DepthBelowSurfaceEnum @@ -1121,8 +1107,4 @@ syn keyword cConstant SmbMeltEnum syn keyword cConstant SmbMonthlytemperaturesEnum -syn keyword cConstant SmbMonthlydsradiationEnum -syn keyword cConstant SmbMonthlydlradiationEnum -syn keyword cConstant SmbMonthlywindspeedEnum -syn keyword cConstant SmbMonthlyairhumidityEnum syn keyword cConstant SmbMSurfEnum syn keyword cConstant SmbNetLWEnum @@ -1134,8 +1116,4 @@ syn keyword cConstant SmbPrecipitationEnum syn keyword cConstant SmbPrecipitationsAnomalyEnum -syn keyword cConstant SmbDsradiationAnomalyEnum -syn keyword cConstant SmbDlradiationAnomalyEnum -syn keyword cConstant SmbWindspeedAnomalyEnum -syn keyword cConstant SmbAirhumidityAnomalyEnum syn keyword cConstant SmbPrecipitationsLgmEnum syn keyword cConstant SmbPrecipitationsPresentdayEnum @@ -1157,7 +1135,4 @@ syn keyword cConstant SmbSmbrefEnum syn keyword cConstant SmbSzaValueEnum -syn keyword cConstant SmbSummerMeltEnum -syn keyword cConstant SmbSummerAlbedoEnum -syn keyword cConstant SmbSnowheightEnum syn keyword cConstant SmbTEnum syn keyword cConstant SmbTaEnum @@ -1679,5 +1654,5 @@ syn keyword cConstant SMBarmaEnum syn keyword cConstant SMBcomponentsEnum -syn keyword cConstant SMBdebrisEvattEnum +syn keyword cConstant SMBdebrisMLEnum syn keyword cConstant SMBd18opddEnum syn keyword cConstant SMBforcingEnum @@ -1805,5 +1780,4 @@ syn keyword cType Cfsurfacesquaretransient syn keyword cType Channel -syn keyword cType classes syn keyword cType Constraint syn keyword cType Constraints @@ -1813,6 +1787,6 @@ syn keyword cType ControlParam syn keyword cType Covertree +syn keyword cType DataSetParam syn keyword cType DatasetInput -syn keyword cType DataSetParam syn keyword cType Definition syn keyword cType DependentObject @@ -1827,6 +1801,6 @@ syn keyword cType ElementInput syn keyword cType ElementMatrix +syn keyword cType ElementVector syn keyword cType Elements -syn keyword cType ElementVector syn keyword cType ExponentialVariogram syn keyword cType ExternalResult @@ -1835,10 +1809,9 @@ syn keyword cType Friction syn keyword cType Gauss -syn keyword cType GaussianVariogram -syn keyword cType gaussobjects syn keyword cType GaussPenta syn keyword cType GaussSeg syn keyword cType GaussTetra syn keyword cType GaussTria +syn keyword cType GaussianVariogram syn keyword cType GenericExternalResult syn keyword cType GenericOption @@ -1857,5 +1830,4 @@ syn keyword cType IssmDirectApplicInterface syn keyword cType IssmParallelDirectApplicInterface -syn keyword cType krigingobjects syn keyword cType Load syn keyword cType Loads @@ -1868,5 +1840,4 @@ syn keyword cType Matice syn keyword cType Matlitho -syn keyword cType matrixobjects syn keyword cType MatrixParam syn keyword cType Misfit @@ -1881,6 +1852,6 @@ syn keyword cType Observations syn keyword cType Option +syn keyword cType OptionUtilities syn keyword cType Options -syn keyword cType OptionUtilities syn keyword cType Param syn keyword cType Parameters @@ -1896,11 +1867,11 @@ syn keyword cType Regionaloutput syn keyword cType Results +syn keyword cType RiftStruct syn keyword cType Riftfront -syn keyword cType RiftStruct syn keyword cType SealevelGeometry syn keyword cType Seg syn keyword cType SegInput +syn keyword cType SegRef syn keyword cType Segment -syn keyword cType SegRef syn keyword cType SpcDynamic syn keyword cType SpcStatic @@ -1921,4 +1892,8 @@ syn keyword cType Vertex syn keyword cType Vertices +syn keyword cType classes +syn keyword cType gaussobjects +syn keyword cType krigingobjects +syn keyword cType matrixobjects syn keyword cType AdjointBalancethickness2Analysis syn keyword cType AdjointBalancethicknessAnalysis Index: /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h =================================================================== --- /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h (revision 27851) +++ /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h (revision 27852) @@ -170,5 +170,4 @@ DebrisRemovalStressThresholdEnum, DebrisPackingFractionEnum, - DebrisMaxdisplacementvelocityEnum, DebugProfilingEnum, DomainDimensionEnum, @@ -531,17 +530,12 @@ SmbAlbedoLandEnum, SmbARMATimestepEnum, - SmbARMAarOrderEnum, - SmbARMAmaOrderEnum, + SmbARMAarOrderEnum, + SmbARMAmaOrderEnum, SmbAveragingEnum, - SmbDebrisalbedoEnum, - SmbDebrisIsMonthlyEnum, - SmbDebrisIsAndersonEnum, - SmbDebrisIsCryokarstEnum, - SmbDebrisAndersonD0Enum, SmbDesfacEnum, SmbDesfacElevEnum, SmbDpermilEnum, SmbDsnowIdxEnum, - SmbElevationBinsEnum, + SmbElevationBinsEnum, SmbCldFracEnum, SmbDelta18oEnum, @@ -552,7 +546,5 @@ SmbEIdxEnum, SmbFEnum, - SmbHumiditylapsEnum, SmbInitDensityScalingEnum, - SmbIcealbedoEnum, SmbIsaccumulationEnum, SmbIsalbedoEnum, @@ -573,9 +565,5 @@ SmbIsturbulentfluxEnum, SmbKEnum, - SmbLapseRatesEnum, - SmbLWlapsEnum, - SmbSWlapsEnum, - SmbSnowalbedoEnum, - SmbMeanElaEnum, + SmbLapseRatesEnum, SmbNumBasinsEnum, SmbNumBreaksEnum, @@ -618,5 +606,4 @@ SmbTemperaturesReconstructedYearsEnum, SmbPrecipitationsReconstructedYearsEnum, - SmbWindspeedlapsEnum, SmoothThicknessMultiplierEnum, SolutionTypeEnum, @@ -788,5 +775,4 @@ DamageFEnum, DebrisThicknessEnum, - DebrisThicknessOldEnum, DegreeOfChannelizationEnum, DepthBelowSurfaceEnum, @@ -1057,5 +1043,5 @@ SmbAdiffiniEnum, SmbAiniEnum, - SmbARMANoiseEnum, + SmbARMANoiseEnum, SmbBasinsIdEnum, SmbBMaxEnum, @@ -1118,8 +1104,4 @@ SmbMInitnum, SmbMonthlytemperaturesEnum, - SmbMonthlydsradiationEnum, - SmbMonthlydlradiationEnum, - SmbMonthlywindspeedEnum, - SmbMonthlyairhumidityEnum, SmbMSurfEnum, SmbNetLWEnum, @@ -1131,8 +1113,4 @@ SmbPrecipitationEnum, SmbPrecipitationsAnomalyEnum, - SmbDsradiationAnomalyEnum, - SmbDlradiationAnomalyEnum, - SmbWindspeedAnomalyEnum, - SmbAirhumidityAnomalyEnum, SmbPrecipitationsLgmEnum, SmbPrecipitationsPresentdayEnum, @@ -1154,7 +1132,4 @@ SmbSmbrefEnum, SmbSzaValueEnum, - SmbSummerMeltEnum, - SmbSummerAlbedoEnum, - SmbSnowheightEnum, SmbTEnum, SmbTaEnum, @@ -1678,5 +1653,5 @@ SMBarmaEnum, SMBcomponentsEnum, - SMBdebrisEvattEnum, + SMBdebrisMLEnum, SMBd18opddEnum, SMBforcingEnum, Index: /issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp =================================================================== --- /issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp (revision 27851) +++ /issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp (revision 27852) @@ -178,5 +178,4 @@ case DebrisRemovalStressThresholdEnum : return "DebrisRemovalStressThreshold"; case DebrisPackingFractionEnum : return "DebrisPackingFraction"; - case DebrisMaxdisplacementvelocityEnum : return "DebrisMaxdisplacementvelocity"; case DebugProfilingEnum : return "DebugProfiling"; case DomainDimensionEnum : return "DomainDimension"; @@ -542,9 +541,4 @@ case SmbARMAmaOrderEnum : return "SmbARMAmaOrder"; case SmbAveragingEnum : return "SmbAveraging"; - case SmbDebrisalbedoEnum : return "SmbDebrisalbedo"; - case SmbDebrisIsMonthlyEnum : return "SmbDebrisIsMonthly"; - case SmbDebrisIsAndersonEnum : return "SmbDebrisIsAnderson"; - case SmbDebrisIsCryokarstEnum : return "SmbDebrisIsCryokarst"; - case SmbDebrisAndersonD0Enum : return "SmbDebrisAndersonD0"; case SmbDesfacEnum : return "SmbDesfac"; case SmbDesfacElevEnum : return "SmbDesfacElev"; @@ -560,7 +554,5 @@ case SmbEIdxEnum : return "SmbEIdx"; case SmbFEnum : return "SmbF"; - case SmbHumiditylapsEnum : return "SmbHumiditylaps"; case SmbInitDensityScalingEnum : return "SmbInitDensityScaling"; - case SmbIcealbedoEnum : return "SmbIcealbedo"; case SmbIsaccumulationEnum : return "SmbIsaccumulation"; case SmbIsalbedoEnum : return "SmbIsalbedo"; @@ -582,8 +574,4 @@ case SmbKEnum : return "SmbK"; case SmbLapseRatesEnum : return "SmbLapseRates"; - case SmbLWlapsEnum : return "SmbLWlaps"; - case SmbSWlapsEnum : return "SmbSWlaps"; - case SmbSnowalbedoEnum : return "SmbSnowalbedo"; - case SmbMeanElaEnum : return "SmbMeanEla"; case SmbNumBasinsEnum : return "SmbNumBasins"; case SmbNumBreaksEnum : return "SmbNumBreaks"; @@ -626,5 +614,4 @@ case SmbTemperaturesReconstructedYearsEnum : return "SmbTemperaturesReconstructedYears"; case SmbPrecipitationsReconstructedYearsEnum : return "SmbPrecipitationsReconstructedYears"; - case SmbWindspeedlapsEnum : return "SmbWindspeedlaps"; case SmoothThicknessMultiplierEnum : return "SmoothThicknessMultiplier"; case SolutionTypeEnum : return "SolutionType"; @@ -794,5 +781,4 @@ case DamageFEnum : return "DamageF"; case DebrisThicknessEnum : return "DebrisThickness"; - case DebrisThicknessOldEnum : return "DebrisThicknessOld"; case DegreeOfChannelizationEnum : return "DegreeOfChannelization"; case DepthBelowSurfaceEnum : return "DepthBelowSurface"; @@ -1123,8 +1109,4 @@ case SmbMeltEnum : return "SmbMelt"; case SmbMonthlytemperaturesEnum : return "SmbMonthlytemperatures"; - case SmbMonthlydsradiationEnum : return "SmbMonthlydsradiation"; - case SmbMonthlydlradiationEnum : return "SmbMonthlydlradiation"; - case SmbMonthlywindspeedEnum : return "SmbMonthlywindspeed"; - case SmbMonthlyairhumidityEnum : return "SmbMonthlyairhumidity"; case SmbMSurfEnum : return "SmbMSurf"; case SmbNetLWEnum : return "SmbNetLW"; @@ -1136,8 +1118,4 @@ case SmbPrecipitationEnum : return "SmbPrecipitation"; case SmbPrecipitationsAnomalyEnum : return "SmbPrecipitationsAnomaly"; - case SmbDsradiationAnomalyEnum : return "SmbDsradiationAnomaly"; - case SmbDlradiationAnomalyEnum : return "SmbDlradiationAnomaly"; - case SmbWindspeedAnomalyEnum : return "SmbWindspeedAnomaly"; - case SmbAirhumidityAnomalyEnum : return "SmbAirhumidityAnomaly"; case SmbPrecipitationsLgmEnum : return "SmbPrecipitationsLgm"; case SmbPrecipitationsPresentdayEnum : return "SmbPrecipitationsPresentday"; @@ -1159,7 +1137,4 @@ case SmbSmbrefEnum : return "SmbSmbref"; case SmbSzaValueEnum : return "SmbSzaValue"; - case SmbSummerMeltEnum : return "SmbSummerMelt"; - case SmbSummerAlbedoEnum : return "SmbSummerAlbedo"; - case SmbSnowheightEnum : return "SmbSnowheight"; case SmbTEnum : return "SmbT"; case SmbTaEnum : return "SmbTa"; @@ -1681,5 +1656,5 @@ case SMBarmaEnum : return "SMBarma"; case SMBcomponentsEnum : return "SMBcomponents"; - case SMBdebrisEvattEnum : return "SMBdebrisEvatt"; + case SMBdebrisMLEnum : return "SMBdebrisML"; case SMBd18opddEnum : return "SMBd18opdd"; case SMBforcingEnum : return "SMBforcing"; Index: /issm/trunk-jpl/src/c/shared/Enum/Enumjl.vim =================================================================== --- /issm/trunk-jpl/src/c/shared/Enum/Enumjl.vim (revision 27851) +++ /issm/trunk-jpl/src/c/shared/Enum/Enumjl.vim (revision 27852) @@ -169,5 +169,4 @@ syn keyword juliaConstC DebrisRemovalStressThresholdEnum syn keyword juliaConstC DebrisPackingFractionEnum -syn keyword juliaConstC DebrisMaxdisplacementvelocityEnum syn keyword juliaConstC DebugProfilingEnum syn keyword juliaConstC DomainDimensionEnum @@ -533,9 +532,4 @@ syn keyword juliaConstC SmbARMAmaOrderEnum syn keyword juliaConstC SmbAveragingEnum -syn keyword juliaConstC SmbDebrisalbedoEnum -syn keyword juliaConstC SmbDebrisIsMonthlyEnum -syn keyword juliaConstC SmbDebrisIsAndersonEnum -syn keyword juliaConstC SmbDebrisIsCryokarstEnum -syn keyword juliaConstC SmbDebrisAndersonD0Enum syn keyword juliaConstC SmbDesfacEnum syn keyword juliaConstC SmbDesfacElevEnum @@ -551,7 +545,5 @@ syn keyword juliaConstC SmbEIdxEnum syn keyword juliaConstC SmbFEnum -syn keyword juliaConstC SmbHumiditylapsEnum syn keyword juliaConstC SmbInitDensityScalingEnum -syn keyword juliaConstC SmbIcealbedoEnum syn keyword juliaConstC SmbIsaccumulationEnum syn keyword juliaConstC SmbIsalbedoEnum @@ -573,8 +565,4 @@ syn keyword juliaConstC SmbKEnum syn keyword juliaConstC SmbLapseRatesEnum -syn keyword juliaConstC SmbLWlapsEnum -syn keyword juliaConstC SmbSWlapsEnum -syn keyword juliaConstC SmbSnowalbedoEnum -syn keyword juliaConstC SmbMeanElaEnum syn keyword juliaConstC SmbNumBasinsEnum syn keyword juliaConstC SmbNumBreaksEnum @@ -617,5 +605,4 @@ syn keyword juliaConstC SmbTemperaturesReconstructedYearsEnum syn keyword juliaConstC SmbPrecipitationsReconstructedYearsEnum -syn keyword juliaConstC SmbWindspeedlapsEnum syn keyword juliaConstC SmoothThicknessMultiplierEnum syn keyword juliaConstC SolutionTypeEnum @@ -785,5 +772,4 @@ syn keyword juliaConstC DamageFEnum syn keyword juliaConstC DebrisThicknessEnum -syn keyword juliaConstC DebrisThicknessOldEnum syn keyword juliaConstC DegreeOfChannelizationEnum syn keyword juliaConstC DepthBelowSurfaceEnum @@ -1114,8 +1100,4 @@ syn keyword juliaConstC SmbMeltEnum syn keyword juliaConstC SmbMonthlytemperaturesEnum -syn keyword juliaConstC SmbMonthlydsradiationEnum -syn keyword juliaConstC SmbMonthlydlradiationEnum -syn keyword juliaConstC SmbMonthlywindspeedEnum -syn keyword juliaConstC SmbMonthlyairhumidityEnum syn keyword juliaConstC SmbMSurfEnum syn keyword juliaConstC SmbNetLWEnum @@ -1127,8 +1109,4 @@ syn keyword juliaConstC SmbPrecipitationEnum syn keyword juliaConstC SmbPrecipitationsAnomalyEnum -syn keyword juliaConstC SmbDsradiationAnomalyEnum -syn keyword juliaConstC SmbDlradiationAnomalyEnum -syn keyword juliaConstC SmbWindspeedAnomalyEnum -syn keyword juliaConstC SmbAirhumidityAnomalyEnum syn keyword juliaConstC SmbPrecipitationsLgmEnum syn keyword juliaConstC SmbPrecipitationsPresentdayEnum @@ -1150,7 +1128,4 @@ syn keyword juliaConstC SmbSmbrefEnum syn keyword juliaConstC SmbSzaValueEnum -syn keyword juliaConstC SmbSummerMeltEnum -syn keyword juliaConstC SmbSummerAlbedoEnum -syn keyword juliaConstC SmbSnowheightEnum syn keyword juliaConstC SmbTEnum syn keyword juliaConstC SmbTaEnum @@ -1672,5 +1647,5 @@ syn keyword juliaConstC SMBarmaEnum syn keyword juliaConstC SMBcomponentsEnum -syn keyword juliaConstC SMBdebrisEvattEnum +syn keyword juliaConstC SMBdebrisMLEnum syn keyword juliaConstC SMBd18opddEnum syn keyword juliaConstC SMBforcingEnum Index: /issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp =================================================================== --- /issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp (revision 27851) +++ /issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp (revision 27852) @@ -181,5 +181,4 @@ else if (strcmp(name,"DebrisRemovalStressThreshold")==0) return DebrisRemovalStressThresholdEnum; else if (strcmp(name,"DebrisPackingFraction")==0) return DebrisPackingFractionEnum; - else if (strcmp(name,"DebrisMaxdisplacementvelocity")==0) return DebrisMaxdisplacementvelocityEnum; else if (strcmp(name,"DebugProfiling")==0) return DebugProfilingEnum; else if (strcmp(name,"DomainDimension")==0) return DomainDimensionEnum; @@ -260,9 +259,9 @@ else if (strcmp(name,"Hydrologyarmadatebreaks")==0) return HydrologyarmadatebreaksEnum; else if (strcmp(name,"Hydrologyarmamalagcoefs")==0) return HydrologyarmamalagcoefsEnum; + else if (strcmp(name,"HydrologyarmamaOrder")==0) return HydrologyarmamaOrderEnum; else stage=3; } if(stage==3){ - if (strcmp(name,"HydrologyarmamaOrder")==0) return HydrologyarmamaOrderEnum; - else if (strcmp(name,"HydrologyarmaMonthlyFactors")==0) return HydrologyarmaMonthlyFactorsEnum; + if (strcmp(name,"HydrologyarmaMonthlyFactors")==0) return HydrologyarmaMonthlyFactorsEnum; else if (strcmp(name,"HydrologyarmaNumBreaks")==0) return HydrologyarmaNumBreaksEnum; else if (strcmp(name,"HydrologyarmaNumParams")==0) return HydrologyarmaNumParamsEnum; @@ -383,9 +382,9 @@ else if (strcmp(name,"MasstransportNumRequestedOutputs")==0) return MasstransportNumRequestedOutputsEnum; else if (strcmp(name,"MasstransportPenaltyFactor")==0) return MasstransportPenaltyFactorEnum; + else if (strcmp(name,"MasstransportRequestedOutputs")==0) return MasstransportRequestedOutputsEnum; else stage=4; } if(stage==4){ - if (strcmp(name,"MasstransportRequestedOutputs")==0) return MasstransportRequestedOutputsEnum; - else if (strcmp(name,"MasstransportStabilization")==0) return MasstransportStabilizationEnum; + if (strcmp(name,"MasstransportStabilization")==0) return MasstransportStabilizationEnum; else if (strcmp(name,"MaterialsBeta")==0) return MaterialsBetaEnum; else if (strcmp(name,"MaterialsEarthDensity")==0) return MaterialsEarthDensityEnum; @@ -506,9 +505,9 @@ else if (strcmp(name,"SolidearthSettingsRunFrequency")==0) return SolidearthSettingsRunFrequencyEnum; else if (strcmp(name,"SolidearthSettingsTimeAcc")==0) return SolidearthSettingsTimeAccEnum; + else if (strcmp(name,"SolidearthSettingsHoriz")==0) return SolidearthSettingsHorizEnum; else stage=5; } if(stage==5){ - if (strcmp(name,"SolidearthSettingsHoriz")==0) return SolidearthSettingsHorizEnum; - else if (strcmp(name,"SolidearthSettingsMaxiter")==0) return SolidearthSettingsMaxiterEnum; + if (strcmp(name,"SolidearthSettingsMaxiter")==0) return SolidearthSettingsMaxiterEnum; else if (strcmp(name,"SolidearthSettingsGrdOcean")==0) return SolidearthSettingsGrdOceanEnum; else if (strcmp(name,"SolidearthSettingsOceanAreaScaling")==0) return SolidearthSettingsOceanAreaScalingEnum; @@ -554,9 +553,4 @@ else if (strcmp(name,"SmbARMAmaOrder")==0) return SmbARMAmaOrderEnum; else if (strcmp(name,"SmbAveraging")==0) return SmbAveragingEnum; - else if (strcmp(name,"SmbDebrisalbedo")==0) return SmbDebrisalbedoEnum; - else if (strcmp(name,"SmbDebrisIsMonthly")==0) return SmbDebrisIsMonthlyEnum; - else if (strcmp(name,"SmbDebrisIsAnderson")==0) return SmbDebrisIsAndersonEnum; - else if (strcmp(name,"SmbDebrisIsCryokarst")==0) return SmbDebrisIsCryokarstEnum; - else if (strcmp(name,"SmbDebrisAndersonD0")==0) return SmbDebrisAndersonD0Enum; else if (strcmp(name,"SmbDesfac")==0) return SmbDesfacEnum; else if (strcmp(name,"SmbDesfacElev")==0) return SmbDesfacElevEnum; @@ -572,7 +566,5 @@ else if (strcmp(name,"SmbEIdx")==0) return SmbEIdxEnum; else if (strcmp(name,"SmbF")==0) return SmbFEnum; - else if (strcmp(name,"SmbHumiditylaps")==0) return SmbHumiditylapsEnum; else if (strcmp(name,"SmbInitDensityScaling")==0) return SmbInitDensityScalingEnum; - else if (strcmp(name,"SmbIcealbedo")==0) return SmbIcealbedoEnum; else if (strcmp(name,"SmbIsaccumulation")==0) return SmbIsaccumulationEnum; else if (strcmp(name,"SmbIsalbedo")==0) return SmbIsalbedoEnum; @@ -594,8 +586,4 @@ else if (strcmp(name,"SmbK")==0) return SmbKEnum; else if (strcmp(name,"SmbLapseRates")==0) return SmbLapseRatesEnum; - else if (strcmp(name,"SmbLWlaps")==0) return SmbLWlapsEnum; - else if (strcmp(name,"SmbSWlaps")==0) return SmbSWlapsEnum; - else if (strcmp(name,"SmbSnowalbedo")==0) return SmbSnowalbedoEnum; - else if (strcmp(name,"SmbMeanEla")==0) return SmbMeanElaEnum; else if (strcmp(name,"SmbNumBasins")==0) return SmbNumBasinsEnum; else if (strcmp(name,"SmbNumBreaks")==0) return SmbNumBreaksEnum; @@ -629,8 +617,5 @@ else if (strcmp(name,"SmbSemicTmax")==0) return SmbSemicTmaxEnum; else if (strcmp(name,"SmbStepsPerStep")==0) return SmbStepsPerStepEnum; - else stage=6; - } - if(stage==6){ - if (strcmp(name,"SmbSwIdx")==0) return SmbSwIdxEnum; + else if (strcmp(name,"SmbSwIdx")==0) return SmbSwIdxEnum; else if (strcmp(name,"SmbT0dry")==0) return SmbT0dryEnum; else if (strcmp(name,"SmbT0wet")==0) return SmbT0wetEnum; @@ -641,9 +626,11 @@ else if (strcmp(name,"SmbTemperaturesReconstructedYears")==0) return SmbTemperaturesReconstructedYearsEnum; else if (strcmp(name,"SmbPrecipitationsReconstructedYears")==0) return SmbPrecipitationsReconstructedYearsEnum; - else if (strcmp(name,"SmbWindspeedlaps")==0) return SmbWindspeedlapsEnum; else if (strcmp(name,"SmoothThicknessMultiplier")==0) return SmoothThicknessMultiplierEnum; else if (strcmp(name,"SolutionType")==0) return SolutionTypeEnum; else if (strcmp(name,"SteadystateMaxiter")==0) return SteadystateMaxiterEnum; - else if (strcmp(name,"SteadystateNumRequestedOutputs")==0) return SteadystateNumRequestedOutputsEnum; + else stage=6; + } + if(stage==6){ + if (strcmp(name,"SteadystateNumRequestedOutputs")==0) return SteadystateNumRequestedOutputsEnum; else if (strcmp(name,"SteadystateReltol")==0) return SteadystateReltolEnum; else if (strcmp(name,"SteadystateRequestedOutputs")==0) return SteadystateRequestedOutputsEnum; @@ -752,8 +739,5 @@ else if (strcmp(name,"BasalforcingsSpatialDeepwaterElevation")==0) return BasalforcingsSpatialDeepwaterElevationEnum; else if (strcmp(name,"BasalforcingsSpatialDeepwaterMeltingRate")==0) return BasalforcingsSpatialDeepwaterMeltingRateEnum; - else stage=7; - } - if(stage==7){ - if (strcmp(name,"BasalforcingsSpatialUpperwaterElevation")==0) return BasalforcingsSpatialUpperwaterElevationEnum; + else if (strcmp(name,"BasalforcingsSpatialUpperwaterElevation")==0) return BasalforcingsSpatialUpperwaterElevationEnum; else if (strcmp(name,"BasalforcingsSpatialUpperwaterMeltingRate")==0) return BasalforcingsSpatialUpperwaterMeltingRateEnum; else if (strcmp(name,"BasalforcingsIsmip6BasinId")==0) return BasalforcingsIsmip6BasinIdEnum; @@ -768,5 +752,8 @@ else if (strcmp(name,"BasalforcingsPicoOverturningCoeff")==0) return BasalforcingsPicoOverturningCoeffEnum; else if (strcmp(name,"BasalforcingsPicoSubShelfOceanOverturning")==0) return BasalforcingsPicoSubShelfOceanOverturningEnum; - else if (strcmp(name,"BasalforcingsPicoSubShelfOceanSalinity")==0) return BasalforcingsPicoSubShelfOceanSalinityEnum; + else stage=7; + } + if(stage==7){ + if (strcmp(name,"BasalforcingsPicoSubShelfOceanSalinity")==0) return BasalforcingsPicoSubShelfOceanSalinityEnum; else if (strcmp(name,"BasalforcingsPicoSubShelfOceanTemp")==0) return BasalforcingsPicoSubShelfOceanTempEnum; else if (strcmp(name,"BasalStressx")==0) return BasalStressxEnum; @@ -812,5 +799,4 @@ else if (strcmp(name,"DamageF")==0) return DamageFEnum; else if (strcmp(name,"DebrisThickness")==0) return DebrisThicknessEnum; - else if (strcmp(name,"DebrisThicknessOld")==0) return DebrisThicknessOldEnum; else if (strcmp(name,"DegreeOfChannelization")==0) return DegreeOfChannelizationEnum; else if (strcmp(name,"DepthBelowSurface")==0) return DepthBelowSurfaceEnum; @@ -875,8 +861,5 @@ else if (strcmp(name,"FrictionPressureAdjustedTemperature")==0) return FrictionPressureAdjustedTemperatureEnum; else if (strcmp(name,"FrictionQ")==0) return FrictionQEnum; - else stage=8; - } - if(stage==8){ - if (strcmp(name,"FrictionSedimentCompressibilityCoefficient")==0) return FrictionSedimentCompressibilityCoefficientEnum; + else if (strcmp(name,"FrictionSedimentCompressibilityCoefficient")==0) return FrictionSedimentCompressibilityCoefficientEnum; else if (strcmp(name,"FrictionTillFrictionAngle")==0) return FrictionTillFrictionAngleEnum; else if (strcmp(name,"FrictionWaterLayer")==0) return FrictionWaterLayerEnum; @@ -892,5 +875,8 @@ else if (strcmp(name,"UGia")==0) return UGiaEnum; else if (strcmp(name,"UGiaRate")==0) return UGiaRateEnum; - else if (strcmp(name,"Gradient")==0) return GradientEnum; + else stage=8; + } + if(stage==8){ + if (strcmp(name,"Gradient")==0) return GradientEnum; else if (strcmp(name,"GroundinglineHeight")==0) return GroundinglineHeightEnum; else if (strcmp(name,"HydraulicPotential")==0) return HydraulicPotentialEnum; @@ -998,8 +984,5 @@ else if (strcmp(name,"SamplingTau")==0) return SamplingTauEnum; else if (strcmp(name,"Sealevel")==0) return SealevelEnum; - else stage=9; - } - if(stage==9){ - if (strcmp(name,"SealevelGRD")==0) return SealevelGRDEnum; + else if (strcmp(name,"SealevelGRD")==0) return SealevelGRDEnum; else if (strcmp(name,"SatGraviGRD")==0) return SatGraviGRDEnum; else if (strcmp(name,"SealevelBarystaticMask")==0) return SealevelBarystaticMaskEnum; @@ -1015,5 +998,8 @@ else if (strcmp(name,"SealevelBarystaticHydroLatbar")==0) return SealevelBarystaticHydroLatbarEnum; else if (strcmp(name,"SealevelBarystaticHydroLongbar")==0) return SealevelBarystaticHydroLongbarEnum; - else if (strcmp(name,"SealevelBarystaticHydroLoad")==0) return SealevelBarystaticHydroLoadEnum; + else stage=9; + } + if(stage==9){ + if (strcmp(name,"SealevelBarystaticHydroLoad")==0) return SealevelBarystaticHydroLoadEnum; else if (strcmp(name,"SealevelBarystaticBpMask")==0) return SealevelBarystaticBpMaskEnum; else if (strcmp(name,"SealevelBarystaticBpWeights")==0) return SealevelBarystaticBpWeightsEnum; @@ -1121,8 +1107,5 @@ else if (strcmp(name,"SmbECDt")==0) return SmbECDtEnum; else if (strcmp(name,"SmbECini")==0) return SmbECiniEnum; - else stage=10; - } - if(stage==10){ - if (strcmp(name,"SmbEla")==0) return SmbElaEnum; + else if (strcmp(name,"SmbEla")==0) return SmbElaEnum; else if (strcmp(name,"SmbEvaporation")==0) return SmbEvaporationEnum; else if (strcmp(name,"SmbFAC")==0) return SmbFACEnum; @@ -1138,5 +1121,8 @@ else if (strcmp(name,"SmbIsInitialized")==0) return SmbIsInitializedEnum; else if (strcmp(name,"SmbMAdd")==0) return SmbMAddEnum; - else if (strcmp(name,"SmbMassBalance")==0) return SmbMassBalanceEnum; + else stage=10; + } + if(stage==10){ + if (strcmp(name,"SmbMassBalance")==0) return SmbMassBalanceEnum; else if (strcmp(name,"SmbMassBalanceSnow")==0) return SmbMassBalanceSnowEnum; else if (strcmp(name,"SmbMassBalanceIce")==0) return SmbMassBalanceIceEnum; @@ -1150,8 +1136,4 @@ else if (strcmp(name,"SmbMelt")==0) return SmbMeltEnum; else if (strcmp(name,"SmbMonthlytemperatures")==0) return SmbMonthlytemperaturesEnum; - else if (strcmp(name,"SmbMonthlydsradiation")==0) return SmbMonthlydsradiationEnum; - else if (strcmp(name,"SmbMonthlydlradiation")==0) return SmbMonthlydlradiationEnum; - else if (strcmp(name,"SmbMonthlywindspeed")==0) return SmbMonthlywindspeedEnum; - else if (strcmp(name,"SmbMonthlyairhumidity")==0) return SmbMonthlyairhumidityEnum; else if (strcmp(name,"SmbMSurf")==0) return SmbMSurfEnum; else if (strcmp(name,"SmbNetLW")==0) return SmbNetLWEnum; @@ -1163,8 +1145,4 @@ else if (strcmp(name,"SmbPrecipitation")==0) return SmbPrecipitationEnum; else if (strcmp(name,"SmbPrecipitationsAnomaly")==0) return SmbPrecipitationsAnomalyEnum; - else if (strcmp(name,"SmbDsradiationAnomaly")==0) return SmbDsradiationAnomalyEnum; - else if (strcmp(name,"SmbDlradiationAnomaly")==0) return SmbDlradiationAnomalyEnum; - else if (strcmp(name,"SmbWindspeedAnomaly")==0) return SmbWindspeedAnomalyEnum; - else if (strcmp(name,"SmbAirhumidityAnomaly")==0) return SmbAirhumidityAnomalyEnum; else if (strcmp(name,"SmbPrecipitationsLgm")==0) return SmbPrecipitationsLgmEnum; else if (strcmp(name,"SmbPrecipitationsPresentday")==0) return SmbPrecipitationsPresentdayEnum; @@ -1186,7 +1164,4 @@ else if (strcmp(name,"SmbSmbref")==0) return SmbSmbrefEnum; else if (strcmp(name,"SmbSzaValue")==0) return SmbSzaValueEnum; - else if (strcmp(name,"SmbSummerMelt")==0) return SmbSummerMeltEnum; - else if (strcmp(name,"SmbSummerAlbedo")==0) return SmbSummerAlbedoEnum; - else if (strcmp(name,"SmbSnowheight")==0) return SmbSnowheightEnum; else if (strcmp(name,"SmbT")==0) return SmbTEnum; else if (strcmp(name,"SmbTa")==0) return SmbTaEnum; @@ -1244,8 +1219,5 @@ else if (strcmp(name,"SurfaceCrevasse")==0) return SurfaceCrevasseEnum; else if (strcmp(name,"Surface")==0) return SurfaceEnum; - else stage=11; - } - if(stage==11){ - if (strcmp(name,"SurfaceOld")==0) return SurfaceOldEnum; + else if (strcmp(name,"SurfaceOld")==0) return SurfaceOldEnum; else if (strcmp(name,"SurfaceLogVelMisfit")==0) return SurfaceLogVelMisfitEnum; else if (strcmp(name,"SurfaceLogVxVyMisfit")==0) return SurfaceLogVxVyMisfitEnum; @@ -1272,5 +1244,8 @@ else if (strcmp(name,"ThicknessResidual")==0) return ThicknessResidualEnum; else if (strcmp(name,"TransientAccumulatedDeltaIceThickness")==0) return TransientAccumulatedDeltaIceThicknessEnum; - else if (strcmp(name,"Vel")==0) return VelEnum; + else stage=11; + } + if(stage==11){ + if (strcmp(name,"Vel")==0) return VelEnum; else if (strcmp(name,"VxAverage")==0) return VxAverageEnum; else if (strcmp(name,"VxBase")==0) return VxBaseEnum; @@ -1367,8 +1342,5 @@ else if (strcmp(name,"Outputdefinition5")==0) return Outputdefinition5Enum; else if (strcmp(name,"Outputdefinition60")==0) return Outputdefinition60Enum; - else stage=12; - } - if(stage==12){ - if (strcmp(name,"Outputdefinition61")==0) return Outputdefinition61Enum; + else if (strcmp(name,"Outputdefinition61")==0) return Outputdefinition61Enum; else if (strcmp(name,"Outputdefinition62")==0) return Outputdefinition62Enum; else if (strcmp(name,"Outputdefinition63")==0) return Outputdefinition63Enum; @@ -1395,5 +1367,8 @@ else if (strcmp(name,"Outputdefinition82")==0) return Outputdefinition82Enum; else if (strcmp(name,"Outputdefinition83")==0) return Outputdefinition83Enum; - else if (strcmp(name,"Outputdefinition84")==0) return Outputdefinition84Enum; + else stage=12; + } + if(stage==12){ + if (strcmp(name,"Outputdefinition84")==0) return Outputdefinition84Enum; else if (strcmp(name,"Outputdefinition85")==0) return Outputdefinition85Enum; else if (strcmp(name,"Outputdefinition86")==0) return Outputdefinition86Enum; @@ -1490,8 +1465,5 @@ else if (strcmp(name,"DefaultAnalysis")==0) return DefaultAnalysisEnum; else if (strcmp(name,"DefaultCalving")==0) return DefaultCalvingEnum; - else stage=13; - } - if(stage==13){ - if (strcmp(name,"Dense")==0) return DenseEnum; + else if (strcmp(name,"Dense")==0) return DenseEnum; else if (strcmp(name,"DependentObject")==0) return DependentObjectEnum; else if (strcmp(name,"DepthAverageAnalysis")==0) return DepthAverageAnalysisEnum; @@ -1518,5 +1490,8 @@ else if (strcmp(name,"ExtrapolationAnalysis")==0) return ExtrapolationAnalysisEnum; else if (strcmp(name,"ExtrudeFromBaseAnalysis")==0) return ExtrudeFromBaseAnalysisEnum; - else if (strcmp(name,"ExtrudeFromTopAnalysis")==0) return ExtrudeFromTopAnalysisEnum; + else stage=13; + } + if(stage==13){ + if (strcmp(name,"ExtrudeFromTopAnalysis")==0) return ExtrudeFromTopAnalysisEnum; else if (strcmp(name,"FSApproximation")==0) return FSApproximationEnum; else if (strcmp(name,"FSSolver")==0) return FSSolverEnum; @@ -1613,8 +1588,5 @@ else if (strcmp(name,"LoveKernelsImag")==0) return LoveKernelsImagEnum; else if (strcmp(name,"LoveKernelsReal")==0) return LoveKernelsRealEnum; - else stage=14; - } - if(stage==14){ - if (strcmp(name,"LoveKf")==0) return LoveKfEnum; + else if (strcmp(name,"LoveKf")==0) return LoveKfEnum; else if (strcmp(name,"LoveKt")==0) return LoveKtEnum; else if (strcmp(name,"LoveLf")==0) return LoveLfEnum; @@ -1641,5 +1613,8 @@ else if (strcmp(name,"Materials")==0) return MaterialsEnum; else if (strcmp(name,"Matestar")==0) return MatestarEnum; - else if (strcmp(name,"Matice")==0) return MaticeEnum; + else stage=14; + } + if(stage==14){ + if (strcmp(name,"Matice")==0) return MaticeEnum; else if (strcmp(name,"Matlitho")==0) return MatlithoEnum; else if (strcmp(name,"Mathydro")==0) return MathydroEnum; @@ -1720,5 +1695,5 @@ else if (strcmp(name,"SMBarma")==0) return SMBarmaEnum; else if (strcmp(name,"SMBcomponents")==0) return SMBcomponentsEnum; - else if (strcmp(name,"SMBdebrisEvatt")==0) return SMBdebrisEvattEnum; + else if (strcmp(name,"SMBdebrisML")==0) return SMBdebrisMLEnum; else if (strcmp(name,"SMBd18opdd")==0) return SMBd18opddEnum; else if (strcmp(name,"SMBforcing")==0) return SMBforcingEnum; @@ -1736,8 +1711,5 @@ else if (strcmp(name,"SSAFSApproximation")==0) return SSAFSApproximationEnum; else if (strcmp(name,"SSAHOApproximation")==0) return SSAHOApproximationEnum; - else stage=15; - } - if(stage==15){ - if (strcmp(name,"Scaled")==0) return ScaledEnum; + else if (strcmp(name,"Scaled")==0) return ScaledEnum; else if (strcmp(name,"SealevelAbsolute")==0) return SealevelAbsoluteEnum; else if (strcmp(name,"SealevelEmotion")==0) return SealevelEmotionEnum; @@ -1764,5 +1736,8 @@ else if (strcmp(name,"SpcTransient")==0) return SpcTransientEnum; else if (strcmp(name,"Sset")==0) return SsetEnum; - else if (strcmp(name,"StatisticsSolution")==0) return StatisticsSolutionEnum; + else stage=15; + } + if(stage==15){ + if (strcmp(name,"StatisticsSolution")==0) return StatisticsSolutionEnum; else if (strcmp(name,"SteadystateSolution")==0) return SteadystateSolutionEnum; else if (strcmp(name,"StressIntensityFactor")==0) return StressIntensityFactorEnum; Index: /issm/trunk-jpl/src/c/shared/io/Marshalling/IoCodeConversions.cpp =================================================================== --- /issm/trunk-jpl/src/c/shared/io/Marshalling/IoCodeConversions.cpp (revision 27851) +++ /issm/trunk-jpl/src/c/shared/io/Marshalling/IoCodeConversions.cpp (revision 27852) @@ -262,5 +262,5 @@ case 12: return SMBsemicEnum; case 13: return SMBarmaEnum; - case 14: return SMBdebrisEvattEnum; + case 14: return SMBdebrisMLEnum; default: _error_("Marshalled SMB code \""<=',0); md = checkfield(md,'fieldname','debris.packingfraction','>=',0); @@ -112,5 +108,4 @@ md = checkfield(md,'fieldname','debris.removal_slope_threshold','>=',0); md = checkfield(md,'fieldname','debris.removal_stress_threshold','>=',0); - md = checkfield(md,'fieldname','debris.max_displacementvelocity','>=',0); md = checkfield(md,'fieldname','debris.requested_outputs','stringrow',1); if ~any(isnan(md.stressbalance.vertex_pairing)), @@ -125,6 +120,5 @@ fielddisplay(self,'stabilization','0: no stabilization, 1: artificial diffusion, 2: streamline upwinding, 3: streamline upwind Petrov-Galerkin (SUPG)'); fielddisplay(self,'removalmodel','frontal removal of debris. 0: no removal, 1: Slope-triggered debris removal, 2: driving-stress triggered debris removal'); - fielddisplay(self,'max_displacementvelocity','maximum velocity of debris transport (v_ice + v_displacement) (m/a)'); - fielddisplay(self,'displacementmodel','debris displacement. 0: no displacement, 1: additional debris velocity above the critical slope/stress threshold'); + fielddisplay(self,'displacementmodel','debris displacement. 0: no displacement, 1: ...'); fielddisplay(self,'removal_slope_threshold','critical slope (degrees) for removalmodel (1)'); fielddisplay(self,'removal_stress_threshold','critical stress (Pa) for removalmodel (2)'); @@ -141,5 +135,4 @@ WriteData(fid,prefix,'object',self,'fieldname','removalmodel','format','Integer'); WriteData(fid,prefix,'object',self,'fieldname','displacementmodel','format','Integer'); - WriteData(fid,prefix,'object',self,'fieldname','max_displacementvelocity','format','Double'); WriteData(fid,prefix,'object',self,'fieldname','removal_slope_threshold','format','Double'); WriteData(fid,prefix,'object',self,'fieldname','removal_stress_threshold','format','Double'); @@ -163,5 +156,4 @@ writejsdouble(fid,[modelname '.debris.removalmodel'],self.removalmodel); writejsdouble(fid,[modelname '.debris.displacementmodel'],self.displacementmodel); - writejsdouble(fid,[modelname '.debris.max_displacementvelocity'],self.displacementmodel); writejsdouble(fid,[modelname '.debris.removal_slope_threshold'],self.removal_slope_threshold); writejsdouble(fid,[modelname '.debris.removal_stress_threshold'],self.removal_stress_threshold); Index: /issm/trunk-jpl/src/m/classes/model.m =================================================================== --- /issm/trunk-jpl/src/m/classes/model.m (revision 27851) +++ /issm/trunk-jpl/src/m/classes/model.m (revision 27852) @@ -1289,5 +1289,4 @@ if isfield(structmd,'cfl_coefficient'), md.timestepping.cfl_coefficient=structmd.cfl_coefficient; end if isfield(structmd,'spcthickness'), md.masstransport.spcthickness=structmd.spcthickness; end - if isfield(structmd,'spcthickness'), md.debris.spcthickness=structmd.spcthickness; end if isfield(structmd,'artificial_diffusivity'), md.masstransport.stabilization=structmd.artificial_diffusivity; end if isfield(structmd,'hydrostatic_adjustment'), md.masstransport.hydrostatic_adjustment=structmd.hydrostatic_adjustment; end Index: /issm/trunk-jpl/src/m/solve/parseresultsfromdisk.m =================================================================== --- /issm/trunk-jpl/src/m/solve/parseresultsfromdisk.m (revision 27851) +++ /issm/trunk-jpl/src/m/solve/parseresultsfromdisk.m (revision 27852) @@ -240,8 +240,4 @@ field = field*yts; elseif strcmp(fieldname,'VyAverage'), - field = field*yts; - elseif strcmp(fieldname,'VxDebris'), - field = field*yts; - elseif strcmp(fieldname,'VyDebris'), field = field*yts; elseif strcmp(fieldname,'BasalforcingsGroundediceMeltingRate'),