source: issm/oecreview/Archive/18296-19100/ISSM-18815-18816.diff@ 19102

Last change on this file since 19102 was 19102, checked in by Mathieu Morlighem, 10 years ago

NEW: added 18296-19100

File size: 21.5 KB
RevLine 
[19102]1Index: ../trunk-jpl/test/NightlyRun/test804.py
2===================================================================
3--- ../trunk-jpl/test/NightlyRun/test804.py (revision 18815)
4+++ ../trunk-jpl/test/NightlyRun/test804.py (revision 18816)
5@@ -24,6 +24,7 @@
6 md.transient.iscalving=True;
7
8 md.calving.calvingrate=1000.*numpy.ones((md.mesh.numberofvertices,1))
9+md.calving.meltingrate=numpy.zeros((md.mesh.numberofvertices,1))
10
11 md=solve(md,TransientSolutionEnum())
12
13Index: ../trunk-jpl/test/NightlyRun/test804.m
14===================================================================
15--- ../trunk-jpl/test/NightlyRun/test804.m (revision 18815)
16+++ ../trunk-jpl/test/NightlyRun/test804.m (revision 18816)
17@@ -14,6 +14,7 @@
18 md.transient.iscalving=1;
19
20 md.calving.calvingrate=1000.*ones(md.mesh.numberofvertices,1);
21+md.calving.meltingrate=zeros(md.mesh.numberofvertices,1);
22
23 md=solve(md,TransientSolutionEnum());
24
25Index: ../trunk-jpl/test/NightlyRun/test805.py
26===================================================================
27--- ../trunk-jpl/test/NightlyRun/test805.py (revision 18815)
28+++ ../trunk-jpl/test/NightlyRun/test805.py (revision 18816)
29@@ -31,6 +31,7 @@
30 md.transient.iscalving=True;
31
32 md.calving.calvingrate=1000.*numpy.ones((md.mesh.numberofvertices,1))
33+md.calving.meltingrate=numpy.zeros((md.mesh.numberofvertices,1))
34
35 md=solve(md,TransientSolutionEnum())
36
37Index: ../trunk-jpl/test/NightlyRun/test805.m
38===================================================================
39--- ../trunk-jpl/test/NightlyRun/test805.m (revision 18815)
40+++ ../trunk-jpl/test/NightlyRun/test805.m (revision 18816)
41@@ -21,6 +21,7 @@
42 md.transient.iscalving=1;
43
44 md.calving.calvingrate=1000.*ones(md.mesh.numberofvertices,1);
45+md.calving.meltingrate=zeros(md.mesh.numberofvertices,1);
46
47 md=solve(md,TransientSolutionEnum());
48
49Index: ../trunk-jpl/test/NightlyRun/test806.py
50===================================================================
51--- ../trunk-jpl/test/NightlyRun/test806.py (revision 18815)
52+++ ../trunk-jpl/test/NightlyRun/test806.py (revision 18816)
53@@ -36,8 +36,9 @@
54 md.transient.isgia=False
55 md.transient.iscalving=True;
56
57-md.calving=calvinglevermann();
58-md.calving.coeff=4.89e13*numpy.ones((md.mesh.numberofvertices,1));
59+md.calving=calvinglevermann()
60+md.calving.coeff=4.89e13*numpy.ones((md.mesh.numberofvertices,1))
61+md.calving.meltingrate=numpy.zeros((md.mesh.numberofvertices,1))
62
63 md.transient.requested_outputs=['default','StrainRateparallel','StrainRateperpendicular','Calvingratex','Calvingratey','CalvingCalvingrate']
64
65Index: ../trunk-jpl/test/NightlyRun/test806.m
66===================================================================
67--- ../trunk-jpl/test/NightlyRun/test806.m (revision 18815)
68+++ ../trunk-jpl/test/NightlyRun/test806.m (revision 18816)
69@@ -25,6 +25,7 @@
70
71 md.calving=calvinglevermann();
72 md.calving.coeff=4.89e13*ones(md.mesh.numberofvertices,1);
73+md.calving.meltingrate=zeros(md.mesh.numberofvertices,1);
74
75 md.transient.requested_outputs={'default','StrainRateparallel','StrainRateperpendicular','Calvingratex','Calvingratey','CalvingCalvingrate'};
76
77Index: ../trunk-jpl/test/NightlyRun/IdToName.m
78===================================================================
79--- ../trunk-jpl/test/NightlyRun/IdToName.m (revision 18815)
80+++ ../trunk-jpl/test/NightlyRun/IdToName.m (revision 18816)
81@@ -194,7 +194,8 @@
82 case 803, name='ValleyGlacierLevelsetEnthalpyHO3d';
83 case 804, name='ValleyGlacierLevelsetCalvingSIA2d';
84 case 805, name='ValleyGlacierLevelsetEnthCalvingHO3d';
85- case 806, name='SquareShelfLevelsetCalvingSIA2dLevermann';
86+ case 806, name='SquareShelfLevelsetCalvingSSA2dLevermann';
87+ case 807, name='SquareShelfLevelsetMeltingSSA2d';
88 case 1101, name='ISMIPAHO';
89 case 1102, name='ISMIPAFS';
90 case 1103, name='ISMIPBHO';
91Index: ../trunk-jpl/test/NightlyRun/IdToName.py
92===================================================================
93--- ../trunk-jpl/test/NightlyRun/IdToName.py (revision 18815)
94+++ ../trunk-jpl/test/NightlyRun/IdToName.py (revision 18816)
95@@ -179,7 +179,8 @@
96 803 : 'ValleyGlacierLevelsetEnthalpyHO3d',
97 804 : 'ValleyGlacierLevelsetCalvingSIA2d',
98 805 : 'ValleyGlacierLevelsetEnthCalvingHO3d',
99- 806 : 'SquareShelfLevelsetCalvingSIA2dLevermann',
100+ 806 : 'SquareShelfLevelsetCalvingSSA2dLevermann',
101+ 807 : 'SquareShelfLevelsetMeltingSSA2d',
102 1101 : 'ISMIPAHO',
103 1102 : 'ISMIPAFS',
104 1103 : 'ISMIPBHO',
105Index: ../trunk-jpl/src/c/analyses/LevelsetAnalysis.cpp
106===================================================================
107--- ../trunk-jpl/src/c/analyses/LevelsetAnalysis.cpp (revision 18815)
108+++ ../trunk-jpl/src/c/analyses/LevelsetAnalysis.cpp (revision 18816)
109@@ -47,9 +47,11 @@
110 switch(calvinglaw){
111 case DefaultCalvingEnum:
112 iomodel->FetchDataToInput(elements,CalvingCalvingrateEnum);
113+ iomodel->FetchDataToInput(elements,CalvingMeltingrateEnum);
114 break;
115 case CalvingLevermannEnum:
116 iomodel->FetchDataToInput(elements,CalvinglevermannCoeffEnum);
117+ iomodel->FetchDataToInput(elements,CalvinglevermannMeltingrateEnum);
118 break;
119 default:
120 _error_("Calving law "<<EnumToStringx(calvinglaw)<<" not supported yet");
121@@ -124,7 +126,7 @@
122 IssmDouble Jdet, dt, D_scalar;
123 IssmDouble h,hx,hy,hz;
124 IssmDouble vel;
125- IssmDouble norm_dlsf, calvingrate;
126+ IssmDouble norm_dlsf, norm_calving, calvingrate, meltingrate;
127 IssmDouble* xyz_list = NULL;
128
129 /*Get problem dimension and whether there is calving or not*/
130@@ -150,6 +152,7 @@
131 IssmDouble* v = xNew<IssmDouble>(dim);
132 IssmDouble* w = xNew<IssmDouble>(dim);
133 IssmDouble* c = xNewZeroInit<IssmDouble>(dim);
134+ IssmDouble* m = xNewZeroInit<IssmDouble>(dim);
135 IssmDouble* dlsf = xNew<IssmDouble>(dim);
136
137 /*Retrieve all inputs and parameters*/
138@@ -162,6 +165,7 @@
139 Input* lsf_slopex_input = NULL;
140 Input* lsf_slopey_input = NULL;
141 Input* calvingrate_input = NULL;
142+ Input* meltingrate_input = NULL;
143
144 /*Load velocities*/
145 if(domaintype==Domain2DhorizontalEnum){
146@@ -185,6 +189,7 @@
147 lsf_slopex_input = basalelement->GetInput(LevelsetfunctionSlopeXEnum); _assert_(lsf_slopex_input);
148 if(dim==2) lsf_slopey_input = basalelement->GetInput(LevelsetfunctionSlopeYEnum); _assert_(lsf_slopey_input);
149 calvingrate_input = basalelement->GetInput(CalvingCalvingrateEnum); _assert_(calvingrate_input);
150+ meltingrate_input = basalelement->GetInput(CalvingMeltingrateEnum); _assert_(meltingrate_input);
151 break;
152 case CalvingLevermannEnum:
153 if(domaintype==Domain2DhorizontalEnum){
154@@ -200,6 +205,7 @@
155 calvingratey_input=basalelement->GetInput(CalvingrateyAverageEnum); _assert_(calvingratey_input);
156 }
157 }
158+ meltingrate_input = basalelement->GetInput(CalvinglevermannMeltingrateEnum); _assert_(meltingrate_input);
159 break;
160 default:
161 _error_("Calving law "<<EnumToStringx(calvinglaw)<<" not supported yet");
162@@ -237,19 +243,35 @@
163 lsf_slopex_input->GetInputValue(&dlsf[0],gauss);
164 if(dim==2) lsf_slopey_input->GetInputValue(&dlsf[1],gauss);
165 calvingrate_input->GetInputValue(&calvingrate,gauss);
166+ meltingrate_input->GetInputValue(&meltingrate,gauss);
167
168 norm_dlsf=0.;
169 for(i=0;i<dim;i++) norm_dlsf+=pow(dlsf[i],2);
170 norm_dlsf=sqrt(norm_dlsf);
171
172 if(norm_dlsf>1.e-10)
173- for(i=0;i<dim;i++) c[i]=calvingrate*dlsf[i]/norm_dlsf;
174+ for(i=0;i<dim;i++){
175+ c[i]=calvingrate*dlsf[i]/norm_dlsf;
176+ m[i]=meltingrate*dlsf[i]/norm_dlsf;
177+ }
178 else
179- for(i=0;i<dim;i++) c[i]=0.;
180+ for(i=0;i<dim;i++){
181+ c[i]=0.;
182+ m[i]=0.;
183+ }
184+
185 break;
186 case CalvingLevermannEnum:
187 calvingratex_input->GetInputValue(&c[0],gauss); // in 3D case, add mesh velocity
188 if(dim==2) calvingratey_input->GetInputValue(&c[1],gauss);
189+ meltingrate_input->GetInputValue(&meltingrate,gauss);
190+
191+ norm_calving=0.;
192+ for(i=0;i<dim;i++) norm_calving+=pow(c[i],2);
193+ norm_calving=sqrt(norm_calving)+1.e-14;
194+
195+ for(i=0;i<dim;i++) m[i]=meltingrate*c[i]/norm_calving;
196+
197 break;
198 default:
199 _error_("Calving law "<<EnumToStringx(calvinglaw)<<" not supported yet");
200@@ -257,7 +279,7 @@
201 }
202
203 /*Levelset speed is ice velocity - calving rate*/
204- for(i=0;i<dim;i++) w[i]=v[i]-c[i];
205+ for(i=0;i<dim;i++) w[i]=v[i]-c[i]-m[i];
206
207 /*Compute D*/
208 for(row=0;row<dim;row++){
209Index: ../trunk-jpl/src/c/shared/Enum/EnumDefinitions.h
210===================================================================
211--- ../trunk-jpl/src/c/shared/Enum/EnumDefinitions.h (revision 18815)
212+++ ../trunk-jpl/src/c/shared/Enum/EnumDefinitions.h (revision 18816)
213@@ -212,10 +212,12 @@
214 StressIntensityFactorEnum,
215 CalvingLawEnum,
216 CalvingCalvingrateEnum,
217+ CalvingMeltingrateEnum,
218 CalvingLevermannEnum,
219 DefaultCalvingEnum,
220 CalvingRequestedOutputsEnum,
221 CalvinglevermannCoeffEnum,
222+ CalvinglevermannMeltingrateEnum,
223 CalvingratexEnum,
224 CalvingrateyEnum,
225 CalvingratexAverageEnum,
226Index: ../trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp
227===================================================================
228--- ../trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp (revision 18815)
229+++ ../trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp (revision 18816)
230@@ -220,10 +220,12 @@
231 case StressIntensityFactorEnum : return "StressIntensityFactor";
232 case CalvingLawEnum : return "CalvingLaw";
233 case CalvingCalvingrateEnum : return "CalvingCalvingrate";
234+ case CalvingMeltingrateEnum : return "CalvingMeltingrate";
235 case CalvingLevermannEnum : return "CalvingLevermann";
236 case DefaultCalvingEnum : return "DefaultCalving";
237 case CalvingRequestedOutputsEnum : return "CalvingRequestedOutputs";
238 case CalvinglevermannCoeffEnum : return "CalvinglevermannCoeff";
239+ case CalvinglevermannMeltingrateEnum : return "CalvinglevermannMeltingrate";
240 case CalvingratexEnum : return "Calvingratex";
241 case CalvingrateyEnum : return "Calvingratey";
242 case CalvingratexAverageEnum : return "CalvingratexAverage";
243Index: ../trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp
244===================================================================
245--- ../trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp (revision 18815)
246+++ ../trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp (revision 18816)
247@@ -223,10 +223,12 @@
248 else if (strcmp(name,"StressIntensityFactor")==0) return StressIntensityFactorEnum;
249 else if (strcmp(name,"CalvingLaw")==0) return CalvingLawEnum;
250 else if (strcmp(name,"CalvingCalvingrate")==0) return CalvingCalvingrateEnum;
251+ else if (strcmp(name,"CalvingMeltingrate")==0) return CalvingMeltingrateEnum;
252 else if (strcmp(name,"CalvingLevermann")==0) return CalvingLevermannEnum;
253 else if (strcmp(name,"DefaultCalving")==0) return DefaultCalvingEnum;
254 else if (strcmp(name,"CalvingRequestedOutputs")==0) return CalvingRequestedOutputsEnum;
255 else if (strcmp(name,"CalvinglevermannCoeff")==0) return CalvinglevermannCoeffEnum;
256+ else if (strcmp(name,"CalvinglevermannMeltingrate")==0) return CalvinglevermannMeltingrateEnum;
257 else if (strcmp(name,"Calvingratex")==0) return CalvingratexEnum;
258 else if (strcmp(name,"Calvingratey")==0) return CalvingrateyEnum;
259 else if (strcmp(name,"CalvingratexAverage")==0) return CalvingratexAverageEnum;
260@@ -257,12 +259,12 @@
261 else if (strcmp(name,"MeshVertexonbase")==0) return MeshVertexonbaseEnum;
262 else if (strcmp(name,"MeshVertexonsurface")==0) return MeshVertexonsurfaceEnum;
263 else if (strcmp(name,"MeshVertexonboundary")==0) return MeshVertexonboundaryEnum;
264- else if (strcmp(name,"MeshX")==0) return MeshXEnum;
265- else if (strcmp(name,"MeshY")==0) return MeshYEnum;
266 else stage=3;
267 }
268 if(stage==3){
269- if (strcmp(name,"MeshZ")==0) return MeshZEnum;
270+ if (strcmp(name,"MeshX")==0) return MeshXEnum;
271+ else if (strcmp(name,"MeshY")==0) return MeshYEnum;
272+ else if (strcmp(name,"MeshZ")==0) return MeshZEnum;
273 else if (strcmp(name,"MeshElementtype")==0) return MeshElementtypeEnum;
274 else if (strcmp(name,"DomainType")==0) return DomainTypeEnum;
275 else if (strcmp(name,"DomainDimension")==0) return DomainDimensionEnum;
276@@ -380,12 +382,12 @@
277 else if (strcmp(name,"SurfaceforcingsIssmbgradients")==0) return SurfaceforcingsIssmbgradientsEnum;
278 else if (strcmp(name,"SolutionType")==0) return SolutionTypeEnum;
279 else if (strcmp(name,"AnalysisType")==0) return AnalysisTypeEnum;
280- else if (strcmp(name,"ConfigurationType")==0) return ConfigurationTypeEnum;
281- else if (strcmp(name,"AdjointBalancethicknessAnalysis")==0) return AdjointBalancethicknessAnalysisEnum;
282 else stage=4;
283 }
284 if(stage==4){
285- if (strcmp(name,"AdjointBalancethickness2Analysis")==0) return AdjointBalancethickness2AnalysisEnum;
286+ if (strcmp(name,"ConfigurationType")==0) return ConfigurationTypeEnum;
287+ else if (strcmp(name,"AdjointBalancethicknessAnalysis")==0) return AdjointBalancethicknessAnalysisEnum;
288+ else if (strcmp(name,"AdjointBalancethickness2Analysis")==0) return AdjointBalancethickness2AnalysisEnum;
289 else if (strcmp(name,"AdjointHorizAnalysis")==0) return AdjointHorizAnalysisEnum;
290 else if (strcmp(name,"AnalysisCounter")==0) return AnalysisCounterEnum;
291 else if (strcmp(name,"DefaultAnalysis")==0) return DefaultAnalysisEnum;
292@@ -503,12 +505,12 @@
293 else if (strcmp(name,"MassconaxpbyName")==0) return MassconaxpbyNameEnum;
294 else if (strcmp(name,"MassconaxpbyNamex")==0) return MassconaxpbyNamexEnum;
295 else if (strcmp(name,"MassconaxpbyNamey")==0) return MassconaxpbyNameyEnum;
296- else if (strcmp(name,"MassconaxpbyAlpha")==0) return MassconaxpbyAlphaEnum;
297- else if (strcmp(name,"MassconaxpbyBeta")==0) return MassconaxpbyBetaEnum;
298 else stage=5;
299 }
300 if(stage==5){
301- if (strcmp(name,"NodeSId")==0) return NodeSIdEnum;
302+ if (strcmp(name,"MassconaxpbyAlpha")==0) return MassconaxpbyAlphaEnum;
303+ else if (strcmp(name,"MassconaxpbyBeta")==0) return MassconaxpbyBetaEnum;
304+ else if (strcmp(name,"NodeSId")==0) return NodeSIdEnum;
305 else if (strcmp(name,"VectorParam")==0) return VectorParamEnum;
306 else if (strcmp(name,"Riftfront")==0) return RiftfrontEnum;
307 else if (strcmp(name,"RiftfrontType")==0) return RiftfrontTypeEnum;
308@@ -626,12 +628,12 @@
309 else if (strcmp(name,"StrainRateyy")==0) return StrainRateyyEnum;
310 else if (strcmp(name,"StrainRateyz")==0) return StrainRateyzEnum;
311 else if (strcmp(name,"StrainRatezz")==0) return StrainRatezzEnum;
312- else if (strcmp(name,"Divergence")==0) return DivergenceEnum;
313- else if (strcmp(name,"MaxDivergence")==0) return MaxDivergenceEnum;
314 else stage=6;
315 }
316 if(stage==6){
317- if (strcmp(name,"GiaCrossSectionShape")==0) return GiaCrossSectionShapeEnum;
318+ if (strcmp(name,"Divergence")==0) return DivergenceEnum;
319+ else if (strcmp(name,"MaxDivergence")==0) return MaxDivergenceEnum;
320+ else if (strcmp(name,"GiaCrossSectionShape")==0) return GiaCrossSectionShapeEnum;
321 else if (strcmp(name,"GiadWdt")==0) return GiadWdtEnum;
322 else if (strcmp(name,"GiaW")==0) return GiaWEnum;
323 else if (strcmp(name,"P0")==0) return P0Enum;
324@@ -749,12 +751,12 @@
325 else if (strcmp(name,"Dense")==0) return DenseEnum;
326 else if (strcmp(name,"MpiDense")==0) return MpiDenseEnum;
327 else if (strcmp(name,"MpiSparse")==0) return MpiSparseEnum;
328- else if (strcmp(name,"Seq")==0) return SeqEnum;
329- else if (strcmp(name,"Mpi")==0) return MpiEnum;
330 else stage=7;
331 }
332 if(stage==7){
333- if (strcmp(name,"Mumps")==0) return MumpsEnum;
334+ if (strcmp(name,"Seq")==0) return SeqEnum;
335+ else if (strcmp(name,"Mpi")==0) return MpiEnum;
336+ else if (strcmp(name,"Mumps")==0) return MumpsEnum;
337 else if (strcmp(name,"Gsl")==0) return GslEnum;
338 else if (strcmp(name,"Option")==0) return OptionEnum;
339 else if (strcmp(name,"GenericOption")==0) return GenericOptionEnum;
340Index: ../trunk-jpl/src/m/classes/calving.py
341===================================================================
342--- ../trunk-jpl/src/m/classes/calving.py (revision 18815)
343+++ ../trunk-jpl/src/m/classes/calving.py (revision 18816)
344@@ -14,6 +14,7 @@
345
346 def __init__(self): # {{{
347 self.calvingrate = float('NaN')
348+ self.meltingrate = float('NaN')
349
350 #set defaults
351 self.setdefaultparameters()
352@@ -22,6 +23,7 @@
353 def __repr__(self): # {{{
354 string=' Calving parameters:'
355 string="%s\n%s"%(string,fielddisplay(self,'calvingrate','calving rate at given location [m/a]'))
356+ string="%s\n%s"%(string,fielddisplay(self,'meltingrate','melting rate at given location [m/a]'))
357
358 return string
359 #}}}
360@@ -36,6 +38,7 @@
361 return md
362
363 md = checkfield(md,'fieldname','calving.calvingrate','NaN',1,'size',[md.mesh.numberofvertices],'>=',0)
364+ md = checkfield(md,'fieldname','calving.meltingrate','NaN',1,'size',[md.mesh.numberofvertices],'>=',0)
365 return md
366 # }}}
367 def marshall(self,md,fid): # {{{
368@@ -44,4 +47,5 @@
369
370 WriteData(fid,'enum',CalvingLawEnum(),'data',DefaultCalvingEnum(),'format','Integer');
371 WriteData(fid,'object',self,'fieldname','calvingrate','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1,'scale',1./yts)
372+ WriteData(fid,'object',self,'fieldname','meltingrate','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1,'scale',1./yts)
373 # }}}
374Index: ../trunk-jpl/src/m/classes/calvinglevermann.py
375===================================================================
376--- ../trunk-jpl/src/m/classes/calvinglevermann.py (revision 18815)
377+++ ../trunk-jpl/src/m/classes/calvinglevermann.py (revision 18816)
378@@ -13,7 +13,8 @@
379 """
380
381 def __init__(self): # {{{
382- self.coeff = float('NaN')
383+ self.coeff = float('NaN')
384+ self.meltingrate = float('NaN')
385
386 #set defaults
387 self.setdefaultparameters()
388@@ -22,6 +23,7 @@
389 def __repr__(self): # {{{
390 string=' Calving Levermann parameters:'
391 string="%s\n%s"%(string,fielddisplay(self,'coeff','proportionality coefficient in Levermann model'))
392+ string="%s\n%s"%(string,fielddisplay(self,'meltingrate','melting rate at given location [m/a]'))
393
394 return string
395 #}}}
396@@ -37,9 +39,12 @@
397 return md
398
399 md = checkfield(md,'fieldname','calving.coeff','size',[md.mesh.numberofvertices],'>',0)
400+ md = checkfield(md,'fieldname','calving.meltingrate','NaN',1,'size',[md.mesh.numberofvertices],'>=',0)
401 return md
402 # }}}
403 def marshall(self,md,fid): # {{{
404+ yts=365.*24.*3600.
405 WriteData(fid,'enum',CalvingLawEnum(),'data',CalvingLevermannEnum(),'format','Integer');
406 WriteData(fid,'enum',CalvinglevermannCoeffEnum(),'data',self.coeff,'format','DoubleMat','mattype',1)
407+ WriteData(fid,'object',self,'fieldname','meltingrate','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1,'scale',1./yts)
408 # }}}
409Index: ../trunk-jpl/src/m/classes/calving.m
410===================================================================
411--- ../trunk-jpl/src/m/classes/calving.m (revision 18815)
412+++ ../trunk-jpl/src/m/classes/calving.m (revision 18816)
413@@ -6,6 +6,7 @@
414 classdef calving
415 properties (SetAccess=public)
416 calvingrate = NaN;
417+ meltingrate = NaN;
418 end
419 methods (Static)
420 function obj = loadobj(obj) % {{{
421@@ -59,15 +60,18 @@
422 %Early return
423 if (solution~=TransientSolutionEnum() | md.transient.iscalving==0), return; end
424 md = checkfield(md,'fieldname','calving.calvingrate','NaN',1,'size',[md.mesh.numberofvertices 1],'>=',0);
425+ md = checkfield(md,'fieldname','calving.meltingrate','NaN',1,'size',[md.mesh.numberofvertices 1],'>=',0);
426 end % }}}
427 function disp(obj) % {{{
428 disp(sprintf(' Calving parameters:'));
429 fielddisplay(obj,'calvingrate','calving rate at given location [m/a]');
430+ fielddisplay(obj,'meltingrate','melting rate at given location [m/a]');
431 end % }}}
432 function marshall(obj,md,fid) % {{{
433 yts=365.0*24.0*3600.0;
434 WriteData(fid,'enum',CalvingLawEnum(),'data',DefaultCalvingEnum(),'format','Integer');
435 WriteData(fid,'object',obj,'fieldname','calvingrate','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1,'scale',1./yts);
436+ WriteData(fid,'object',obj,'fieldname','meltingrate','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1,'scale',1./yts);
437 end % }}}
438 end
439 end
440Index: ../trunk-jpl/src/m/classes/calvinglevermann.m
441===================================================================
442--- ../trunk-jpl/src/m/classes/calvinglevermann.m (revision 18815)
443+++ ../trunk-jpl/src/m/classes/calvinglevermann.m (revision 18816)
444@@ -5,7 +5,8 @@
445
446 classdef calvinglevermann
447 properties (SetAccess=public)
448- coeff = NaN;
449+ coeff = NaN;
450+ meltingrate = NaN;
451 end
452 methods
453 function obj = calvinglevermann(varargin) % {{{
454@@ -34,15 +35,19 @@
455 %Early return
456 if (solution~=TransientSolutionEnum() | md.transient.iscalving==0), return; end
457 md = checkfield(md,'fieldname','calving.coeff','>',0,'size',[md.mesh.numberofvertices 1]);
458+ md = checkfield(md,'fieldname','calving.meltingrate','NaN',1,'size',[md.mesh.numberofvertices 1],'>=',0);
459 end % }}}
460 function disp(obj) % {{{
461 disp(sprintf(' Calving Levermann parameters:'));
462 fielddisplay(obj,'coeff','proportionality coefficient in Levermann model');
463+ fielddisplay(obj,'meltingrate','melting rate at given location [m/a]');
464
465 end % }}}
466 function marshall(obj,md,fid) % {{{
467+ yts=365.0*24.0*3600.0;
468 WriteData(fid,'enum',CalvingLawEnum(),'data',CalvingLevermannEnum(),'format','Integer');
469 WriteData(fid,'enum',CalvinglevermannCoeffEnum(),'data',obj.coeff,'format','DoubleMat','mattype',1);
470+ WriteData(fid,'object',obj,'fieldname','meltingrate','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1,'scale',1./yts);
471 end % }}}
472 end
473 end
Note: See TracBrowser for help on using the repository browser.