source: issm/oecreview/Archive/25834-26739/ISSM-26664-26665.diff@ 26740

Last change on this file since 26740 was 26740, checked in by Mathieu Morlighem, 3 years ago

CHG: added 25834-26739

File size: 15.7 KB
RevLine 
[26740]1Index: ../trunk-jpl/test/Archives/Archive543.arch
2===================================================================
3Cannot display: file marked as a binary type.
4svn:mime-type = application/octet-stream
5Index: ../trunk-jpl/test/NightlyRun/test543.m
6===================================================================
7--- ../trunk-jpl/test/NightlyRun/test543.m (revision 26664)
8+++ ../trunk-jpl/test/NightlyRun/test543.m (revision 26665)
9@@ -1,4 +1,4 @@
10-%Test Name: PigTranFrontalforcingsrignotautoregression
11+%Test Name: PigTranStochasticforcing
12 md=triangle(model(),'../Exp/Pig.exp',10000.);
13 md=setmask(md,'../Exp/PigShelves.exp','../Exp/PigIslands.exp');
14 md=parameterize(md,'../Par/Pig.par');
15@@ -5,28 +5,48 @@
16 md=setflowequation(md,'SSA','all');
17 md.timestepping.start_time = 0;
18 md.timestepping.time_step = 1;
19-md.timestepping.final_time = 25;
20+md.timestepping.final_time = 10;
21
22-%separate domain in 2 basins
23-idbasin = zeros(md.mesh.numberofelements,1);
24+%Basin separation TF
25+idb_tf = zeros(md.mesh.numberofelements,1);
26 iid1 = find(md.mesh.x<=-1.6e6);
27 for ii=1:md.mesh.numberofelements
28 for vertex=1:3
29 if any(iid1==md.mesh.elements(ii,vertex)) %one vertex in basin 1
30- idbasin(ii) = 1;
31+ idb_tf(ii) = 1;
32 end
33 end
34- if idbasin(ii)==0 %no vertex was found in basin 1
35- idbasin(ii) = 2;
36+ if idb_tf(ii)==0 %no vertex was found in basin 1
37+ idb_tf(ii) = 2;
38 end
39 end
40+% Basin separation default
41+idb_df = zeros(md.mesh.numberofelements,1);
42+iid1 = find(md.mesh.x<=-1.62e6);
43+for ii=1:md.mesh.numberofelements
44+ for vertex=1:3
45+ if any(iid1==md.mesh.elements(ii,vertex)) %one vertex in basin 1
46+ idb_df(ii) = 1;
47+ end
48+ end
49+ if idb_df(ii)==0 %no vertex was found in basin 1
50+ idb_df(ii) = 2;
51+ end
52+end
53+% Dimensionalities
54+nb_tf = 2;
55+nb_clv = 2;
56+nb_flmlt = 2;
57
58-%calving and frontalforcings parameters
59+%Calving parameters
60 md.mask.ice_levelset = 1e4*(md.mask.ice_levelset + 0.5);
61-md.calving=calvingvonmises();
62+md.calving.calvingrate = 0.3*ones(md.mesh.numberofvertices,1);
63+md.levelset.spclevelset = NaN(md.mesh.numberofvertices,1);
64+md.levelset.migration_max = 10.0; %avoid fast advance/retreat of the front
65+%Frontal forcing parameters
66 md.frontalforcings=frontalforcingsrignotautoregression();
67-md.frontalforcings.num_basins = 2;
68-md.frontalforcings.basin_id = idbasin;
69+md.frontalforcings.num_basins = nb_tf;
70+md.frontalforcings.basin_id = idb_tf;
71 md.frontalforcings.subglacial_discharge = 0.1*ones(md.mesh.numberofvertices,1);
72 md.frontalforcings.beta0 = [0.05,0.01]; %intercept values of TF in basins [C]
73 md.frontalforcings.beta1 = [0.001,0.0001]; %trend values of TF in basins [C/yr]
74@@ -35,57 +55,65 @@
75 md.frontalforcings.ar_timestep = 2; %timestep of the autoregressive model [yr]
76 md.frontalforcings.phi = [[0.1,-0.1,0.01,-0.01];[0.2,-0.2,0.1,0.0]]; %autoregressive parameters
77
78+% Floating Ice Melt parameters
79+md.basalforcings.floatingice_melting_rate = 0.1*ones(md.mesh.numberofvertices,1);
80+
81+% Covariance matrix
82+covtf = 1e-4*eye(nb_tf);
83+covclv = 1e-1*eye(nb_clv);
84+covclv(1,1) = 1/10*covclv(1,1);
85+covflmlt = 0.05*eye(nb_flmlt);
86+covglob = blkdiag(covtf,covclv,covflmlt);
87+
88 %Stochastic forcing
89 md.stochasticforcing.isstochasticforcing = 1;
90-md.stochasticforcing.fields = [{'FrontalForcingsRignotAutoregression'}];
91-md.stochasticforcing.covariance = 1e-4*[[1.5,0.5];[0.5,0.4]]; %global covariance among- and between-fields
92+md.stochasticforcing.fields = [{'FrontalForcingsRignotAutoregression'},{'DefaultCalving'},{'FloatingMeltRate'}];
93+md.stochasticforcing.defaultdimension = 2;
94+md.stochasticforcing.default_id = idb_df;
95+md.stochasticforcing.covariance = covglob; %global covariance among- and between-fields
96 md.stochasticforcing.randomflag = 0; %determines true/false randomness
97
98-md.transient.ismovingfront = 1;
99-md.levelset.spclevelset = NaN(md.mesh.numberofvertices,1);
100-md.levelset.migration_max = 1e10;
101+md.transient.ismovingfront = 1;
102+md.transient.isgroundingline = 1;
103
104-md.transient.requested_outputs = {'default','CalvingMeltingrate'};
105+md.transient.requested_outputs = {'default','CalvingCalvingrate','CalvingMeltingrate','BasalforcingsFloatingiceMeltingRate'};
106 md.cluster=generic('name',oshostname(),'np',2);
107 md=solve(md,'Transient');
108
109 %Fields and tolerances to track changes
110 field_names ={...
111- 'Vx1' ,'Vy1' ,'Vel1' ,'Pressure1' ,'Bed1' ,'Surface1' ,'Thickness1' ,'MaskIceLevelset1', 'CalvingMeltingrate1',...
112- 'Vx2' ,'Vy2' ,'Vel2' ,'Pressure2' ,'Bed2' ,'Surface2' ,'Thickness2' ,'MaskIceLevelset2', 'CalvingMeltingrate2',...
113- 'Vx10','Vy10','Vel10','Pressure10','Bed10','Surface10','Thickness10','MaskIceLevelset10', 'CalvingMeltingrate10',...
114+ 'Vx1' ,'Vy1' ,'Vel1' ,'Thickness1' ,'MaskIceLevelset1' ,'CalvingCalvingrate1' ,'CalvingMeltingrate1' ,'BasalforcingsFloatingiceMeltingRate1',...
115+ 'Vx2' ,'Vy2' ,'Vel2' ,'Thickness2' ,'MaskIceLevelset2' ,'CalvingCalvingrate1' ,'CalvingMeltingrate2' ,'BasalforcingsFloatingiceMeltingRate2',...
116+ 'Vx10','Vy10','Vel10','Thickness10','MaskIceLevelset10','CalvingCalvingrate10','CalvingMeltingrate10','BasalforcingsFloatingiceMeltingRate10',...
117 };
118 field_tolerances={...
119- 1e-11,2e-11,2e-11,1e-12,2e-11,6e-12,9e-12,1e-12,1e-9,...
120- 2e-11,1e-11,1e-11,9e-12,2e-1,2e-11,2e-11,1e-11,1e-9,...
121- 2e-10,1e-10,1e-10,1e-10,1e-10,1e-10,1e-10,1e-9,1e-9,...
122+ 1e-11,2e-11,2e-11,1e-11,1e-9,1e-10,1e-10,1e-10,...
123+ 2e-11,1e-11,1e-11,9e-11,2e-9,1e-10,1e-10,1e-10,...
124+ 2e-10,1e-10,1e-10,1e-10,5e-9,1e-10,1e-10,1e-10,...
125 };
126 field_values={...
127 (md.results.TransientSolution(1).Vx),...
128 (md.results.TransientSolution(1).Vy),...
129 (md.results.TransientSolution(1).Vel),...
130- (md.results.TransientSolution(1).Pressure),...
131- (md.results.TransientSolution(1).Base),...
132- (md.results.TransientSolution(1).Surface),...
133 (md.results.TransientSolution(1).Thickness),...
134 (md.results.TransientSolution(1).MaskIceLevelset),...
135+ (md.results.TransientSolution(1).CalvingCalvingrate),...
136 (md.results.TransientSolution(1).CalvingMeltingrate),...
137- (md.results.TransientSolution(2).Vx),...
138- (md.results.TransientSolution(2).Vy),...
139- (md.results.TransientSolution(2).Vel),...
140- (md.results.TransientSolution(2).Pressure),...
141- (md.results.TransientSolution(2).Base),...
142- (md.results.TransientSolution(2).Surface),...
143- (md.results.TransientSolution(2).Thickness),...
144- (md.results.TransientSolution(2).MaskIceLevelset),...
145- (md.results.TransientSolution(2).CalvingMeltingrate),...
146- (md.results.TransientSolution(10).Vx),...
147- (md.results.TransientSolution(10).Vy),...
148- (md.results.TransientSolution(10).Vel),...
149- (md.results.TransientSolution(10).Pressure),...
150- (md.results.TransientSolution(10).Base),...
151- (md.results.TransientSolution(10).Surface),...
152- (md.results.TransientSolution(10).Thickness),...
153- (md.results.TransientSolution(10).MaskIceLevelset),...
154- (md.results.TransientSolution(10).CalvingMeltingrate),...
155- };
156+ (md.results.TransientSolution(1).BasalforcingsFloatingiceMeltingRate),...
157+ (md.results.TransientSolution(5).Vx),...
158+ (md.results.TransientSolution(5).Vy),...
159+ (md.results.TransientSolution(5).Vel),...
160+ (md.results.TransientSolution(5).Thickness),...
161+ (md.results.TransientSolution(5).MaskIceLevelset),...
162+ (md.results.TransientSolution(5).CalvingCalvingrate),...
163+ (md.results.TransientSolution(5).CalvingMeltingrate),...
164+ (md.results.TransientSolution(5).BasalforcingsFloatingiceMeltingRate),...
165+ (md.results.TransientSolution(10).Vx),...
166+ (md.results.TransientSolution(10).Vy),...
167+ (md.results.TransientSolution(10).Vel),...
168+ (md.results.TransientSolution(10).Thickness),...
169+ (md.results.TransientSolution(10).MaskIceLevelset),...
170+ (md.results.TransientSolution(10).CalvingCalvingrate),...
171+ (md.results.TransientSolution(10).CalvingMeltingrate),...
172+ (md.results.TransientSolution(10).BasalforcingsFloatingiceMeltingRate),...
173+ };
174Index: ../trunk-jpl/test/NightlyRun/test543.py
175===================================================================
176--- ../trunk-jpl/test/NightlyRun/test543.py (revision 26664)
177+++ ../trunk-jpl/test/NightlyRun/test543.py (revision 26665)
178@@ -1,7 +1,6 @@
179-#Test Name: PigTranFrontalforcingsrignotautoregression
180+#Test Name: PigTranStochasticforcing
181 import numpy as np
182
183-from calvingvonmises import *
184 from frontalforcingsrignotautoregression import *
185 from socket import gethostname
186 from model import *
187@@ -18,25 +17,42 @@
188 md = setflowequation(md, 'SSA', 'all')
189 md.timestepping.start_time = 0
190 md.timestepping.time_step = 1
191-md.timestepping.final_time = 25
192+md.timestepping.final_time = 10
193
194-# Separate domain in 2 basins
195-idbasin = np.zeros((md.mesh.numberofelements,))
196+# Basin separation TF
197+idb_tf = np.zeros((md.mesh.numberofelements,))
198 iid1 = np.where(md.mesh.x <= -1.6e6)[0]
199 for ii in range(md.mesh.numberofelements):
200 for vertex in range(3):
201 if md.mesh.elements[ii][vertex] - 1 in iid1: # one vertex in basin 1; NOTE: offset because of 1-based vertex indexing
202- idbasin[ii] = 1
203+ idb_tf[ii] = 1
204 if idbasin[ii] == 0: # no vertex was found in basin 1
205 for vertex in range(3):
206- idbasin[ii] = 2
207+ idb_tf[ii] = 2
208+# Basin separation default
209+idb_df = np.zeros((md.mesh.numberofelements,))
210+iid1 = np.where(md.mesh.x <= -1.62e6)[0]
211+for ii in range(md.mesh.numberofelements):
212+ for vertex in range(3):
213+ if md.mesh.elements[ii][vertex] - 1 in iid1: # one vertex in basin 1; NOTE: offset because of 1-based vertex indexing
214+ idb_df[ii] = 1
215+ if idbasin[ii] == 0: # no vertex was found in basin 1
216+ for vertex in range(3):
217+ idb_df[ii] = 2
218+#Dimensionalities
219+nb_tf = 2
220+nb_clv = 2
221+nb_flmlt = 2
222
223-# Calving and frontalforcings parameters
224+# Calving parameters
225 md.mask.ice_levelset = 1e4 * (md.mask.ice_levelset + 0.5)
226-md.calving = calvingvonmises()
227+md.calving.calvingrate = 0.3*np.ones((md.mesh.numberofvertices,))
228+md.levelset.spclevelset = np.full((md.mesh.numberofvertices,), np.nan)
229+md.levelset.migration_max = 10.0
230+#Frontal forcing parameters
231 md.frontalforcings = frontalforcingsrignotautoregression()
232-md.frontalforcings.num_basins = 2
233-md.frontalforcings.basin_id = idbasin - 1 # NOTE: offset because of 1-based vertex indexing
234+md.frontalforcings.num_basins = nb_tf
235+md.frontalforcings.basin_id = idb_tf - 1 # NOTE: offset because of 1-based vertex indexing
236 md.frontalforcings.subglacial_discharge = 0.1 * np.ones((md.mesh.numberofvertices,))
237 md.frontalforcings.beta0 = np.array([[0.05, 0.01]]) # intercept values of TF in basins [C]
238 md.frontalforcings.beta1 = np.array([[0.001, 0.0001]]) # trend values of TF in basins [C/yr]
239@@ -44,60 +60,69 @@
240 md.frontalforcings.ar_order = 4
241 md.frontalforcings.ar_timestep = 2 # timestep of the autoregressive model [yr]
242 md.frontalforcings.phi = np.array([[0.1, -0.1, 0.01, -0.01], [0.2, -0.2, 0.1, 0.0]]) # autoregressive parameters
243+#Floating Ice Melt parameters
244+md.basalforcings.floatingice_melting_rate = 0.1*np.ones((md.mesh.numberofvertices,))
245
246+#Covariance matrix
247+covtf = 1e-4*np.identity(nb_tf)
248+covclv = 1e-1*np.identity(nb_clv)
249+covclv[0,0] = 1/10*covclv[0,0]
250+covflmlt = 0.05*np.identity(nb_flmlt)
251+covglob = np.zeros([6,6])
252+covglob[0:2,0:2] = covtf
253+covglob[2:4,2:4] = covclv
254+covglob[4:6,4:6] = covflmlt
255+
256 # Stochastic forcing
257 md.stochasticforcing.isstochasticforcing = 1
258-md.stochasticforcing.fields = ['FrontalForcingsRignotAutoregression']
259+md.stochasticforcing.fields = ['FrontalForcingsRignotAutoregression','DefaultCalving','FloatingMeltRate']
260+md.stochasticforcing.defauldimension = 2
261+md.stochasticforcing.default_id = idb_df
262+md.stochasticforcing.covariance = covglob # global covariance among- and between-fields
263+md.stochasticforcing.randomflag = 0 # determines true/false randomness
264
265-md.stochasticforcing.covariance = 1e-4 * np.array([[1.5, 0.5], [0.5, 0.4]]) # global covariance among- and between-fields
266-md.stochasticforcing.randomflag = 0 # determines true/false randomness
267+md.transient.ismovingfront = 1
268+md.transient.isgroundingline = 1
269
270-
271-md.transient.ismovingfront = 1
272-md.levelset.spclevelset = np.full((md.mesh.numberofvertices,), np.nan)
273-md.levelset.migration_max = 1e10
274-
275-md.transient.requested_outputs = ['default', 'CalvingMeltingrate']
276+md.transient.requested_outputs = ['default', 'CalvingCalvingrate','CalvingMeltingrate','BasalforcingsFloatingiceMeltingRate']
277 md.cluster = generic('name', gethostname(), 'np', 2)
278 md = solve(md, 'Transient')
279
280 # Fields and tolerances to track changes
281 field_names = [
282- 'Vx1', 'Vy1', 'Vel1', 'Pressure1', 'Bed1', 'Surface1', 'Thickness1', 'MaskIceLevelset1', 'CalvingMetlingRate1',
283- 'Vx2', 'Vy2', 'Vel2', 'Pressure2', 'Bed2', 'Surface2', 'Thickness2', 'MaskIceLevelset2', 'CalvingMetlingRate2',
284- 'Vx10', 'Vy10', 'Vel10', 'Pressure10', 'Bed10', 'Surface10', 'Thickness10', 'MaskIceLevelset10', 'CalvingMetlingRate10'
285+ 'Vx1' ,'Vy1' ,'Vel1' ,'Thickness1' ,'MaskIceLevelset1' ,'CalvingCalvingrate1' ,'CalvingMeltingrate1' ,'BasalforcingsFloatingiceMeltingRate1',
286+ 'Vx2' ,'Vy2' ,'Vel2' ,'Thickness2' ,'MaskIceLevelset2' ,'CalvingCalvingrate1' ,'CalvingMeltingrate2' ,'BasalforcingsFloatingiceMeltingRate2',
287+ 'Vx10','Vy10','Vel10','Thickness10','MaskIceLevelset10','CalvingCalvingrate10','CalvingMeltingrate10','BasalforcingsFloatingiceMeltingRate10'
288 ]
289+
290 field_tolerances = [
291- 1e-11, 2e-11, 2e-11, 1e-12, 2e-11, 6e-12, 9e-12, 1e-12, 1e-9,
292- 2e-11, 1e-11, 1e-11, 9e-12, 2e-1, 2e-11, 2e-11, 1e-11, 1e-9,
293- 2e-10, 1e-10, 1e-10, 1e-10, 1e-10, 1e-10, 1e-10, 1e-9, 1e-9
294+ 1e-11,2e-11,2e-11,1e-11,1e-9,1e-10,1e-10,1e-10,
295+ 2e-11,1e-11,1e-11,9e-11,2e-9,1e-10,1e-10,1e-10,
296+ 2e-10,1e-10,1e-10,1e-10,5e-9,1e-10,1e-10,1e-10
297 ]
298 field_values = [
299 md.results.TransientSolution[0].Vx,
300 md.results.TransientSolution[0].Vy,
301 md.results.TransientSolution[0].Vel,
302- md.results.TransientSolution[0].Pressure,
303- md.results.TransientSolution[0].Base,
304- md.results.TransientSolution[0].Surface,
305 md.results.TransientSolution[0].Thickness,
306 md.results.TransientSolution[0].MaskIceLevelset,
307+ md.results.TransientSolution[0].CalvingCalvingrate,
308 md.results.TransientSolution[0].CalvingMeltingrate,
309- md.results.TransientSolution[1].Vx,
310- md.results.TransientSolution[1].Vy,
311- md.results.TransientSolution[1].Vel,
312- md.results.TransientSolution[1].Pressure,
313- md.results.TransientSolution[1].Base,
314- md.results.TransientSolution[1].Surface,
315- md.results.TransientSolution[1].Thickness,
316- md.results.TransientSolution[1].MaskIceLevelset,
317- md.results.TransientSolution[1].CalvingMeltingrate,
318+ md.results.TransientSolution[0].BasalforcingsFloatingiceMeltingRate,
319+ md.results.TransientSolution[4].Vx,
320+ md.results.TransientSolution[4].Vy,
321+ md.results.TransientSolution[4].Vel,
322+ md.results.TransientSolution[4].Thickness,
323+ md.results.TransientSolution[4].MaskIceLevelset,
324+ md.results.TransientSolution[4].CalvingCalvingrate,
325+ md.results.TransientSolution[4].CalvingMeltingrate,
326+ md.results.TransientSolution[4].BasalforcingsFloatingiceMeltingRate,
327 md.results.TransientSolution[9].Vx,
328 md.results.TransientSolution[9].Vy,
329 md.results.TransientSolution[9].Vel,
330- md.results.TransientSolution[9].Pressure,
331- md.results.TransientSolution[9].Base,
332- md.results.TransientSolution[9].Surface,
333 md.results.TransientSolution[9].Thickness,
334 md.results.TransientSolution[9].MaskIceLevelset,
335- md.results.TransientSolution[9].CalvingMeltingrate
336+ md.results.TransientSolution[9].CalvingCalvingrate,
337+ md.results.TransientSolution[9].CalvingMeltingrate,
338+ md.results.TransientSolution[9].BasalforcingsFloatingiceMeltingRate
339 ]
Note: See TracBrowser for help on using the repository browser.