Index: /issm/trunk/examples/Bumps/Bump1_surface_bed/Square.par =================================================================== --- /issm/trunk/examples/Bumps/Bump1_surface_bed/Square.par (revision 8398) +++ /issm/trunk/examples/Bumps/Bump1_surface_bed/Square.par (revision 8399) @@ -68,6 +68,6 @@ md.rheology_n=3*ones(md.numberofelements,1); - disp(' creating accumulation rates'); - md.accumulation_rate=ones(md.numberofgrids,1)/md.yts; %1m/a + disp(' creating surface mass balance'); + md.surface_mass_balance=ones(md.numberofgrids,1)/md.yts; %1m/a md.basal_melting_rate=0*ones(md.numberofgrids,1)/md.yts; %1m/a Index: /issm/trunk/examples/Bumps/Bump1_surface_bed/runme.m =================================================================== --- /issm/trunk/examples/Bumps/Bump1_surface_bed/runme.m (revision 8398) +++ /issm/trunk/examples/Bumps/Bump1_surface_bed/runme.m (revision 8399) @@ -12,5 +12,5 @@ md.vx=PatchToVec(md.results.DiagnosticSolution.Vx); md.vy=PatchToVec(md.results.DiagnosticSolution.Vy); - md.accumulation_rate(:)=0; + md.surface_mass_balance(:)=0; md.basal_melting_rate(:)=0; md.thickness(:)=1; Index: /issm/trunk/examples/SquareIceShelf/runme.m =================================================================== --- /issm/trunk/examples/SquareIceShelf/runme.m (revision 8398) +++ /issm/trunk/examples/SquareIceShelf/runme.m (revision 8399) @@ -6,6 +6,6 @@ md.dt=1; md.ndt=3; -md.forcings.AccumulationRate=[ones(md.numberofnodes,1),2*ones(md.numberofnodes,1),3*ones(md.numberofnodes,1)]; -md.forcings.AccumulationRate(end+1,:)=[1.5 2.5 5.5]*md.yts; -md.forcings.MeltingRate=md.forcings.AccumulationRate; +md.forcings.surface_mass_balance=[ones(md.numberofnodes,1),2*ones(md.numberofnodes,1),3*ones(md.numberofnodes,1)]; +md.forcings.surface_mass_balance(end+1,:)=[1.5 2.5 5.5]*md.yts; +md.forcings.basal_melting_rate=md.forcings.surface_mass_balance; md=solve(md,Transient2DSolutionEnum); Index: /issm/trunk/src/c/EnumDefinitions/EnumDefinitions.h =================================================================== --- /issm/trunk/src/c/EnumDefinitions/EnumDefinitions.h (revision 8398) +++ /issm/trunk/src/c/EnumDefinitions/EnumDefinitions.h (revision 8399) @@ -175,5 +175,4 @@ /*}}}*/ /*Inputs {{{1*/ - AccumulationRateEnum, AdjointxEnum, AdjointyEnum, @@ -247,6 +246,9 @@ StabilizeConstraintsEnum, StokesReconditioningEnum, + SurfaceAccumulationRateEnum, + SurfaceAblationRateEnum, SurfaceAreaEnum, SurfaceEnum, + SurfaceMassBalanceEnum, SurfaceSlopeXEnum, SurfaceSlopeYEnum, Index: /issm/trunk/src/c/EnumDefinitions/EnumToModelField.cpp =================================================================== --- /issm/trunk/src/c/EnumDefinitions/EnumToModelField.cpp (revision 8398) +++ /issm/trunk/src/c/EnumDefinitions/EnumToModelField.cpp (revision 8399) @@ -23,6 +23,8 @@ case VyObsEnum : return "vy_obs"; case GroundingLineMigrationEnum : return "gl_migration"; - case AccumulationRateEnum: return "accumulation_rate"; case BasalMeltingRateEnum: return "basal_melting_rate"; + case SurfaceAccumulationRateEnum: return "surface_accumulation_rate"; + case SurfaceAblationRateEnum: return "surface_ablation_rate"; + case SurfaceMassBalanceEnum: return "surface_mass_balance"; default : _error_("No model field is associated to enum %s",EnumToStringx(en)); } Index: /issm/trunk/src/c/modules/EnumToStringx/EnumToStringx.cpp =================================================================== --- /issm/trunk/src/c/modules/EnumToStringx/EnumToStringx.cpp (revision 8398) +++ /issm/trunk/src/c/modules/EnumToStringx/EnumToStringx.cpp (revision 8399) @@ -143,5 +143,4 @@ case OpenEnum : return "Open"; case ClosedEnum : return "Closed"; - case AccumulationRateEnum : return "AccumulationRate"; case AdjointxEnum : return "Adjointx"; case AdjointyEnum : return "Adjointy"; @@ -215,6 +214,9 @@ case StabilizeConstraintsEnum : return "StabilizeConstraints"; case StokesReconditioningEnum : return "StokesReconditioning"; + case SurfaceAccumulationRateEnum : return "SurfaceAccumulationRate"; + case SurfaceAblationRateEnum : return "SurfaceAblationRate"; case SurfaceAreaEnum : return "SurfaceArea"; case SurfaceEnum : return "Surface"; + case SurfaceMassBalanceEnum : return "SurfaceMassBalance"; case SurfaceSlopeXEnum : return "SurfaceSlopeX"; case SurfaceSlopeYEnum : return "SurfaceSlopeY"; Index: /issm/trunk/src/c/modules/ModelProcessorx/Balancethickness/UpdateElementsBalancethickness.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/Balancethickness/UpdateElementsBalancethickness.cpp (revision 8398) +++ /issm/trunk/src/c/modules/ModelProcessorx/Balancethickness/UpdateElementsBalancethickness.cpp (revision 8399) @@ -30,5 +30,7 @@ IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy"); IoModelFetchData(&iomodel->basal_melting_rate,NULL,NULL,iomodel_handle,"basal_melting_rate"); - IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"accumulation_rate"); + IoModelFetchData(&iomodel->surface_accumulation_rate,NULL,NULL,iomodel_handle,"surface_accumulation_rate"); + IoModelFetchData(&iomodel->surface_ablation_rate,NULL,NULL,iomodel_handle,"surface_ablation_rate"); + IoModelFetchData(&iomodel->surface_mass_balance,NULL,NULL,iomodel_handle,"surface_mass_balance"); IoModelFetchData(&iomodel->dhdt,NULL,NULL,iomodel_handle,"dhdt"); if (iomodel->dim==3){ @@ -59,4 +61,6 @@ xfree((void**)&iomodel->vy); xfree((void**)&iomodel->basal_melting_rate); - xfree((void**)&iomodel->accumulation_rate); + xfree((void**)&iomodel->surface_accumulation_rate); + xfree((void**)&iomodel->surface_ablation_rate); + xfree((void**)&iomodel->surface_mass_balance); } Index: /issm/trunk/src/c/modules/ModelProcessorx/Balancevelocities/UpdateElementsBalancevelocities.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/Balancevelocities/UpdateElementsBalancevelocities.cpp (revision 8398) +++ /issm/trunk/src/c/modules/ModelProcessorx/Balancevelocities/UpdateElementsBalancevelocities.cpp (revision 8399) @@ -30,5 +30,7 @@ IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy"); IoModelFetchData(&iomodel->basal_melting_rate,NULL,NULL,iomodel_handle,"basal_melting_rate"); - IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"accumulation_rate"); + IoModelFetchData(&iomodel->surface_accumulation_rate,NULL,NULL,iomodel_handle,"surface_accumulation_rate"); + IoModelFetchData(&iomodel->surface_ablation_rate,NULL,NULL,iomodel_handle,"surface_ablation_rate"); + IoModelFetchData(&iomodel->surface_mass_balance,NULL,NULL,iomodel_handle,"surface_mass_balance"); if (iomodel->dim==3){ @@ -60,4 +62,6 @@ xfree((void**)&iomodel->vz); xfree((void**)&iomodel->basal_melting_rate); - xfree((void**)&iomodel->accumulation_rate); + xfree((void**)&iomodel->surface_accumulation_rate); + xfree((void**)&iomodel->surface_ablation_rate); + xfree((void**)&iomodel->surface_mass_balance); } Index: /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticHoriz/UpdateElementsDiagnosticHoriz.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticHoriz/UpdateElementsDiagnosticHoriz.cpp (revision 8398) +++ /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticHoriz/UpdateElementsDiagnosticHoriz.cpp (revision 8399) @@ -45,5 +45,7 @@ IoModelFetchData(&iomodel->elementonsurface,NULL,NULL,iomodel_handle,"elementonsurface"); IoModelFetchData(&iomodel->vz,NULL,NULL,iomodel_handle,"vz"); - IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"accumulation_rate"); + IoModelFetchData(&iomodel->surface_accumulation_rate,NULL,NULL,iomodel_handle,"surface_accumulation_rate"); + IoModelFetchData(&iomodel->surface_ablation_rate,NULL,NULL,iomodel_handle,"surface_ablation_rate"); + IoModelFetchData(&iomodel->surface_mass_balance,NULL,NULL,iomodel_handle,"surface_mass_balance"); IoModelFetchData(&iomodel->basal_melting_rate,NULL,NULL,iomodel_handle,"basal_melting_rate"); IoModelFetchData(&iomodel->nodeonstokes,NULL,NULL,iomodel_handle,"nodeonstokes"); @@ -81,5 +83,7 @@ xfree((void**)&iomodel->elementonsurface); xfree((void**)&iomodel->elementonwater); - xfree((void**)&iomodel->accumulation_rate); + xfree((void**)&iomodel->surface_accumulation_rate); + xfree((void**)&iomodel->surface_ablation_rate); + xfree((void**)&iomodel->surface_mass_balance); xfree((void**)&iomodel->basal_melting_rate); xfree((void**)&iomodel->nodeonstokes); Index: /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticVert/UpdateElementsDiagnosticVert.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticVert/UpdateElementsDiagnosticVert.cpp (revision 8398) +++ /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticVert/UpdateElementsDiagnosticVert.cpp (revision 8399) @@ -34,5 +34,7 @@ IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); IoModelFetchData(&iomodel->basal_melting_rate,NULL,NULL,iomodel_handle,"basal_melting_rate"); - IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"accumulation_rate"); + IoModelFetchData(&iomodel->surface_accumulation_rate,NULL,NULL,iomodel_handle,"surface_accumulation_rate"); + IoModelFetchData(&iomodel->surface_ablation_rate,NULL,NULL,iomodel_handle,"surface_ablation_rate"); + IoModelFetchData(&iomodel->surface_mass_balance,NULL,NULL,iomodel_handle,"surface_mass_balance"); IoModelFetchData(&iomodel->vx,NULL,NULL,iomodel_handle,"vx"); IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy"); @@ -60,5 +62,7 @@ xfree((void**)&iomodel->elementonwater); xfree((void**)&iomodel->basal_melting_rate); - xfree((void**)&iomodel->accumulation_rate); + xfree((void**)&iomodel->surface_accumulation_rate); + xfree((void**)&iomodel->surface_ablation_rate); + xfree((void**)&iomodel->surface_mass_balance); xfree((void**)&iomodel->vx); xfree((void**)&iomodel->vy); Index: /issm/trunk/src/c/modules/ModelProcessorx/Melting/UpdateElementsMelting.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/Melting/UpdateElementsMelting.cpp (revision 8398) +++ /issm/trunk/src/c/modules/ModelProcessorx/Melting/UpdateElementsMelting.cpp (revision 8399) @@ -38,5 +38,7 @@ IoModelFetchData(&iomodel->rheology_B,NULL,NULL,iomodel_handle,"rheology_B"); IoModelFetchData(&iomodel->rheology_n,NULL,NULL,iomodel_handle,"rheology_n"); - IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"accumulation_rate"); + IoModelFetchData(&iomodel->surface_accumulation_rate,NULL,NULL,iomodel_handle,"surface_accumulation_rate"); + IoModelFetchData(&iomodel->surface_ablation_rate,NULL,NULL,iomodel_handle,"surface_ablation_rate"); + IoModelFetchData(&iomodel->surface_mass_balance,NULL,NULL,iomodel_handle,"surface_mass_balance"); IoModelFetchData(&iomodel->basal_melting_rate,NULL,NULL,iomodel_handle,"basal_melting_rate"); IoModelFetchData(&iomodel->pressure,NULL,NULL,iomodel_handle,"pressure"); @@ -68,6 +70,8 @@ xfree((void**)&iomodel->rheology_B); xfree((void**)&iomodel->rheology_n); - xfree((void**)&iomodel->accumulation_rate); xfree((void**)&iomodel->basal_melting_rate); + xfree((void**)&iomodel->surface_accumulation_rate); + xfree((void**)&iomodel->surface_ablation_rate); + xfree((void**)&iomodel->surface_mass_balance); xfree((void**)&iomodel->pressure); } Index: /issm/trunk/src/c/modules/ModelProcessorx/Prognostic/UpdateElementsPrognostic.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/Prognostic/UpdateElementsPrognostic.cpp (revision 8398) +++ /issm/trunk/src/c/modules/ModelProcessorx/Prognostic/UpdateElementsPrognostic.cpp (revision 8399) @@ -27,5 +27,7 @@ IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf"); IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); - IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"accumulation_rate"); + IoModelFetchData(&iomodel->surface_accumulation_rate,NULL,NULL,iomodel_handle,"surface_accumulation_rate"); + IoModelFetchData(&iomodel->surface_ablation_rate,NULL,NULL,iomodel_handle,"surface_ablation_rate"); + IoModelFetchData(&iomodel->surface_mass_balance,NULL,NULL,iomodel_handle,"surface_mass_balance"); IoModelFetchData(&iomodel->basal_melting_rate,NULL,NULL,iomodel_handle,"basal_melting_rate"); if(iomodel->basal_melting_rate_correction_apply)IoModelFetchData(&iomodel->basal_melting_rate_correction,NULL,NULL,iomodel_handle,"basal_melting_rate_correction"); @@ -59,5 +61,7 @@ xfree((void**)&iomodel->elementonsurface); xfree((void**)&iomodel->elementonwater); - xfree((void**)&iomodel->accumulation_rate); + xfree((void**)&iomodel->surface_accumulation_rate); + xfree((void**)&iomodel->surface_ablation_rate); + xfree((void**)&iomodel->surface_mass_balance); xfree((void**)&iomodel->basal_melting_rate); xfree((void**)&iomodel->basal_melting_rate_correction); Index: /issm/trunk/src/c/modules/ModelProcessorx/Transient/UpdateElementsTransient.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/Transient/UpdateElementsTransient.cpp (revision 8398) +++ /issm/trunk/src/c/modules/ModelProcessorx/Transient/UpdateElementsTransient.cpp (revision 8399) @@ -49,5 +49,5 @@ /*we now have the forcing for the corresponding time step, for all the nodes. Use this - *to write over the existing accumulation input: */ + *to write over the existing smb input: */ counter=0; for (k=0;knumberofelements;k++){ Index: /issm/trunk/src/c/modules/StringToEnumx/StringToEnumx.cpp =================================================================== --- /issm/trunk/src/c/modules/StringToEnumx/StringToEnumx.cpp (revision 8398) +++ /issm/trunk/src/c/modules/StringToEnumx/StringToEnumx.cpp (revision 8399) @@ -141,5 +141,4 @@ else if (strcmp(name,"Open")==0) return OpenEnum; else if (strcmp(name,"Closed")==0) return ClosedEnum; - else if (strcmp(name,"AccumulationRate")==0) return AccumulationRateEnum; else if (strcmp(name,"Adjointx")==0) return AdjointxEnum; else if (strcmp(name,"Adjointy")==0) return AdjointyEnum; @@ -213,6 +212,9 @@ else if (strcmp(name,"StabilizeConstraints")==0) return StabilizeConstraintsEnum; else if (strcmp(name,"StokesReconditioning")==0) return StokesReconditioningEnum; + else if (strcmp(name,"SurfaceAccumulationRate")==0) return SurfaceAccumulationRateEnum; + else if (strcmp(name,"SurfaceAblationRate")==0) return SurfaceAblationRateEnum; else if (strcmp(name,"SurfaceArea")==0) return SurfaceAreaEnum; else if (strcmp(name,"Surface")==0) return SurfaceEnum; + else if (strcmp(name,"SurfaceMassBalance")==0) return SurfaceMassBalanceEnum; else if (strcmp(name,"SurfaceSlopeX")==0) return SurfaceSlopeXEnum; else if (strcmp(name,"SurfaceSlopeY")==0) return SurfaceSlopeYEnum; Index: /issm/trunk/src/c/objects/Elements/Penta.cpp =================================================================== --- /issm/trunk/src/c/objects/Elements/Penta.cpp (revision 8398) +++ /issm/trunk/src/c/objects/Elements/Penta.cpp (revision 8399) @@ -4007,7 +4007,15 @@ this->inputs->AddInput(new PentaVertexInput(BasalMeltingRateEnum,nodeinputs)); } - if (iomodel->accumulation_rate) { - for(i=0;i<6;i++)nodeinputs[i]=iomodel->accumulation_rate[penta_vertex_ids[i]-1]/iomodel->yts; - this->inputs->AddInput(new PentaVertexInput(AccumulationRateEnum,nodeinputs)); + if (iomodel->surface_accumulation_rate) { + for(i=0;i<6;i++)nodeinputs[i]=iomodel->surface_accumulation_rate[penta_vertex_ids[i]-1]/iomodel->yts; + this->inputs->AddInput(new PentaVertexInput(SurfaceAccumulationRateEnum,nodeinputs)); + } + if (iomodel->surface_ablation_rate) { + for(i=0;i<6;i++)nodeinputs[i]=iomodel->surface_ablation_rate[penta_vertex_ids[i]-1]/iomodel->yts; + this->inputs->AddInput(new PentaVertexInput(SurfaceAblationRateEnum,nodeinputs)); + } + if (iomodel->surface_mass_balance) { + for(i=0;i<6;i++)nodeinputs[i]=iomodel->surface_mass_balance[penta_vertex_ids[i]-1]/iomodel->yts; + this->inputs->AddInput(new PentaVertexInput(SurfaceMassBalanceEnum,nodeinputs)); } if (iomodel->geothermalflux) { @@ -5307,5 +5315,5 @@ name==SurfaceSlopeYEnum || name==BasalMeltingRateEnum || - name==AccumulationRateEnum || + name==SurfaceMassBalanceEnum || name==GeothermalFluxEnum || name==SurfaceAreaEnum|| @@ -6285,6 +6293,6 @@ if(step==0){ - /*This is the first time we are trying to replace the accumulation penta vertex input by an accumulation pentav vertex - *forcing, which is essentially a penta vertex input that can vary with time. So firt, if there is already + /*This is the first time we are trying to replace the smb penta vertex input by a smb pentav vertex + *forcing, which is essentially a penta vertex input that can vary with time. So first, if there is already an input with enum FieldEnum, squash it: */ this->inputs->AddInput(new PentaVertexForcing(FieldEnum,nodeinputs,time,iomodel->forcing_numtimesteps,this->parameters)); Index: /issm/trunk/src/c/objects/Elements/Tria.cpp =================================================================== --- /issm/trunk/src/c/objects/Elements/Tria.cpp (revision 8398) +++ /issm/trunk/src/c/objects/Elements/Tria.cpp (revision 8399) @@ -1661,5 +1661,5 @@ int i,j,ig; double xyz_list[NUMVERTICES][3]; - double dhdt_g,melting_g,accumulation_g,Jdettria; + double dhdt_g,basal_melting_g,surface_mass_balance_g,Jdettria; double L[NUMVERTICES]; GaussTria* gauss=NULL; @@ -1670,6 +1670,6 @@ /*Retrieve all inputs and parameters*/ GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES); - Input* accumulation_input=inputs->GetInput(AccumulationRateEnum); _assert_(accumulation_input); - Input* melting_input=inputs->GetInput(BasalMeltingRateEnum); _assert_(melting_input); + Input* surface_mass_balance_input=inputs->GetInput(SurfaceMassBalanceEnum); _assert_(surface_mass_balance_input); + Input* basal_melting_input=inputs->GetInput(BasalMeltingRateEnum); _assert_(basal_melting_input); Input* dhdt_input=inputs->GetInput(DhDtEnum); _assert_(dhdt_input); @@ -1680,6 +1680,6 @@ gauss->GaussPoint(ig); - accumulation_input->GetParameterValue(&accumulation_g,gauss); - melting_input->GetParameterValue(&melting_g,gauss); + surface_mass_balance_input->GetParameterValue(&surface_mass_balance_g,gauss); + basal_melting_input->GetParameterValue(&basal_melting_g,gauss); dhdt_input->GetParameterValue(&dhdt_g,gauss); @@ -1687,5 +1687,5 @@ GetL(&L[0], &xyz_list[0][0], gauss,NDOF1); - for(i=0;ivalues[i]+=Jdettria*gauss->weight*(accumulation_g-melting_g-dhdt_g)*L[i]; + for(i=0;ivalues[i]+=Jdettria*gauss->weight*(surface_mass_balance_g-basal_melting_g-dhdt_g)*L[i]; } @@ -1704,5 +1704,5 @@ int i,j,ig; double xyz_list[NUMVERTICES][3]; - double melting_g,accumulation_g,dhdt_g,Jdettria; + double basal_melting_g,surface_mass_balance_g,dhdt_g,Jdettria; double L[NUMVERTICES]; GaussTria* gauss=NULL; @@ -1713,6 +1713,6 @@ /*Retrieve all inputs and parameters*/ GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES); - Input* accumulation_input=inputs->GetInput(AccumulationRateEnum); _assert_(accumulation_input); - Input* melting_input=inputs->GetInput(BasalMeltingRateEnum); _assert_(melting_input); + Input* surface_mass_balance_input=inputs->GetInput(SurfaceMassBalanceEnum); _assert_(surface_mass_balance_input); + Input* basal_melting_input=inputs->GetInput(BasalMeltingRateEnum); _assert_(basal_melting_input); Input* dhdt_input=inputs->GetInput(DhDtEnum); _assert_(dhdt_input); @@ -1723,6 +1723,6 @@ gauss->GaussPoint(ig); - accumulation_input->GetParameterValue(&accumulation_g,gauss); - melting_input->GetParameterValue(&melting_g,gauss); + surface_mass_balance_input->GetParameterValue(&surface_mass_balance_g,gauss); + basal_melting_input->GetParameterValue(&basal_melting_g,gauss); dhdt_input->GetParameterValue(&dhdt_g,gauss); @@ -1730,5 +1730,5 @@ GetL(&L[0], &xyz_list[0][0], gauss,NDOF1); - for(i=0;ivalues[i]+=Jdettria*gauss->weight*(accumulation_g-melting_g-dhdt_g)*L[i]; + for(i=0;ivalues[i]+=Jdettria*gauss->weight*(surface_mass_balance_g-basal_melting_g-dhdt_g)*L[i]; } @@ -1747,5 +1747,5 @@ int i,j,ig; double xyz_list[NUMVERTICES][3]; - double Jdettria,accumulation_g,melting_g; + double Jdettria,surface_mass_balance_g,basal_melting_g; double L[NUMVERTICES]; GaussTria* gauss=NULL; @@ -1756,6 +1756,6 @@ /*Retrieve all inputs and parameters*/ GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES); - Input* accumulation_input=inputs->GetInput(AccumulationRateEnum); _assert_(accumulation_input); - Input* melting_input=inputs->GetInput(BasalMeltingRateEnum); _assert_(melting_input); + Input* surface_mass_balance_input=inputs->GetInput(SurfaceMassBalanceEnum); _assert_(surface_mass_balance_input); + Input* basal_melting_input=inputs->GetInput(BasalMeltingRateEnum); _assert_(basal_melting_input); /* Start looping on the number of gaussian points: */ @@ -1768,8 +1768,8 @@ GetL(&L[0], &xyz_list[0][0], gauss,NDOF1); - accumulation_input->GetParameterValue(&accumulation_g,gauss); - melting_input->GetParameterValue(&melting_g,gauss); - - for(i=0;ivalues[i]+=Jdettria*gauss->weight*(accumulation_g-melting_g)*L[i]; + surface_mass_balance_input->GetParameterValue(&surface_mass_balance_g,gauss); + basal_melting_input->GetParameterValue(&basal_melting_g,gauss); + + for(i=0;ivalues[i]+=Jdettria*gauss->weight*(surface_mass_balance_g-basal_melting_g)*L[i]; } @@ -1790,5 +1790,5 @@ double xyz_list[NUMVERTICES][3]; double Jdet; - double vx,vy,vz,dbdx,dbdy,meltingvalue; + double vx,vy,vz,dbdx,dbdy,basalmeltingvalue; double slope[2]; double L[NUMVERTICES]; @@ -1802,5 +1802,5 @@ inputs->GetParameterValue(&approximation,ApproximationEnum); Input* bed_input=inputs->GetInput(BedEnum); _assert_(bed_input); - Input* melting_input=inputs->GetInput(BasalMeltingRateEnum); _assert_(melting_input); + Input* basal_melting_input=inputs->GetInput(BasalMeltingRateEnum); _assert_(basal_melting_input); Input* vx_input=inputs->GetInput(VxEnum); _assert_(vx_input); Input* vy_input=inputs->GetInput(VyEnum); _assert_(vy_input); @@ -1816,5 +1816,5 @@ gauss->GaussPoint(ig); - melting_input->GetParameterValue(&meltingvalue, gauss); + basal_melting_input->GetParameterValue(&basalmeltingvalue, gauss); bed_input->GetParameterDerivativeValue(&slope[0],&xyz_list[0][0],gauss); vx_input->GetParameterValue(&vx, gauss); @@ -1831,5 +1831,5 @@ GetL(&L[0], &xyz_list[0][0], gauss,NDOF1); - for(i=0;ivalues[i]+=-Jdet*gauss->weight*(vx*dbdx+vy*dbdy-vz-meltingvalue)*L[i]; + for(i=0;ivalues[i]+=-Jdet*gauss->weight*(vx*dbdx+vy*dbdy-vz-basalmeltingvalue)*L[i]; } @@ -2389,5 +2389,5 @@ int i,j,ig; double Jdettria,dt; - double melting_g; + double basal_melting_g; double old_watercolumn_g; double xyz_list[NUMVERTICES][3]; @@ -2401,8 +2401,8 @@ GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES); this->parameters->FindParam(&dt,DtEnum); - Input* melting_input=inputs->GetInput(BasalMeltingRateEnum); _assert_(melting_input); + Input* basal_melting_input=inputs->GetInput(BasalMeltingRateEnum); _assert_(basal_melting_input); Input* old_watercolumn_input=inputs->GetInput(WaterColumnOldEnum);_assert_(old_watercolumn_input); - /*Initialize melting_correction_g to 0, do not forget!:*/ + /*Initialize basal_melting_correction_g to 0, do not forget!:*/ /* Start looping on the number of gaussian points: */ gauss=new GaussTria(2); @@ -2414,9 +2414,9 @@ GetL(&L[0], &xyz_list[0][0], gauss,NDOF1); - melting_input->GetParameterValue(&melting_g,gauss); + basal_melting_input->GetParameterValue(&basal_melting_g,gauss); old_watercolumn_input->GetParameterValue(&old_watercolumn_g,gauss); - if(dt)for(i=0;ivalues[i]+=Jdettria*gauss->weight*(old_watercolumn_g+dt*melting_g)*L[i]; - else for(i=0;ivalues[i]+=Jdettria*gauss->weight*melting_g*L[i]; + if(dt)for(i=0;ivalues[i]+=Jdettria*gauss->weight*(old_watercolumn_g+dt*basal_melting_g)*L[i]; + else for(i=0;ivalues[i]+=Jdettria*gauss->weight*basal_melting_g*L[i]; } @@ -2435,5 +2435,5 @@ int i,j,ig; double Jdettria,dt; - double accumulation_g,melting_g,melting_correction_g,thickness_g; + double surface_mass_balance_g,basal_melting_g,basal_melting_correction_g,thickness_g; double xyz_list[NUMVERTICES][3]; double L[NUMVERTICES]; @@ -2446,10 +2446,10 @@ GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES); this->parameters->FindParam(&dt,DtEnum); - Input* accumulation_input=inputs->GetInput(AccumulationRateEnum); _assert_(accumulation_input); - Input* melting_input=inputs->GetInput(BasalMeltingRateEnum); _assert_(melting_input); - Input* melting_correction_input=inputs->GetInput(BasalMeltingRateCorrectionEnum); + Input* surface_mass_balance_input=inputs->GetInput(SurfaceMassBalanceEnum); _assert_(surface_mass_balance_input); + Input* basal_melting_input=inputs->GetInput(BasalMeltingRateEnum); _assert_(basal_melting_input); + Input* basal_melting_correction_input=inputs->GetInput(BasalMeltingRateCorrectionEnum); Input* thickness_input=inputs->GetInput(ThicknessEnum); _assert_(thickness_input); - /*Initialize melting_correction_g to 0, do not forget!:*/ + /*Initialize basal_melting_correction_g to 0, do not forget!:*/ /* Start looping on the number of gaussian points: */ gauss=new GaussTria(2); @@ -2461,10 +2461,10 @@ GetL(&L[0], &xyz_list[0][0], gauss,NDOF1); - accumulation_input->GetParameterValue(&accumulation_g,gauss); - melting_input->GetParameterValue(&melting_g,gauss); + surface_mass_balance_input->GetParameterValue(&surface_mass_balance_g,gauss); + basal_melting_input->GetParameterValue(&basal_melting_g,gauss); thickness_input->GetParameterValue(&thickness_g,gauss); - if(melting_correction_input) melting_correction_input->GetParameterValue(&melting_correction_g,gauss); - - for(i=0;ivalues[i]+=Jdettria*gauss->weight*(thickness_g+dt*(accumulation_g-melting_g-melting_correction_g))*L[i]; + if(basal_melting_correction_input) basal_melting_correction_input->GetParameterValue(&basal_melting_correction_g,gauss); + + for(i=0;ivalues[i]+=Jdettria*gauss->weight*(thickness_g+dt*(surface_mass_balance_g-basal_melting_g-basal_melting_correction_g))*L[i]; } @@ -2483,5 +2483,5 @@ int i,j,ig; double Jdettria,dt; - double accumulation_g,melting_g,thickness_g; + double surface_mass_balance_g,basal_melting_g,thickness_g; double xyz_list[NUMVERTICES][3]; double L[NUMVERTICES]; @@ -2494,6 +2494,6 @@ this->parameters->FindParam(&dt,DtEnum); GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES); - Input* accumulation_input=inputs->GetInput(AccumulationRateEnum); _assert_(accumulation_input); - Input* melting_input=inputs->GetInput(BasalMeltingRateEnum); _assert_(melting_input); + Input* surface_mass_balance_input=inputs->GetInput(SurfaceMassBalanceEnum); _assert_(surface_mass_balance_input); + Input* basal_melting_input=inputs->GetInput(BasalMeltingRateEnum); _assert_(basal_melting_input); Input* thickness_input=inputs->GetInput(ThicknessEnum); _assert_(thickness_input); @@ -2507,9 +2507,9 @@ GetL(&L[0], &xyz_list[0][0], gauss,NDOF1); - accumulation_input->GetParameterValue(&accumulation_g,gauss); - melting_input->GetParameterValue(&melting_g,gauss); + surface_mass_balance_input->GetParameterValue(&surface_mass_balance_g,gauss); + basal_melting_input->GetParameterValue(&basal_melting_g,gauss); thickness_input->GetParameterValue(&thickness_g,gauss); - for(i=0;ivalues[i]+=Jdettria*gauss->weight*(thickness_g+dt*(accumulation_g-melting_g))*L[i]; + for(i=0;ivalues[i]+=Jdettria*gauss->weight*(thickness_g+dt*(surface_mass_balance_g-basal_melting_g))*L[i]; } @@ -3885,8 +3885,15 @@ this->inputs->AddInput(new TriaVertexInput(WaterColumnOldEnum,nodeinputs)); } - - if (iomodel->accumulation_rate) { - for(i=0;i<3;i++)nodeinputs[i]=iomodel->accumulation_rate[tria_vertex_ids[i]-1]/iomodel->yts; - this->inputs->AddInput(new TriaVertexInput(AccumulationRateEnum,nodeinputs)); + if (iomodel->surface_accumulation_rate) { + for(i=0;i<3;i++)nodeinputs[i]=iomodel->surface_accumulation_rate[tria_vertex_ids[i]-1]/iomodel->yts; + this->inputs->AddInput(new TriaVertexInput(SurfaceAccumulationRateEnum,nodeinputs)); + } + if (iomodel->surface_ablation_rate) { + for(i=0;i<3;i++)nodeinputs[i]=iomodel->surface_ablation_rate[tria_vertex_ids[i]-1]/iomodel->yts; + this->inputs->AddInput(new TriaVertexInput(SurfaceAblationRateEnum,nodeinputs)); + } + if (iomodel->surface_mass_balance) { + for(i=0;i<3;i++)nodeinputs[i]=iomodel->surface_mass_balance[tria_vertex_ids[i]-1]/iomodel->yts; + this->inputs->AddInput(new TriaVertexInput(SurfaceMassBalanceEnum,nodeinputs)); } if (iomodel->geothermalflux) { @@ -4530,5 +4537,5 @@ name==BasalMeltingRateEnum || name==WaterColumnEnum || - name==AccumulationRateEnum || + name==SurfaceMassBalanceEnum || name==SurfaceAreaEnum|| name==VxEnum || @@ -5146,5 +5153,5 @@ this->inputs->AddInput(new BoolInput(ElementOnIceShelfEnum,elementonshelf)); - /*If this element just became ungrounded, set its melting rate at 50 m/yr:*/ + /*If this element just became ungrounded, set its basal melting rate at 50 m/yr:*/ if(swap){ Input* basal_melting_rate_input =inputs->GetInput(BasalMeltingRateEnum); _assert_(basal_melting_rate_input); @@ -5886,6 +5893,6 @@ if(step==0){ - /*This is the first time we are trying to replace the accumulation tria vertex input by an accumulation triav vertex - *forcing, which is essentially a tria vertex input that can vary with time. So firt, if there is already + /*This is the first time we are trying to replace the smb tria vertex input by a smb triav vertex + *forcing, which is essentially a tria vertex input that can vary with time. So first, if there is already an input with enum FieldEnum, squash it: */ this->inputs->AddInput(new TriaVertexForcing(FieldEnum,nodeinputs,time,iomodel->forcing_numtimesteps,this->parameters)); Index: /issm/trunk/src/c/objects/IoModel.cpp =================================================================== --- /issm/trunk/src/c/objects/IoModel.cpp (revision 8398) +++ /issm/trunk/src/c/objects/IoModel.cpp (revision 8399) @@ -88,5 +88,7 @@ xfree((void**)&this->watercolumn); xfree((void**)&this->basal_melting_rate_correction); - xfree((void**)&this->accumulation_rate); + xfree((void**)&this->surface_accumulation_rate); + xfree((void**)&this->surface_ablation_rate); + xfree((void**)&this->surface_mass_balance); xfree((void**)&this->forcingtypes); xfree((void**)&this->forcing); @@ -398,6 +400,8 @@ this->penalties=NULL; - /*!basal: */ - this->accumulation_rate=NULL; + /*!surface: */ + this->surface_mass_balance=NULL; + this->surface_accumulation_rate=NULL; + this->surface_ablation_rate=NULL; this->dhdt=NULL; Index: /issm/trunk/src/c/objects/IoModel.h =================================================================== --- /issm/trunk/src/c/objects/IoModel.h (revision 8398) +++ /issm/trunk/src/c/objects/IoModel.h (revision 8399) @@ -194,5 +194,8 @@ double* basal_melting_rate_correction; int basal_melting_rate_correction_apply; - double* accumulation_rate; + double* surface_mass_balance; + double* surface_accumulation_rate; + double* surface_ablation_rate; + double* dhdt; Index: /issm/trunk/src/c/shared/Numerics/UnitConversion.cpp =================================================================== --- /issm/trunk/src/c/shared/Numerics/UnitConversion.cpp (revision 8398) +++ /issm/trunk/src/c/shared/Numerics/UnitConversion.cpp (revision 8399) @@ -68,6 +68,8 @@ case VelObsEnum: scale=yts;break; //m/yr case DhDtEnum: scale=yts;break; //m/yr - case BasalMeltingRateEnum: scale=yts;break; //m/yr - case AccumulationRateEnum: scale=yts;break; //m/yr + case BasalMeltingRateEnum: scale=yts;break; //m/yr + case SurfaceAccumulationRateEnum: scale=yts;break; //m/yr + case SurfaceAblationRateEnum: scale=yts;break; //m/yr + case SurfaceMassBalanceEnum: scale=yts;break; //m/yr case MisfitEnum: scale=pow(yts,2);break; //(m/yr)^2 case MassFluxEnum: scale=pow(10,-12)*yts;break; // (GigaTon/year) Index: /issm/trunk/src/c/solutions/transient2d_core.cpp =================================================================== --- /issm/trunk/src/c/solutions/transient2d_core.cpp (revision 8398) +++ /issm/trunk/src/c/solutions/transient2d_core.cpp (revision 8399) @@ -86,5 +86,5 @@ InputToResultx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,SurfaceEnum,step,time); InputToResultx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,BedEnum,step,time); - InputToResultx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,AccumulationRateEnum,step,time); + InputToResultx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,SurfaceMassBalanceEnum,step,time); InputToResultx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,BasalMeltingRateEnum,step,time); if(gl_migration!=NoneEnum)InputToResultx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,ElementOnIceShelfEnum,step,time); Index: /issm/trunk/src/m/classes/model.m =================================================================== --- /issm/trunk/src/m/classes/model.m (revision 8398) +++ /issm/trunk/src/m/classes/model.m (revision 8399) @@ -176,5 +176,7 @@ vel_bal=NaN; vel_obs_raw=NaN; - accumulation_rate=NaN; + surface_accumulation_rate=NaN; + surface_ablation_rate=NaN; + surface_mass_balance=NaN; dhdt=NaN; geothermalflux=NaN; @@ -373,4 +375,6 @@ n; melting; + melting_rate; + accumulation_rate; accumulation; type; @@ -426,5 +430,7 @@ if isfield(structmd,'n'), md.rheology_n=structmd.n; end if isfield(structmd,'melting'), md.basal_melting_rate=structmd.melting; end - if isfield(structmd,'accumulation'), md.accumulation_rate=structmd.accumulation; end + if isfield(structmd,'melting_rate'), md.basal_melting_rate=structmd.melting_rate; end + if isfield(structmd,'accumulation'), md.surface_mass_balance=structmd.accumulation; end + if isfield(structmd,'accumulation_rate'), md.surface_mass_balance=structmd.accumulation_rate; end if isfield(structmd,'numberofgrids'), md.numberofnodes=structmd.numberofgrids; end if isfield(structmd,'numberofgrids2d'), md.numberofnodes2d=structmd.numberofgrids2d; end Index: sm/trunk/src/m/enum/AccumulationRateEnum.m =================================================================== --- /issm/trunk/src/m/enum/AccumulationRateEnum.m (revision 8398) +++ (revision ) @@ -1,11 +1,0 @@ -function macro=AccumulationRateEnum() -%ACCUMULATIONRATEENUM - Enum of AccumulationRate -% -% WARNING: DO NOT MODIFY THIS FILE -% this file has been automatically generated by src/c/EnumDefinitions/Synchronize.sh -% Please read src/c/EnumDefinitions/README for more information -% -% Usage: -% macro=AccumulationRateEnum() - -macro=StringToEnum('AccumulationRate'); Index: /issm/trunk/src/m/enum/EnumToModelField.m =================================================================== --- /issm/trunk/src/m/enum/EnumToModelField.m (revision 8398) +++ /issm/trunk/src/m/enum/EnumToModelField.m (revision 8399) @@ -21,6 +21,8 @@ case VyObsEnum(), string='vy_obs'; return case GroundingLineMigrationEnum(), string='gl_migration'; return - case AccumulationRateEnum(), string='accumulation_rate'; return case BasalMeltingRateEnum(), string='basal_melting_rate'; return + case SurfaceAccumulationRateEnum(), string='surface_accumulation_rate'; return + case SurfaceAblationRateEnum(), string='surface_ablation_rate'; return + case SurfaceMassBalanceEnum(), string='surface_mass_balance'; return otherwise, error(['Enum ' num2str(enum) ' not found associated to any model field']); Index: /issm/trunk/src/m/enum/SurfaceAblationRateEnum.m =================================================================== --- /issm/trunk/src/m/enum/SurfaceAblationRateEnum.m (revision 8399) +++ /issm/trunk/src/m/enum/SurfaceAblationRateEnum.m (revision 8399) @@ -0,0 +1,11 @@ +function macro=SurfaceAblationRateEnum() +%SURFACEABLATIONRATEENUM - Enum of SurfaceAblationRate +% +% WARNING: DO NOT MODIFY THIS FILE +% this file has been automatically generated by src/c/EnumDefinitions/Synchronize.sh +% Please read src/c/EnumDefinitions/README for more information +% +% Usage: +% macro=SurfaceAblationRateEnum() + +macro=StringToEnum('SurfaceAblationRate'); Index: /issm/trunk/src/m/enum/SurfaceAccumulationRateEnum.m =================================================================== --- /issm/trunk/src/m/enum/SurfaceAccumulationRateEnum.m (revision 8399) +++ /issm/trunk/src/m/enum/SurfaceAccumulationRateEnum.m (revision 8399) @@ -0,0 +1,11 @@ +function macro=SurfaceAccumulationRateEnum() +%SURFACEACCUMULATIONRATEENUM - Enum of SurfaceAccumulationRate +% +% WARNING: DO NOT MODIFY THIS FILE +% this file has been automatically generated by src/c/EnumDefinitions/Synchronize.sh +% Please read src/c/EnumDefinitions/README for more information +% +% Usage: +% macro=SurfaceAccumulationRateEnum() + +macro=StringToEnum('SurfaceAccumulationRate'); Index: /issm/trunk/src/m/enum/SurfaceMassBalanceEnum.m =================================================================== --- /issm/trunk/src/m/enum/SurfaceMassBalanceEnum.m (revision 8399) +++ /issm/trunk/src/m/enum/SurfaceMassBalanceEnum.m (revision 8399) @@ -0,0 +1,11 @@ +function macro=SurfaceMassBalanceEnum() +%SURFACEMASSBALANCEENUM - Enum of SurfaceMassBalance +% +% WARNING: DO NOT MODIFY THIS FILE +% this file has been automatically generated by src/c/EnumDefinitions/Synchronize.sh +% Please read src/c/EnumDefinitions/README for more information +% +% Usage: +% macro=SurfaceMassBalanceEnum() + +macro=StringToEnum('SurfaceMassBalance'); Index: /issm/trunk/src/m/model/collapse.m =================================================================== --- /issm/trunk/src/m/model/collapse.m (revision 8398) +++ /issm/trunk/src/m/model/collapse.m (revision 8399) @@ -29,5 +29,5 @@ if ~isnan(md.vy_obs), md.vy_obs=project2d(md,md.vy_obs,md.numlayers); end; if ~isnan(md.vel_obs), md.vel_obs=project2d(md,md.vel_obs,md.numlayers); end; -if ~isnan(md.accumulation_rate), md.accumulation_rate=project2d(md,md.accumulation_rate,md.numlayers); end; +if ~isnan(md.surface_mass_balance), md.surface_mass_balance=project2d(md,md.surface_mass_balance,md.numlayers); end; if ~isnan(md.dhdt), md.dhdt=project2d(md,md.dhdt,md.numlayers); end; if ~isnan(md.firn_layer), md.firn_layer=project2d(md,md.firn_layer,md.numlayers); end; Index: /issm/trunk/src/m/model/display/displayobservations.m =================================================================== --- /issm/trunk/src/m/model/display/displayobservations.m (revision 8398) +++ /issm/trunk/src/m/model/display/displayobservations.m (revision 8399) @@ -17,5 +17,5 @@ fielddisplay(md,'vy_obs_raw','raw observed velocity y component [m/a]'); fielddisplay(md,'vel_obs_raw','raw observed magnitude [m/a]'); -fielddisplay(md,'accumulation_rate','surface accumulation rate [m/a]'); +fielddisplay(md,'surface_mass_balance','surface mass balance [m/a]'); fielddisplay(md,'dhdt','surface dhdt rate [m/a]'); fielddisplay(md,'observed_temperature','observed temperature [K]'); Index: /issm/trunk/src/m/model/extrude.m =================================================================== --- /issm/trunk/src/m/model/extrude.m (revision 8398) +++ /issm/trunk/src/m/model/extrude.m (revision 8399) @@ -150,5 +150,5 @@ md.vel_bal=project3d(md,md.vel_bal,'node'); md.vel_obs_raw=project3d(md,md.vel_obs_raw,'node'); -md.accumulation_rate=project3d(md,md.accumulation_rate,'node'); +md.surface_mass_balance=project3d(md,md.surface_mass_balance,'node'); md.dhdt=project3d(md,md.dhdt,'node'); md.firn_layer=project3d(md,md.firn_layer,'node',md.numlayers); Index: /issm/trunk/src/m/model/ismodelselfconsistent.m =================================================================== --- /issm/trunk/src/m/model/ismodelselfconsistent.m (revision 8398) +++ /issm/trunk/src/m/model/ismodelselfconsistent.m (revision 8399) @@ -92,5 +92,5 @@ %}}} %SIZE NUMBEROFNODES {{{1 -fields={'x','y','z','rheology_B','drag_coefficient','basal_melting_rate','accumulation_rate','surface','thickness','bed','nodeonbed','nodeonsurface'}; +fields={'x','y','z','rheology_B','drag_coefficient','basal_melting_rate','surface_mass_balance','surface','thickness','bed','nodeonbed','nodeonsurface'}; checksize(md,fields,[md.numberofnodes 1]); %}}} @@ -505,5 +505,5 @@ %INITIAL TEMPERATURE, MELTING AND ACCUMULATION - fields={'temperature','accumulation_rate','basal_melting_rate'}; + fields={'temperature','surface_mass_balance','basal_melting_rate'}; checksize(md,fields,[md.numberofnodes 1]); checknan(md,fields); @@ -521,5 +521,5 @@ % {{{2 %VELOCITIES MELTING AND ACCUMULATION - fields={'vx','vy','accumulation_rate','basal_melting_rate','dhdt'}; + fields={'vx','vy','surface_mass_balance','basal_melting_rate','dhdt'}; checksize(md,fields,[md.numberofnodes 1]); checknan(md,fields); @@ -540,5 +540,5 @@ % {{{2 %VELOCITIES MELTING AND ACCUMULATION - fields={'vx','vy','accumulation_rate','basal_melting_rate'}; + fields={'vx','vy','surface_mass_balance','basal_melting_rate'}; checksize(md,fields,[md.numberofnodes 1]); checknan(md,fields); Index: /issm/trunk/src/m/model/marshall.m =================================================================== --- /issm/trunk/src/m/model/marshall.m (revision 8398) +++ /issm/trunk/src/m/model/marshall.m (revision 8399) @@ -89,5 +89,7 @@ WriteData(fid,md.geothermalflux,'Mat','geothermalflux'); -WriteData(fid,md.accumulation_rate,'Mat','accumulation_rate'); +WriteData(fid,md.surface_accumulation_rate,'Mat','surface_accumulation_rate'); +WriteData(fid,md.surface_ablation_rate,'Mat','surface_ablation_rate'); +WriteData(fid,md.surface_mass_balance,'Mat','surface_mass_balance'); WriteData(fid,md.gl_melting_rate,'Scalar','gl_melting_rate'); WriteData(fid,md.basal_melting_rate,'Mat','basal_melting_rate'); @@ -101,6 +103,10 @@ for i=1:numforcings, switch (forcingnames{i}) - case 'accumulation_rate' - forcingtypes(i)=AccumulationRateEnum; + case 'surface_accumulation_rate' + forcingtypes(i)=SurfaceAccumulationRateEnum; + case 'surface_ablation_rate' + forcingtypes(i)=SurfaceAblationRateEnum; + case 'surface_mass_balance' + forcingtypes(i)=SurfaceMassBalanceEnum; case 'basal_melting_rate' forcingtypes(i)=BasalMeltingRateEnum; Index: /issm/trunk/src/m/model/modeldefault/defaultparams.m =================================================================== --- /issm/trunk/src/m/model/modeldefault/defaultparams.m (revision 8398) +++ /issm/trunk/src/m/model/modeldefault/defaultparams.m (revision 8399) @@ -88,7 +88,9 @@ %end - disp(' creating accumulation rates'); - md.accumulation_rate=ones(md.numberofnodes,1); %1m/a - + disp(' creating surface balance rates'); + md.surface_accumulation_rate=ones(md.numberofnodes,1); %1m/a + md.surface_ablation_rate=ones(md.numberofnodes,0); %0m/a + md.surface_mass_balance=ones(md.numberofnodes,1); %1m/a + %Deal with boundary conditions: Index: /issm/trunk/src/m/model/structtomodel.m =================================================================== --- /issm/trunk/src/m/model/structtomodel.m (revision 8398) +++ /issm/trunk/src/m/model/structtomodel.m (revision 8399) @@ -26,5 +26,7 @@ if isfield(structmd,'n'), md.rheology_n=structmd.n; end if isfield(structmd,'melting'), md.basal_melting_rate=structmd.melting; end -if isfield(structmd,'accumulation'), md.accumulation_rate=structmd.accumulation; end +if isfield(structmd,'melting_rate'), md.basal_melting_rate=structmd.melting_rate; end +if isfield(structmd,'accumulation'), md.surface_mass_balance=structmd.accumulation; end +if isfield(structmd,'accumulation_rate'), md.surface_mass_balance=structmd.accumulation_rate; end if isfield(structmd,'numberofgrids'), md.numberofnodes=structmd.numberofgrids; end if isfield(structmd,'numberofgrids2d'), md.numberofnodes2d=structmd.numberofgrids2d; end Index: /issm/trunk/src/m/utils/BC/SetIceSheetBC.m =================================================================== --- /issm/trunk/src/m/utils/BC/SetIceSheetBC.m (revision 8398) +++ /issm/trunk/src/m/utils/BC/SetIceSheetBC.m (revision 8399) @@ -28,9 +28,19 @@ end -%Create zeros basal_melting_rate and accumulation_rate if not specified -if isnan(md.accumulation_rate), - md.accumulation_rate=zeros(md.numberofnodes,1); - md.forcings.accumulation_rate=[zeros(md.numberofnodes+1,1)]; - disp(' no accumulation_rate specified: values set as zero'); +%Create zeros basal_melting_rate and surface mass balance if not specified +if isnan(md.surface_accumulation_rate), + md.surface_accumulation_rate=zeros(md.numberofnodes,1); + md.forcings.surface_accumulation_rate=zeros(md.numberofnodes+1,1); + disp(' no surface_accumulation_rate specified: values set as zero'); +end +if isnan(md.surface_ablation_rate), + md.surface_ablation_rate=zeros(md.numberofnodes,1); + md.forcings.surface_ablation_rate=zeros(md.numberofnodes+1,1); + disp(' no surface_ablation_rate specified: values set as zero'); +end +if isnan(md.surface_mass_balance), + md.surface_mass_balance=zeros(md.numberofnodes,1); + md.forcings.surface_mass_balance=zeros(md.numberofnodes+1,1); + disp(' no surface_mass_balance specified: values set as zero'); end if isnan(md.basal_melting_rate), Index: /issm/trunk/src/m/utils/BC/SetIceShelfBC.m =================================================================== --- /issm/trunk/src/m/utils/BC/SetIceShelfBC.m (revision 8398) +++ /issm/trunk/src/m/utils/BC/SetIceShelfBC.m (revision 8399) @@ -59,9 +59,20 @@ md.pressureload=pressureload; -%Create zeros basal_melting_rate and accumulation_rate if not specified -if isnan(md.accumulation_rate), - md.accumulation_rate=zeros(md.numberofnodes,1); - md.forcings.accumulation_rate=zeros(md.numberofnodes+1,1); - disp(' no accumulation_rate specified: values set as zero'); +%Create zeros basal_melting_rate, surface_ablation_rate, surface_accumulation_rate +% and surface_mass_balance if not specified +if isnan(md.surface_accumulation_rate), + md.surface_accumulation_rate=zeros(md.numberofnodes,1); + md.forcings.surface_accumulation_rate=zeros(md.numberofnodes+1,1); + disp(' no surface_accumulation_rate specified: values set as zero'); +end +if isnan(md.surface_ablation_rate), + md.surface_ablation_rate=zeros(md.numberofnodes,1); + md.forcings.surface_ablation_rate=zeros(md.numberofnodes+1,1); + disp(' no surface_ablation_rate specified: values set as zero'); +end +if isnan(md.surface_mass_balance), + md.surface_mass_balance=zeros(md.numberofnodes,1); + md.forcings.surface_mass_balance=zeros(md.numberofnodes+1,1); + disp(' no surface_mass_balance specified: values set as zero'); end if isnan(md.basal_melting_rate), Index: /issm/trunk/src/m/utils/BC/SetMarineIceSheetBC.m =================================================================== --- /issm/trunk/src/m/utils/BC/SetMarineIceSheetBC.m (revision 8398) +++ /issm/trunk/src/m/utils/BC/SetMarineIceSheetBC.m (revision 8399) @@ -69,9 +69,21 @@ md.pressureload=pressureload; -%Create zeros basal_melting_rate and accumulation_rate if not specified -if isnan(md.accumulation_rate), - md.accumulation_rate=zeros(md.numberofnodes,1); - md.forcings.accumulation_rate=zeros(md.numberofnodes+1,1); - disp(' no accumulation_rate specified: values set as zero'); + +%Create zeros basal_melting_rate, surface_ablation_rate, surface_accumulation_rate +% and surface_mass_balance if not specified +if isnan(md.surface_accumulation_rate), + md.surface_accumulation_rate=zeros(md.numberofnodes,1); + md.forcings.surface_accumulation_rate=zeros(md.numberofnodes+1,1); + disp(' no surface_accumulation_rate specified: values set as zero'); +end +if isnan(md.surface_ablation_rate), + md.surface_ablation_rate=zeros(md.numberofnodes,1); + md.forcings.surface_ablation_rate=zeros(md.numberofnodes+1,1); + disp(' no surface_ablation_rate specified: values set as zero'); +end +if isnan(md.surface_mass_balance), + md.surface_mass_balance=zeros(md.numberofnodes,1); + md.forcings.surface_mass_balance=zeros(md.numberofnodes+1,1); + disp(' no surface_mass_balance specified: values set as zero'); end if isnan(md.basal_melting_rate), Index: /issm/trunk/test/Par/79North.par =================================================================== --- /issm/trunk/test/Par/79North.par (revision 8398) +++ /issm/trunk/test/Par/79North.par (revision 8399) @@ -23,11 +23,11 @@ md.drag_q=ones(md.numberofelements,1); -%Ice shelf melting and accumulation +%Ice shelf melting and surface mass balance md.basal_melting_rate=zeros(md.numberofnodes,1); pos=zeros(md.numberofnodes,1); pos(md.elements(find(md.elementoniceshelf),:))=1; md.basal_melting_rate(find(pos))=10; -md.accumulation_rate=15*ones(md.numberofnodes,1); -md.forcings.accumulation_rate=[15*ones(md.numberofnodes,1);1]; +md.surface_mass_balance=15*ones(md.numberofnodes,1); +md.forcings.surface_mass_balance=[15*ones(md.numberofnodes,1);1]; %Numerical parameters Index: /issm/trunk/test/Par/RoundSheetEISMINT.par =================================================================== --- /issm/trunk/test/Par/RoundSheetEISMINT.par (revision 8398) +++ /issm/trunk/test/Par/RoundSheetEISMINT.par (revision 8399) @@ -23,10 +23,10 @@ md.rheology_n=3*ones(md.numberofelements,1); -disp(' creating accumulation rates'); -acc_max=0.5; %m/yr +disp(' creating surface mass balance'); +smb_max=0.5; %m/yr sb=10^-2/1000; %m/yr/m rel=450*1000; %m -md.accumulation_rate=min(acc_max,sb*(rel-radius)); -md.forcings.accumulation_rate=[min(acc_max,sb*(rel-radius));1]; +md.surface_mass_balance=min(smb_max,sb*(rel-radius)); +md.forcings.surface_mass_balance=[min(smb_max,sb*(rel-radius));1]; disp(' creating velocities'); Index: /issm/trunk/test/Par/RoundSheetShelf.par =================================================================== --- /issm/trunk/test/Par/RoundSheetShelf.par (revision 8398) +++ /issm/trunk/test/Par/RoundSheetShelf.par (revision 8399) @@ -45,7 +45,7 @@ md.temperature=md.observed_temperature; -%Accumulation and melting -md.accumulation_rate=-10*ones(md.numberofnodes,1); -md.forcings.accumulation_rate=[-10*ones(md.numberofnodes,1);1]; +%Surface mass balance and basal melting +md.surface_mass_balance=-10*ones(md.numberofnodes,1); +md.forcings.surface_mass_balance=[-10*ones(md.numberofnodes,1);1]; md.basal_melting_rate=zeros(md.numberofnodes,1); pos=find(md.nodeoniceshelf);md.basal_melting_rate(pos)=10; Index: /issm/trunk/test/Par/RoundSheetStaticEISMINT.par =================================================================== --- /issm/trunk/test/Par/RoundSheetStaticEISMINT.par (revision 8398) +++ /issm/trunk/test/Par/RoundSheetStaticEISMINT.par (revision 8399) @@ -29,10 +29,10 @@ md.rheology_n=3*ones(md.numberofelements,1); -disp(' creating accumulation rates'); -acc_max=0.5; %m/yr +disp(' creating surface mass balance'); +smb_max=0.5; %m/yr sb=10^-2/1000; %m/yr/m rel=450*1000; %m -md.accumulation_rate=min(acc_max,sb*(rel-radius)); -md.forcings.accumulation_rate=[min(acc_max,sb*(rel-radius));1]; +md.surface_mass_balance=min(smb_max,sb*(rel-radius)); +md.forcings.surface_mass_balance=[min(smb_max,sb*(rel-radius));1]; disp(' creating velocities'); Index: /issm/trunk/test/Par/SquareEISMINT.par =================================================================== --- /issm/trunk/test/Par/SquareEISMINT.par (revision 8398) +++ /issm/trunk/test/Par/SquareEISMINT.par (revision 8399) @@ -26,7 +26,7 @@ md.rheology_n=3*ones(md.numberofelements,1); -disp(' creating accumulation rates'); -md.accumulation_rate=0.2*ones(md.numberofnodes,1); %0m/a -md.forcings.accumulation_rate=[0.2*ones(md.numberofnodes,1);1]; %0m/a +disp(' creating surface mass balance'); +md.surface_mass_balance=0.2*ones(md.numberofnodes,1); %0m/a +md.forcings.surface_mass_balance=[0.2*ones(md.numberofnodes,1);1]; %0m/a md.basal_melting_rate=0*ones(md.numberofnodes,1); %0m/a Index: /issm/trunk/test/Par/SquareSheetShelf.par =================================================================== --- /issm/trunk/test/Par/SquareSheetShelf.par (revision 8398) +++ /issm/trunk/test/Par/SquareSheetShelf.par (revision 8399) @@ -28,6 +28,6 @@ %Accumulation and melting -md.accumulation_rate=10*ones(md.numberofnodes,1); -md.forcings.accumulation_rate=[10*ones(md.numberofnodes,1);1]; +md.surface_mass_balance=10*ones(md.numberofnodes,1); +md.forcings.surface_mass_balance=[10*ones(md.numberofnodes,1);1]; md.basal_melting_rate=5*ones(md.numberofnodes,1); Index: /issm/trunk/test/Par/SquareShelfConstrained.par =================================================================== --- /issm/trunk/test/Par/SquareShelfConstrained.par (revision 8398) +++ /issm/trunk/test/Par/SquareShelfConstrained.par (revision 8399) @@ -24,7 +24,7 @@ md.temperature=md.observed_temperature; -%Accumulation and melting -md.accumulation_rate=10*ones(md.numberofnodes,1); -md.forcings.accumulation_rate=[10*ones(md.numberofnodes,1);1]; +%Surface mass balance and basal melting +md.surface_mass_balance=10*ones(md.numberofnodes,1); +md.forcings.surface_mass_balance=[10*ones(md.numberofnodes,1);1]; md.basal_melting_rate=5*ones(md.numberofnodes,1); Index: /issm/trunk/test/Par/SquareThermal.par =================================================================== --- /issm/trunk/test/Par/SquareThermal.par (revision 8398) +++ /issm/trunk/test/Par/SquareThermal.par (revision 8399) @@ -33,7 +33,7 @@ md.rheology_n=3*ones(md.numberofelements,1); -disp(' creating accumulation rates'); -md.accumulation_rate=ones(md.numberofnodes,1)/md.yts; %1m/a -md.forcings.accumulation_rate=[ones(md.numberofnodes,1)/md.yts;1]; %1m/a +disp(' creating surface_mass_balance'); +md.surface_mass_balance=ones(md.numberofnodes,1)/md.yts; %1m/a +md.forcings.surface_mass_balance=[ones(md.numberofnodes,1)/md.yts;1]; %1m/a md.basal_melting_rate=0*ones(md.numberofnodes,1)/md.yts; %1m/a